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
46 #include "gmx_fatal.h"
49 int ifunc_index(directive d,int type)
76 gmx_fatal(FARGS,"Invalid bond type %d",type);
88 return F_CROSS_BOND_BONDS;
90 return F_CROSS_BOND_ANGLES;
92 return F_UREY_BRADLEY;
94 return F_QUARTIC_ANGLES;
98 return F_LINEAR_ANGLES;
100 gmx_fatal(FARGS,"Invalid angle type %d",type);
106 if (type == 1 || (d == d_pairtypes && type == 2))
111 gmx_fatal(FARGS,"Invalid pairs type %d",type);
114 return F_LJC_PAIRS_NB;
116 case d_dihedraltypes:
131 return F_PDIHS; /* proper dihedrals where we allow multiple terms over single bond */
133 gmx_fatal(FARGS,"Invalid dihedral type %d",type);
140 case d_nonbond_params:
158 gmx_fatal(FARGS,"Invalid vsites3 type %d",type);
168 gmx_fatal(FARGS,"Invalid vsites4 type %d",type);
174 case d_constrainttypes:
181 gmx_fatal(FARGS,"Invalid constraints type %d",type);
186 case d_position_restraints:
191 gmx_fatal(FARGS,"Water polarization should now be listed under [ water_polarization ]\n");
194 gmx_fatal(FARGS,"Invalid position restraint type %d",type);
200 return F_POLARIZATION;
204 gmx_fatal(FARGS,"Invalid polarization type %d",type);
207 case d_thole_polarization:
209 case d_water_polarization:
211 case d_angle_restraints:
213 case d_angle_restraints_z:
215 case d_distance_restraints:
217 case d_orientation_restraints:
219 case d_dihedral_restraints:
222 gmx_fatal(FARGS,"invalid directive %s in ifunc_index (%s:%s)",
223 dir2str(d),__FILE__,__LINE__);
228 const char *dir2str (directive d)
236 directive str2dir (char *dstr)
239 char buf[STRLEN],*ptr;
241 /* Hack to be able to read old topologies */
242 if (gmx_strncasecmp_min(dstr,"dummies",7) == 0) {
243 sprintf(buf,"virtual_sites%s",dstr+7);
249 for (d=0; (d<d_maxdir); d++)
250 if (gmx_strcasecmp_min(ptr,dir2str((directive)d)) == 0)
256 static directive **necessary = NULL;
258 static void set_nec(directive **n, ...)
259 /* Must always have at least one extra argument */
267 d = (directive)va_arg(ap,int);
270 } while (d != d_none);
274 void DS_Init(DirStack **DS)
276 if (necessary==NULL) {
279 snew(necessary,d_maxdir);
280 set_nec(&(necessary[d_defaults]),d_none);
281 set_nec(&(necessary[d_atomtypes]),d_defaults,d_none);
282 set_nec(&(necessary[d_bondtypes]),d_atomtypes,d_none);
283 set_nec(&(necessary[d_constrainttypes]),d_atomtypes,d_none);
284 set_nec(&(necessary[d_pairtypes]),d_atomtypes,d_none);
285 set_nec(&(necessary[d_angletypes]),d_atomtypes,d_none);
286 set_nec(&(necessary[d_dihedraltypes]),d_atomtypes,d_none);
287 set_nec(&(necessary[d_nonbond_params]),d_atomtypes,d_none);
288 set_nec(&(necessary[d_implicit_genborn_params]),d_atomtypes,d_none);
289 set_nec(&(necessary[d_implicit_surface_params]),d_atomtypes,d_none);
290 set_nec(&(necessary[d_cmaptypes]),d_atomtypes,d_none);
291 set_nec(&(necessary[d_moleculetype]),d_atomtypes,d_none);
292 set_nec(&(necessary[d_atoms]),d_moleculetype,d_none);
293 set_nec(&(necessary[d_vsites2]),d_atoms,d_none);
294 set_nec(&(necessary[d_vsites3]),d_atoms,d_none);
295 set_nec(&(necessary[d_vsites4]),d_atoms,d_none);
296 set_nec(&(necessary[d_vsitesn]),d_atoms,d_none);
297 set_nec(&(necessary[d_bonds]),d_atoms,d_none);
298 set_nec(&(necessary[d_exclusions]),d_bonds,d_constraints,d_settles,d_none);
299 set_nec(&(necessary[d_pairs]),d_atoms,d_none);
300 set_nec(&(necessary[d_pairs_nb]),d_atoms,d_none);
301 set_nec(&(necessary[d_angles]),d_atoms,d_none);
302 set_nec(&(necessary[d_polarization]),d_atoms,d_none);
303 set_nec(&(necessary[d_water_polarization]),d_atoms,d_none);
304 set_nec(&(necessary[d_thole_polarization]),d_atoms,d_none);
305 set_nec(&(necessary[d_dihedrals]),d_atoms,d_none);
306 set_nec(&(necessary[d_constraints]),d_atoms,d_none);
307 set_nec(&(necessary[d_settles]),d_atoms,d_none);
308 set_nec(&(necessary[d_system]),d_moleculetype,d_none);
309 set_nec(&(necessary[d_molecules]),d_system,d_none);
310 set_nec(&(necessary[d_position_restraints]),d_atoms,d_none);
311 set_nec(&(necessary[d_angle_restraints]),d_atoms,d_none);
312 set_nec(&(necessary[d_angle_restraints_z]),d_atoms,d_none);
313 set_nec(&(necessary[d_distance_restraints]),d_atoms,d_none);
314 set_nec(&(necessary[d_orientation_restraints]),d_atoms,d_none);
315 set_nec(&(necessary[d_dihedral_restraints]),d_atoms,d_none);
316 set_nec(&(necessary[d_cmap]), d_atoms, d_none);
318 for(i=0; (i<d_maxdir); i++) {
320 fprintf(debug,"%20s: ",dir2str((directive)i));
327 fprintf(debug,"%20s ",dir2str(d));
338 void DS_Done (DirStack **DS)
342 while (*DS != NULL) {
349 void DS_Push (DirStack **DS, directive d)
359 int DS_Search(DirStack *DS, directive d)
364 while ((D != NULL) && (D->d != d))
370 int DS_Check_Order(DirStack *DS,directive d)
375 /* Check if parameter definitions appear after a moleculetype directive */
376 if (d<d_moleculetype && DS_Search(DS,d_moleculetype))
379 /* Check if all the necessary directives have appeared before directive d */
380 if (necessary[d][0] == d_none)
384 d0=necessary[d][i++];
385 if (DS_Search(DS,d0))