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)
78 gmx_fatal(FARGS, "Invalid bond type %d", type);
91 return F_CROSS_BOND_BONDS;
93 return F_CROSS_BOND_ANGLES;
95 return F_UREY_BRADLEY;
97 return F_QUARTIC_ANGLES;
101 return F_LINEAR_ANGLES;
103 gmx_fatal(FARGS, "Invalid angle type %d", type);
109 if (type == 1 || (d == d_pairtypes && type == 2))
119 gmx_fatal(FARGS, "Invalid pairs type %d", type);
123 return F_LJC_PAIRS_NB;
125 case d_dihedraltypes:
141 return F_PDIHS; /* proper dihedrals where we allow multiple terms over single bond */
143 gmx_fatal(FARGS, "Invalid dihedral type %d", type);
150 case d_nonbond_params:
173 gmx_fatal(FARGS, "Invalid vsites3 type %d", type);
184 gmx_fatal(FARGS, "Invalid vsites4 type %d", type);
190 case d_constrainttypes:
198 gmx_fatal(FARGS, "Invalid constraints type %d", type);
203 case d_position_restraints:
212 gmx_fatal(FARGS, "Invalid position restraint type %d", type);
219 return F_POLARIZATION;
223 gmx_fatal(FARGS, "Invalid polarization type %d", type);
226 case d_thole_polarization:
228 case d_water_polarization:
230 case d_angle_restraints:
232 case d_angle_restraints_z:
234 case d_distance_restraints:
236 case d_orientation_restraints:
238 case d_dihedral_restraints:
241 gmx_fatal(FARGS, "invalid directive %s in ifunc_index (%s:%s)",
242 dir2str(d), __FILE__, __LINE__);
247 const char *dir2str (directive d)
259 directive str2dir (char *dstr)
262 char buf[STRLEN], *ptr;
264 /* Hack to be able to read old topologies */
265 if (gmx_strncasecmp_min(dstr, "dummies", 7) == 0)
267 sprintf(buf, "virtual_sites%s", dstr+7);
275 for (d = 0; (d < d_maxdir); d++)
277 if (gmx_strcasecmp_min(ptr, dir2str((directive)d)) == 0)
286 static directive **necessary = NULL;
288 static void set_nec(directive **n, ...)
289 /* Must always have at least one extra argument */
298 d = (directive)va_arg(ap, int);
306 void DS_Init(DirStack **DS)
308 if (necessary == NULL)
312 snew(necessary, d_maxdir);
313 set_nec(&(necessary[d_defaults]), d_none);
314 set_nec(&(necessary[d_atomtypes]), d_defaults, d_none);
315 set_nec(&(necessary[d_bondtypes]), d_atomtypes, d_none);
316 set_nec(&(necessary[d_constrainttypes]), d_atomtypes, d_none);
317 set_nec(&(necessary[d_pairtypes]), d_atomtypes, d_none);
318 set_nec(&(necessary[d_angletypes]), d_atomtypes, d_none);
319 set_nec(&(necessary[d_dihedraltypes]), d_atomtypes, d_none);
320 set_nec(&(necessary[d_nonbond_params]), d_atomtypes, d_none);
321 set_nec(&(necessary[d_implicit_genborn_params]), d_atomtypes, d_none);
322 set_nec(&(necessary[d_implicit_surface_params]), d_atomtypes, d_none);
323 set_nec(&(necessary[d_cmaptypes]), d_atomtypes, d_none);
324 set_nec(&(necessary[d_moleculetype]), d_atomtypes, d_none);
325 set_nec(&(necessary[d_atoms]), d_moleculetype, d_none);
326 set_nec(&(necessary[d_vsites2]), d_atoms, d_none);
327 set_nec(&(necessary[d_vsites3]), d_atoms, d_none);
328 set_nec(&(necessary[d_vsites4]), d_atoms, d_none);
329 set_nec(&(necessary[d_vsitesn]), d_atoms, d_none);
330 set_nec(&(necessary[d_bonds]), d_atoms, d_none);
331 set_nec(&(necessary[d_exclusions]), d_bonds, d_constraints, d_settles, d_none);
332 set_nec(&(necessary[d_pairs]), d_atoms, d_none);
333 set_nec(&(necessary[d_pairs_nb]), d_atoms, d_none);
334 set_nec(&(necessary[d_angles]), d_atoms, d_none);
335 set_nec(&(necessary[d_polarization]), d_atoms, d_none);
336 set_nec(&(necessary[d_water_polarization]), d_atoms, d_none);
337 set_nec(&(necessary[d_thole_polarization]), d_atoms, d_none);
338 set_nec(&(necessary[d_dihedrals]), d_atoms, d_none);
339 set_nec(&(necessary[d_constraints]), d_atoms, d_none);
340 set_nec(&(necessary[d_settles]), d_atoms, d_none);
341 set_nec(&(necessary[d_system]), d_moleculetype, d_none);
342 set_nec(&(necessary[d_molecules]), d_system, d_none);
343 set_nec(&(necessary[d_position_restraints]), d_atoms, d_none);
344 set_nec(&(necessary[d_angle_restraints]), d_atoms, d_none);
345 set_nec(&(necessary[d_angle_restraints_z]), d_atoms, d_none);
346 set_nec(&(necessary[d_distance_restraints]), d_atoms, d_none);
347 set_nec(&(necessary[d_orientation_restraints]), d_atoms, d_none);
348 set_nec(&(necessary[d_dihedral_restraints]), d_atoms, d_none);
349 set_nec(&(necessary[d_cmap]), d_atoms, d_none);
351 for (i = 0; (i < d_maxdir); i++)
355 fprintf(debug, "%20s: ", dir2str((directive)i));
363 d = necessary[i][j++];
366 fprintf(debug, "%20s ", dir2str(d));
373 fprintf(debug, "\n");
381 void DS_Done (DirStack **DS)
393 void DS_Push (DirStack **DS, directive d)
403 int DS_Search(DirStack *DS, directive d)
408 while ((D != NULL) && (D->d != d))
416 int DS_Check_Order(DirStack *DS, directive d)
421 /* Check if parameter definitions appear after a moleculetype directive */
422 if (d < d_moleculetype && DS_Search(DS, d_moleculetype))
427 /* Check if all the necessary directives have appeared before directive d */
428 if (necessary[d][0] == d_none)
436 d0 = necessary[d][i++];
437 if (DS_Search(DS, d0))
442 while (d0 != d_none);