3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
10 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
11 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
12 * Copyright (c) 2001-2004, The GROMACS development team,
13 * check out http://www.gromacs.org for more information.
15 * This program is free software; you can redistribute it and/or
16 * modify it under the terms of the GNU General Public License
17 * as published by the Free Software Foundation; either version 2
18 * of the License, or (at your option) any later version.
20 * If you want to redistribute modifications, please consider that
21 * scientific software is very special. Version control is crucial -
22 * bugs must be traceable. We will be happy to consider code for
23 * inclusion in the official distribution, but derived work must not
24 * be called official GROMACS. Details are found in the README & COPYING
25 * files - if they are missing, get the official version at www.gromacs.org.
27 * To help us fund GROMACS development, we humbly ask that you cite
28 * the papers on the package - you can find them in the top README file.
30 * For more info, check our website at http://www.gromacs.org
33 * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
35 /* This file is completely threadsafe - keep it that way! */
46 #include "hackblock.h"
48 #include "gmx_fatal.h"
50 static void clear_atom_list(int i0, atom_id a[])
54 for (i=i0; i < MAXATOMLIST; i++)
58 static void clear_force_param(int i0, real c[])
62 for (i=i0; i < MAXFORCEPARAM; i++)
66 void add_param(t_params *ps,int ai,int aj, real *c, char *s)
70 if ((ai < 0) || (aj < 0))
71 gmx_fatal(FARGS,"Trying to add impossible atoms: ai=%d, aj=%d",ai,aj);
73 ps->param[ps->nr].AI=ai;
74 ps->param[ps->nr].AJ=aj;
75 clear_atom_list(2, ps->param[ps->nr].a);
77 clear_force_param(0, ps->param[ps->nr].c);
79 for(i=0; (i < MAXFORCEPARAM); i++)
80 ps->param[ps->nr].c[i]=c[i];
81 set_p_string(&(ps->param[ps->nr]),s);
85 void add_imp_param(t_params *ps,int ai,int aj,int ak,int al,real c0, real c1,
89 ps->param[ps->nr].AI=ai;
90 ps->param[ps->nr].AJ=aj;
91 ps->param[ps->nr].AK=ak;
92 ps->param[ps->nr].AL=al;
93 clear_atom_list (4, ps->param[ps->nr].a);
94 ps->param[ps->nr].C0=c0;
95 ps->param[ps->nr].C1=c1;
96 clear_force_param(2, ps->param[ps->nr].c);
97 set_p_string(&(ps->param[ps->nr]),s);
101 void add_dih_param(t_params *ps,int ai,int aj,int ak,int al,real c0, real c1,
105 ps->param[ps->nr].AI=ai;
106 ps->param[ps->nr].AJ=aj;
107 ps->param[ps->nr].AK=ak;
108 ps->param[ps->nr].AL=al;
109 clear_atom_list (4, ps->param[ps->nr].a);
110 ps->param[ps->nr].C0=c0;
111 ps->param[ps->nr].C1=c1;
112 ps->param[ps->nr].C2=c2;
113 clear_force_param(3, ps->param[ps->nr].c);
114 set_p_string(&(ps->param[ps->nr]),s);
118 void add_cmap_param(t_params *ps, int ai, int aj, int ak, int al, int am, char *s)
121 ps->param[ps->nr].AI=ai;
122 ps->param[ps->nr].AJ=aj;
123 ps->param[ps->nr].AK=ak;
124 ps->param[ps->nr].AL=al;
125 ps->param[ps->nr].AM=am;
126 clear_atom_list(5,ps->param[ps->nr].a);
127 clear_force_param(0,ps->param[ps->nr].c);
128 set_p_string(&(ps->param[ps->nr]),s);
132 void add_vsite2_atoms(t_params *ps,int ai,int aj,int ak)
135 ps->param[ps->nr].AI=ai;
136 ps->param[ps->nr].AJ=aj;
137 ps->param[ps->nr].AK=ak;
138 clear_atom_list (3, ps->param[ps->nr].a);
139 clear_force_param(0, ps->param[ps->nr].c);
140 set_p_string(&(ps->param[ps->nr]),"");
144 void add_vsite2_param(t_params *ps,int ai,int aj,int ak,real c0)
147 ps->param[ps->nr].AI=ai;
148 ps->param[ps->nr].AJ=aj;
149 ps->param[ps->nr].AK=ak;
150 clear_atom_list (3, ps->param[ps->nr].a);
151 ps->param[ps->nr].C0=c0;
152 clear_force_param(1, ps->param[ps->nr].c);
153 set_p_string(&(ps->param[ps->nr]),"");
157 void add_vsite3_param(t_params *ps,int ai,int aj,int ak,int al,
161 ps->param[ps->nr].AI=ai;
162 ps->param[ps->nr].AJ=aj;
163 ps->param[ps->nr].AK=ak;
164 ps->param[ps->nr].AL=al;
165 clear_atom_list (4, ps->param[ps->nr].a);
166 ps->param[ps->nr].C0=c0;
167 ps->param[ps->nr].C1=c1;
168 clear_force_param(2, ps->param[ps->nr].c);
169 set_p_string(&(ps->param[ps->nr]),"");
173 void add_vsite3_atoms(t_params *ps,int ai,int aj,int ak,int al, gmx_bool bSwapParity)
176 ps->param[ps->nr].AI=ai;
177 ps->param[ps->nr].AJ=aj;
178 ps->param[ps->nr].AK=ak;
179 ps->param[ps->nr].AL=al;
180 clear_atom_list (4, ps->param[ps->nr].a);
181 clear_force_param(0, ps->param[ps->nr].c);
183 ps->param[ps->nr].C1=-1;
184 set_p_string(&(ps->param[ps->nr]),"");
188 void add_vsite4_atoms(t_params *ps,int ai,int aj,int ak,int al,int am)
191 ps->param[ps->nr].AI=ai;
192 ps->param[ps->nr].AJ=aj;
193 ps->param[ps->nr].AK=ak;
194 ps->param[ps->nr].AL=al;
195 ps->param[ps->nr].AM=am;
196 clear_atom_list (5, ps->param[ps->nr].a);
197 clear_force_param(0, ps->param[ps->nr].c);
198 set_p_string(&(ps->param[ps->nr]),"");
202 int search_jtype(t_restp *rtp,char *name,gmx_bool bNterm)
204 int niter,iter,j,k,kmax,jmax,minstrlen;
205 char *rtpname,searchname[12];
207 strcpy(searchname,name);
209 /* Do a best match comparison */
210 /* for protein N-terminus, allow renaming of H1, H2 and H3 to H */
211 if ( bNterm && (strlen(searchname)==2) && (searchname[0] == 'H') &&
212 ( (searchname[1] == '1') || (searchname[1] == '2') ||
213 (searchname[1] == '3') ) ) {
220 for(iter=0; (iter<niter && jmax==-1); iter++) {
222 /* Try without the hydrogen number in the N-terminus */
223 searchname[1] = '\0';
225 for(j=0; (j<rtp->natom); j++) {
226 rtpname=*(rtp->atomname[j]);
227 if (gmx_strcasecmp(searchname,rtpname) == 0) {
229 kmax=strlen(searchname);
232 if (iter == niter - 1) {
233 minstrlen = min(strlen(searchname), strlen(rtpname));
234 for(k=0; k < minstrlen; k++)
235 if (searchname[k] != rtpname[k])
245 gmx_fatal(FARGS,"Atom %s not found in rtp database in residue %s",
246 searchname,rtp->resname);
247 if (kmax != strlen(searchname))
248 gmx_fatal(FARGS,"Atom %s not found in rtp database in residue %s, "
249 "it looks a bit like %s",
250 searchname,rtp->resname,*(rtp->atomname[jmax]));