7 |Gromacs| must know on which atoms and combinations of atoms the various
8 contributions to the potential functions (see chapter :ref:`ff`) must act.
9 It must also know what parameters must be applied to the various
10 functions. All this is described in the *topology* file :ref:`top`, which
11 lists the *constant attributes* of each atom. There are many more atom
12 types than elements, but only atom types present in biological systems
13 are parameterized in the force field, plus some metals, ions and
14 silicon. The bonded and special interactions are determined by fixed
15 lists that are included in the topology file. Certain non-bonded
16 interactions must be excluded (first and second neighbors), as these are
17 already treated in bonded interactions. In addition, there are *dynamic
18 attributes* of atoms - their positions, velocities and forces. These do
19 not strictly belong to the molecular topology, and are stored in the
20 coordinate file :ref:`gro` (positions and velocities), or
21 trajectory file :ref:`trr` (positions, velocities, forces).
23 This chapter describes the setup of the topology file, the :ref:`top` file and
24 the database files: what the parameters stand for and how/where to
25 change them if needed. First, all file formats are explained. Section
26 :ref:`fffiles` describes the organization of the files in each force
29 **Note:** if you construct your own topologies, we encourage you to
30 upload them to our topology archive at our `webpage`_! Just imagine how thankful
31 you’d have been if your topology had been available there before you
32 started. The same goes for new force fields or modified versions of the
33 standard force fields - contribute them to the force field archive!
35 .. _homepage: `webpage`_
40 In |Gromacs|, there are three types of
42 , see :numref:`Table %s <tab-ptype>`. Only regular atoms and virtual
43 interaction sites are used in |Gromacs|; shells are necessary for
44 polarizable models like the Shell-Water models \ :ref:`45 <refMaaren2001a>`.
48 .. table:: Particle types in |Gromacs|
50 +--------------+----------+
52 +==============+==========+
54 +--------------+----------+
56 +--------------+----------+
57 | virtual side | V (or D) |
58 +--------------+----------+
66 Each force field defines a set of atom
68 which have a characteristic name or number, and mass (in a.m.u.). These
69 listings are found in the ``atomtypes.atp`` file (:ref:`atp` =
70 **a**\ tom **t**\ ype **p**\ arameter file). Therefore, it is in this
71 file that you can begin to change and/or add an atom type. A sample from
72 the ``gromos43a1.ff`` force field is listed below.
76 | O 15.99940 ; carbonyl oxygen (C=O)
77 | OM 15.99940 ; carboxyl oxygen (CO-)
78 | OA 15.99940 ; hydroxyl, sugar or ester oxygen
79 | OW 15.99940 ; water oxygen
80 | N 14.00670 ; peptide nitrogen (N or NH)
81 | NT 14.00670 ; terminal nitrogen (NH2)
82 | NL 14.00670 ; terminal nitrogen (NH3)
83 | NR 14.00670 ; aromatic nitrogen
84 | NZ 14.00670 ; Arg NH (NH2)
85 | NE 14.00670 ; Arg NE (NH)
86 | C 12.01100 ; bare carbon
87 |CH1 13.01900 ; aliphatic or sugar CH-group
88 |CH2 14.02700 ; aliphatic or sugar CH2-group
89 |CH3 15.03500 ; aliphatic CH3-group
91 **Note:** |Gromacs| makes use of the atom types as a name, *not* as a
92 number (as *e.g.* in GROMOS).
99 Some force fields use virtual interaction sites (interaction sites that
100 are constructed from other particle positions) on which certain
101 interactions are located (*e.g.* on benzene rings, to reproduce the
102 correct quadrupole). This is described in sec. :ref:`virtualsites`.
104 To make virtual sites in your system, you should include a section
105 ``[ virtual_sites? ]`` (for backward compatibility the old
106 name ``[ dummies? ]`` can also be used) in your topology
107 file, where the ``?`` stands for the number constructing
108 particles for the virtual site. This will be `:ref:`2`` for
109 type 2, `:ref:`3`` for types 3, 3fd, 3fad and 3out and
110 `:ref:`4`` for type 4fdn. The last of these replace an older
111 4fd type (with the ‘type’ value 1) that could occasionally be unstable;
112 while it is still supported internally in the code, the old 4fd type
113 should not be used in new input files. The different types are explained
114 in sec. :ref:`virtualsites`.
116 Parameters for type 2 should look like this:
124 for type 3 like this:
129 ; Site from funct a b
130 5 1 2 3 1 0.7439756 0.128012
132 for type 3fd like this:
137 ; Site from funct a d
140 for type 3fad like this:
145 ; Site from funct theta d
148 for type 3out like this:
153 ; Site from funct a b c
154 5 1 2 3 4 -0.4 -0.4 6.9281
156 for type 4fdn like this:
161 ; Site from funct a b c
162 5 1 2 3 4 2 1.0 0.9 0.105
164 This will result in the construction of a virtual site, number 5 (first
165 column ``Site``), based on the positions of the atoms
166 whose indices are 1 and 2 or 1, 2 and 3 or 1, 2, 3 and 4 (next two,
167 three or four columns ``from``) following the rules
168 determined by the function number (next column ``funct``)
169 with the parameters specified (last one, two or three columns
170 ``a b . .``). Obviously, the atom numbers (including
171 virtual site number) depend on the molecule. It may be instructive to
172 study the topologies for TIP4P or TIP5P water models that are included
173 with the |Gromacs| distribution.
175 **Note** that if any constant bonded interactions are defined between
176 virtual sites and/or normal atoms, they will be removed by
177 :ref:`grompp <gmx grompp>` (unless the option ``-normvsbds`` is used). This
178 removal of bonded interactions is done after generating exclusions, as
179 the generation of exclusions is based on “chemically” bonded
182 Virtual sites can be constructed in a more generic way using basic
183 geometric parameters. The directive that can be used is ``[ virtual_sitesn ]``. Required
184 parameters are listed in :numref:`Table %s <tab-topfile2>`. An example entry for
185 defining a virtual site at the center of geometry of a given set of
200 The *static* properties (see :numref:`Table %s <tab-statprop>`) assigned to the atom
201 types are assigned based on data in several places. The mass is listed
202 in ``atomtypes.atp`` (see :ref:`atomtype`), whereas the charge is listed
203 in :ref:`rtp` (:ref:`rtp` = **r**\ esidue **t**\ opology **p**\ arameter file,
204 see :ref:`rtp`). This implies that the charges are only defined in the
205 building blocks of amino acids, nucleic acids or otherwise, as defined
206 by the user. When generating a :ref:`topology <top>` using the :ref:`pdb2gmx <gmx pdb2gmx>`
207 program, the information from these files is combined.
211 .. table:: Static atom type properties in |Gromacs|
213 +----------+------------------+----------+
214 | Property | Symbol | Unit |
215 +==========+==================+==========+
217 +----------+------------------+----------+
218 | Mass | m | a.m.u. |
219 +----------+------------------+----------+
220 | Charge | q | electron |
221 +----------+------------------+----------+
222 | epsilon | :math:`\epsilon` | kJ/mol |
223 +----------+------------------+----------+
224 | sigma | :math:`\sigma` | nm |
225 +----------+------------------+----------+
230 Non-bonded parameters
231 ~~~~~~~~~~~~~~~~~~~~~
233 The non-bonded parameters consist of the van der Waals parameters V (``c6``
234 or :math:`\sigma`, depending on the combination rule) and W (``c12`` or
235 :math:`\epsilon`), as listed in the file ``ffnonbonded.itp``, where ``ptype`` is
236 the particle type (see :numref:`Table %s <tab-ptype>`). As with the bonded
237 parameters, entries in ``[ *type ]`` directives are applied to their counterparts in
238 the topology file. Missing parameters generate warnings, except as noted
239 below in section :ref:`pairinteractions`.
244 ;name at.num mass charge ptype V(c6) W(c12)
245 O 8 15.99940 0.000 A 0.22617E-02 0.74158E-06
246 OM 8 15.99940 0.000 A 0.22617E-02 0.74158E-06
250 ; i j func V(c6) W(c12)
251 O O 1 0.22617E-02 0.74158E-06
252 O OA 1 0.22617E-02 0.13807E-05
255 **Note** that most of the included force fields also include the ``at.num.``
256 column, but this same information is implied in the OPLS-AA ``bond_type``
257 column. The interpretation of the parameters V and W depends on the
258 combination rule that was chosen in the ``[ defaults ]`` section of the topology file
259 (see :ref:`topfile`):
264 \mbox{for combination rule 1}: & &
265 \begin{array}{llllll}
266 \mbox{V}_{ii} & = & C^{(6)}_{i} & = & 4\,\epsilon_i\sigma_i^{6} &
267 \mbox{[ kJ mol$^{-1}$ nm$^{6}$ ]}\\
268 \mbox{W}_{ii} & = & C^{(12)}_{i} & = & 4\,\epsilon_i\sigma_i^{12} &
269 \mbox{[ kJ mol$^{-1}$ nm$^{12}$ ]}\\
272 \mbox{for combination rules 2 and 3}: & &
274 \mbox{V}_{ii} & = & \sigma_i & \mbox{[ nm ]} \\
275 \mbox{W}_{ii} & = & \epsilon_i & \mbox{[ kJ mol$^{-1}$ ]}
276 \end{array}\end{aligned}
278 Some or all combinations for different atom types can be given in the
279 ``[ nonbond_params ]`` section, again with parameters V and
280 W as defined above. Any combination that is not given will be computed
281 from the parameters for the corresponding atom types, according to the
287 \mbox{for combination rules 1 and 3}: & &
289 C^{(6)}_{ij} & = & \left(C^{(6)}_i\,C^{(6)}_j\right)^{\frac{1}{2}} \\
290 C^{(12)}_{ij} & = & \left(C^{(12)}_i\,C^{(12)}_j\right)^{\frac{1}{2}}
293 \mbox{for combination rule 2}: & &
295 \sigma_{ij} & = & \frac{1}{2}(\sigma_i+\sigma_j) \\
296 \epsilon_{ij} & = & \sqrt{\epsilon_i\,\epsilon_j}
297 \end{array}\end{aligned}
299 When :math:`\sigma` and :math:`\epsilon` need to be supplied (rules 2
300 and 3), it would seem it is impossible to have a non-zero :math:`C^{12}`
301 combined with a zero :math:`C^6` parameter. However, providing a
302 negative :math:`\sigma` will do exactly that, such that :math:`C^6` is
303 set to zero and :math:`C^{12}` is calculated normally. This situation
304 represents a special case in reading the value of :math:`\sigma`, and
307 There is only one set of combination rules for Buckingham potentials:
312 A_{ij} &=& \left(A_{ii} \, A_{jj}\right)^{1/2} \\
313 B_{ij} &=& 2 / \left(\frac{1}{B_{ii}} + \frac{1}{B_{jj}}\right) \\
314 C_{ij} &=& \left(C_{ii} \, C_{jj}\right)^{1/2}
322 (*i.e.* bonds, bond angles, improper and proper dihedrals) are listed in
323 ``ffbonded.itp``. The entries in this database describe,
324 respectively, the atom types in the interactions, the type of the
325 interaction, and the parameters associated with that interaction. These
326 parameters are then read by
327 :ref:`grompp <gmx grompp>` when processing a
328 topology and applied to the relevant bonded parameters, *i.e.*
329 ``bondtypes`` are applied to entries in the
330 ``[ bonds ]`` directive, etc. Any bonded parameter that is
331 missing from the relevant :``[ *type ]`` directive generates
332 a fatal error. The types of interactions are listed in
333 :numref:`Table %s <tab-topfile2>`. Example excerpts from such files
340 C O 1 0.12300 502080.
341 C OM 1 0.12500 418400.
346 HO OA C 1 109.500 397.480
347 HO OA CH1 1 109.500 397.480
352 NR5* NR5 2 0.000 167.360
353 NR5* NR5* 2 0.000 167.360
357 ; j k func phi0 cp mult
358 C OA 1 180.000 16.736 2
359 C N 1 180.000 33.472 2
364 ; Ryckaert-Bellemans Dihedrals
367 CP2 CP2 3 9.2789 12.156 -13.120 -3.0597 26.240 -31.495
369 In the ``ffbonded.itp`` file, you can add bonded parameters.
370 If you want to include parameters for new atom types, make sure you
371 define them in ``atomtypes.atp`` as well.
373 For most interaction types, bonded parameters are searched and assigned
374 using an exact match for all type names and allowing only a single set
375 of parameters. The exception to this rule are
378 ``[ dihedraltypes ]`` wildcard atom type names can be
379 specified with the letter ``X`` in one or more of the four
380 positions. Thus one can for example assign proper dihedral parameters
381 based on the types of the middle two atoms. The parameters for the entry
382 with the most exact matches, i.e. the least wildcard matches, will be
383 used. Note that |Gromacs| versions older than 5.1.3 used the first match,
384 which means that a full match would be ignored if it is preceded by an
385 entry that matches on wildcards. Thus it is suggested to put wildcard
386 entries at the end, in case someone might use a forcefield with older
387 versions of |Gromacs|. In addition there is a dihedral type 9 which adds
388 the possibility of assigning multiple dihedral potentials, useful for
389 combining terms with different multiplicities. The different dihedral
390 potential parameter sets should be on directly adjacent lines in the
391 ``[ dihedraltypes ]`` section.
399 An organizational structure that usually corresponds to molecules is the
400 ``[ moleculetype ]`` entry. This entry serves two main
401 purposes. One is to give structure to the topology file(s), usually
402 corresponding to real molecules. This makes the topology easier to read
403 and writing it less labor intensive. A second purpose is computational
404 efficiency. The system definition that is kept in memory is proportional
405 in size of the ``moleculetype`` definitions. If a molecule
406 is present in 100000 copies, this saves a factor of 100000 in memory,
407 which means the system usually fits in cache, which can improve
408 performance tremendously. Interactions that correspond to chemical
409 bonds, that generate exclusions, can only be defined between atoms
410 within a ``moleculetype``. It is allowed to have multiple
411 molecules which are not covalently bonded in one
412 ``moleculetype`` definition. Molecules can be made
413 infinitely long by connecting to themselves over periodic boundaries.
414 When such periodic molecules are present, an option in the
415 :ref:`mdp` file needs to be set to tell |Gromacs| not to attempt
416 to make molecules that are broken over periodic boundaries whole again.
418 Intermolecular interactions
419 ~~~~~~~~~~~~~~~~~~~~~~~~~~~
421 In some cases, one would like atoms in different molecules to also
422 interact with other interactions than the usual non-bonded interactions.
423 This is often the case in binding studies. When the molecules are
424 covalently bound, e.g. a ligand binding covalently to a protein, they
425 are effectively one molecule and they should be defined in one
426 ``[ moleculetype ]`` entry. Note that
427 :ref:`pdb2gmx <gmx pdb2gmx>` has an option to put two or more molecules in
428 one ``[ moleculetype ]`` entry. When molecules are not
429 covalently bound, it is much more convenient to use separate
430 ``moleculetype`` definitions and specify the intermolecular
431 interactions in the ``[ intermolecular_interactions]``
432 section. In this section, which is placed at the end of the topology
433 (see :numref:`Table %s <tab-topfile1>`), normal bonded interactions
434 can be specified using global atom indices. The only restrictions are
435 that no interactions can be used that generates exclusions and no
436 constraints can be used.
438 .. _pairinteractions:
440 Intramolecular pair interactions
441 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
443 Extra Lennard-Jones and electrostatic interactions between pairs of
444 atoms in a molecule can be added in the ``[ pairs ]`` section of a molecule
445 definition. The parameters for these interactions can be set
446 independently from the non-bonded interaction parameters. In the GROMOS
447 force fields, pairs are only used to modify the 1-4 interactions
448 (interactions of atoms separated by three bonds). In these force fields
449 the 1-4 interactions are excluded from the non-bonded interactions (see
456 ; i j func cs6 cs12 ; THESE ARE 1-4 INTERACTIONS
457 O O 1 0.22617E-02 0.74158E-06
458 O OM 1 0.22617E-02 0.74158E-06
461 The pair interaction parameters for the atom types in ``ffnonbonded.itp``
462 are listed in the ``[ pairtypes ]`` section. The GROMOS force fields list all these
463 interaction parameters explicitly, but this section might be empty for
464 force fields like OPLS that calculate the 1-4 interactions by uniformly
465 scaling the parameters. Pair parameters that are not present in the ``[ pairtypes ]``
466 section are only generated when ``gen-pairs`` is set to ``yes`` in the
467 ``[ defaults ]`` directive of ``forcefield.itp`` (see :ref:`topfile`). When ``gen-pairs`` is
468 set to ``no``, :ref:`grompp <gmx grompp>` will give a warning for each pair type for which no
469 parameters are given.
471 The normal pair interactions, intended for 1-4 interactions, have
472 function type 1. Function type 2 and the ``[ pairs_nb ]`` are intended for free-energy
473 simulations. When determining hydration free energies, the solute needs
474 to be decoupled from the solvent. This can be done by adding a B-state
475 topology (see sec. :ref:`fecalc`) that uses zero for all solute
476 non-bonded parameters, *i.e.* charges and LJ parameters. However, the
477 free energy difference between the A and B states is not the total
478 hydration free energy. One has to add the free energy for reintroducing
479 the internal Coulomb and LJ interactions in the solute when in vacuum.
480 This second step can be combined with the first step when the Coulomb
481 and LJ interactions within the solute are not modified. For this
482 purpose, there is a pairs function type 2, which is identical to
483 function type 1, except that the B-state parameters are always identical
484 to the A-state parameters. For searching the parameters in the ``[ pairtypes ]`` section,
485 no distinction is made between function type 1 and 2. The pairs section
486 ``[ pairs_nb ]`` is intended to replace the non-bonded interaction. It uses the unscaled
487 charges and the non-bonded LJ parameters; it also only uses the A-state
488 parameters. **Note** that one should add exclusions for all atom pairs
489 listed in ``[ pairs_nb ]``, otherwise such pairs will also end up in the normal neighbor
492 Alternatively, this same behavior can be achieved without ever touching
493 the topology, by using the ``couple-moltype``, ``couple-lambda0``,
494 ``couple-lambda1``, and ``couple-intramol`` keywords. See sections
495 sec. :ref:`fecalc` and sec. :ref:`dgimplement` for more information.
497 All three pair types always use plain Coulomb interactions, even when
498 Reaction-field, PME, Ewald or shifted Coulomb interactions are selected
499 for the non-bonded interactions. Energies for types 1 and 2 are written
500 to the energy and log file in separate “LJ-14” and “Coulomb-14” entries
501 per energy group pair. Energies for ``[ pairs_nb ]`` are added to the “LJ-(SR)” and
502 “Coulomb-(SR)” terms.
509 The exclusions for non-bonded interactions are generated by :ref:`grompp <gmx grompp>` for
510 neighboring atoms up to a certain number of bonds away, as defined in
511 the ``[ moleculetype ]`` section in the topology file (see :ref:`topfile`). Particles are
512 considered bonded when they are connected by “chemical” bonds (``[ bonds ]`` types 1
513 to 5, 7 or 8) or constraints (``[ constraints ]`` type 1). Type 5 ``[ bonds ]`` can be used to create a
514 connection between two atoms without creating an interaction. There is a
515 harmonic interaction (``[ bonds ]`` type 6) that does not connect the atoms by a
516 chemical bond. There is also a second constraint type (``[ constraints ]`` type 2) that
517 fixes the distance, but does not connect the atoms by a chemical bond.
518 For a complete list of all these interactions, see :numref:`Table %s <tab-topfile2>`.
520 Extra exclusions within a molecule can be added manually in a
521 ``[ exclusions ]`` section. Each line should start with one
522 atom index, followed by one or more atom indices. All non-bonded
523 interactions between the first atom and the other atoms will be
526 When all non-bonded interactions within or between groups of atoms need
527 to be excluded, is it more convenient and much more efficient to use
528 energy monitor group exclusions (see sec. :ref:`groupconcept`).
532 Constraint algorithms
533 ---------------------
535 Constraints are defined in the ``[ constraints ]`` section. The format is two atom numbers
536 followed by the function type, which can be 1 or 2, and the constraint
537 distance. The only difference between the two types is that type 1 is
538 used for generating exclusions and type 2 is not (see sec. :ref:`excl`).
539 The distances are constrained using the LINCS or the SHAKE algorithm,
540 which can be selected in the :ref:`mdp` file. Both types of constraints can be
541 perturbed in free-energy calculations by adding a second constraint
542 distance (see :ref:`constraintforce`). Several types of bonds and
543 angles (see :numref:`Table %s <tab-topfile2>`) can be converted automatically to
544 constraints by :ref:`grompp <gmx grompp>`. There are several options for this in the :ref:`mdp`
547 We have also implemented the SETTLE
548 algorithm \ :ref:`47 <refMiyamoto92>`, which is an analytical solution of SHAKE, specifically for
549 water. SETTLE can be selected in the topology file. See, for instance,
550 the SPC molecule definition:
559 ; nr at type res nr ren nm at nm cg nr charge
560 1 OW 1 SOL OW1 1 -0.82
561 2 HW 1 SOL HW2 1 0.41
562 3 HW 1 SOL HW3 1 0.41
573 The ``[ settles ]`` directive defines the first atom of the
574 water molecule. The settle funct is always 1, and the distance between
575 O-H and H-H distances must be given. **Note** that the algorithm can
576 also be used for TIP3P and TIP4P \ :ref:`128 <refJorgensen83>`. TIP3P just has
577 another geometry. TIP4P has a virtual site, but since that is generated
578 it does not need to be shaken (nor stirred).
582 :ref:`pdb2gmx <gmx pdb2gmx>` input files
583 ----------------------------------------
585 The |Gromacs| program :ref:`pdb2gmx <gmx pdb2gmx>` generates a topology for the input
586 coordinate file. Several formats are supported for that coordinate file,
587 but :ref:`pdb` is the most commonly-used format (hence the name :ref:`pdb2gmx <gmx pdb2gmx>`).
588 :ref:`pdb2gmx <gmx pdb2gmx>` searches for force fields in sub-directories of the |Gromacs|
589 ``share/top`` directory and your working directory. Force fields are
590 recognized from the file ``forcefield.itp`` in a directory with the
591 extension ``.ff``. The file ``forcefield.doc`` may be present, and if so, its
592 first line will be used by :ref:`pdb2gmx <gmx pdb2gmx>` to present a short description to the
593 user to help in choosing a force field. Otherwise, the user can choose a
594 force field with the ``-ff xxx`` command-line argument to :ref:`pdb2gmx <gmx pdb2gmx>`, which
595 indicates that a force field in a ``xxx.ff`` directory is desired. :ref:`pdb2gmx <gmx pdb2gmx>`
596 will search first in the working directory, then in the |Gromacs|
597 ``share/top`` directory, and use the first matching ``xxx.ff`` directory found.
599 Two general files are read by :ref:`pdb2gmx <gmx pdb2gmx>`: an atom type file (extension
600 :ref:`atp`, see :ref:`atomtype`) from the force-field directory, and a file
601 called ``residuetypes.dat`` from either the working directory, or the
602 |Gromacs| ``share/top`` directory. ``residuetypes.dat`` determines which residue
603 names are considered protein, DNA, RNA, water, and ions.
605 :ref:`pdb2gmx <gmx pdb2gmx>` can read one or multiple databases with topological information
606 for different types of molecules. A set of files belonging to one
607 database should have the same basename, preferably telling something
608 about the type of molecules (*e.g.* aminoacids, rna, dna). The possible
613 - ``<basename>.r2b (optional)``
615 - ``<basename>.arn (optional)``
617 - ``<basename>.hdb (optional)``
619 - ``<basename>.n.tdb (optional)``
621 - ``<basename>.c.tdb (optional)``
623 Only the :ref:`rtp` file, which contains the topologies of the building
624 blocks, is mandatory. Information from other files will only be used for
625 building blocks that come from an :ref:`rtp` file with the same base name. The
626 user can add building blocks to a force field by having additional files
627 with the same base name in their working directory. By default, only
628 extra building blocks can be defined, but calling :ref:`pdb2gmx <gmx pdb2gmx>` with the ``-rtpo``
629 option will allow building blocks in a local file to replace the default
630 ones in the force field.
635 The files holding the residue databases have the extension :ref:`rtp`.
636 Originally this file contained building blocks (amino acids) for
637 proteins, and is the |Gromacs| interpretation of the ``rt37c4.dat`` file of
638 GROMOS. So the residue database file contains information (bonds,
639 charges, charge groups, and improper dihedrals) for a frequently-used
640 building block. It is better *not* to change this file because it is
641 standard input for :ref:`pdb2gmx <gmx pdb2gmx>`, but if changes are needed make them in the
642 :ref:`top` file (see :ref:`topfile`), or in a :ref:`rtp` file in the working
643 directory as explained in sec. :ref:`pdb2gmxfiles`. Defining topologies
644 of new small molecules is probably easier by writing an include topology
645 file :ref:`itp` directly. This will be discussed in section :ref:`molitp`.
646 When adding a new protein residue to the database, don’t forget to add
647 the residue name to the residuetypes.dat file, so that :ref:`grompp <gmx grompp>`, :ref:`make_ndx <gmx make_ndx>`
648 and analysis tools can recognize the residue as a protein residue (see
649 :ref:`defaultgroups`).
651 The :ref:`rtp` files are only used by :ref:`pdb2gmx <gmx pdb2gmx>`. As mentioned before, the only
652 extra information this program needs from the :ref:`rtp` database is bonds,
653 charges of atoms, charge groups, and improper dihedrals, because the
654 rest is read from the coordinate input file. Some proteins contain
655 residues that are not standard, but are listed in the coordinate file.
656 You have to construct a building block for this “strange” residue,
657 otherwise you will not obtain a :ref:`top` file. This also holds for molecules
658 in the coordinate file such as ligands, polyatomic ions, crystallization
659 co-solvents, etc. The residue database is constructed in the following
664 [ bondedtypes ] ; mandatory
665 ; bonds angles dihedrals impropers
670 [ atoms ] ; mandatory
671 ; name type charge chargegroup
686 [ exclusions ] ; optional
689 [ angles ] ; optional
690 ;atom1 atom2 atom3 th0 cth
692 [ dihedrals ] ; optional
693 ;atom1 atom2 atom3 atom4 phi0 cp mult
695 [ impropers ] ; optional
696 ;atom1 atom2 atom3 atom4 q0 cq
705 The file is free format; the only restriction is that there can be at
706 most one entry on a line. The first field in the file is the ``[ bondedtypes ]`` field,
707 which is followed by four numbers, indicating the interaction type for
708 bonds, angles, dihedrals, and improper dihedrals. The file contains
709 residue entries, which consist of atoms and (optionally) bonds, angles,
710 dihedrals, and impropers. The charge group codes denote the charge group
711 numbers. Atoms in the same charge group should always be ordered
712 consecutively. When using the hydrogen database with :ref:`pdb2gmx <gmx pdb2gmx>` for adding
713 missing hydrogens (see :ref:`hdb`), the atom names defined in the :ref:`rtp`
714 entry should correspond exactly to the naming convention used in the
715 hydrogen database. The atom names in the bonded interaction can be
716 preceded by a minus or a plus, indicating that the atom is in the
717 preceding or following residue respectively. Explicit parameters added
718 to bonds, angles, dihedrals, and impropers override the standard
719 parameters in the :ref:`itp` files. This should only be used in special cases.
720 Instead of parameters, a string can be added for each bonded
721 interaction. This is used in GROMOS-96 :ref:`rtp` files. These strings are
722 copied to the topology file and can be replaced by force-field
723 parameters by the C-preprocessor in :ref:`grompp <gmx grompp>` using ``#define`` statements.
725 :ref:`pdb2gmx <gmx pdb2gmx>` automatically generates all angles. This means
726 that for most force fields the ``[ angles ]`` field is only
727 useful for overriding :ref:`itp` parameters. For the GROMOS-96
728 force field the interaction number of all angles needs to be specified.
730 :ref:`pdb2gmx <gmx pdb2gmx>` automatically generates one proper dihedral for every rotatable
731 bond, preferably on heavy atoms. When the ``[ dihedrals ]`` field is used, no other
732 dihedrals will be generated for the bonds corresponding to the specified
733 dihedrals. It is possible to put more than one dihedral function on a
734 rotatable bond. In the case of CHARMM27 FF :ref:`pdb2gmx <gmx pdb2gmx>` can add correction
735 maps to the dihedrals using the default ``-cmap`` option. Please refer to
736 :ref:`charmmff` for more information.
738 :ref:`pdb2gmx <gmx pdb2gmx>` sets the number of exclusions to 3, which means
739 that interactions between atoms connected by at most 3 bonds are
740 excluded. Pair interactions are generated for all pairs of atoms that
741 are separated by 3 bonds (except pairs of hydrogens). When more
742 interactions need to be excluded, or some pair interactions should not
743 be generated, an ``[ exclusions ]`` field can be added,
744 followed by pairs of atom names on separate lines. All non-bonded and
745 pair interactions between these atoms will be excluded.
747 Residue to building block database
748 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
750 Each force field has its own naming convention for residues. Most
751 residues have consistent naming, but some, especially those with
752 different protonation states, can have many different names. The
753 :ref:`r2b` files are used to convert standard residue names to
754 the force-field build block names. If no :ref:`r2b` is present
755 in the force-field directory or a residue is not listed, the building
756 block name is assumed to be identical to the residue name. The
757 :ref:`r2b` can contain 2 or 5 columns. The 2-column format has
758 the residue name in the first column and the building block name in the
759 second. The 5-column format has 3 additional columns with the building
760 block for the residue occurring in the N-terminus, C-terminus and both
761 termini at the same time (single residue molecule). This is useful for,
762 for instance, the AMBER force fields. If one or more of the terminal
763 versions are not present, a dash should be entered in the corresponding
766 There is a |Gromacs| naming convention for residues which is only apparent
767 (except for the :ref:`pdb2gmx <gmx pdb2gmx>` code) through the
768 :ref:`r2b` file and ``specbond.dat`` files. This
769 convention is only of importance when you are adding residue types to an
770 :ref:`rtp` file. The convention is listed in :numref:`Table %s <tab-r2b>`.
771 For special bonds with, for instance,
772 a heme group, the |Gromacs| naming convention is introduced through
773 ``specbond.dat`` (see :ref:`specbond`),
774 which can subsequently be translated by the :ref:`r2b` file,
777 .. |NDEL| replace:: N\ :math:`_\delta`
778 .. |NEPS| replace:: N\ :math:`_\epsilon`
782 .. table:: Internal |Gromacs| residue naming convention.
784 +--------------+-----------------------------------------------------------+
785 | |Gromacs| ID | Residue |
786 +==============+===========================================================+
787 | ARG | protonated arginine |
788 +--------------+-----------------------------------------------------------+
789 | ARGN | neutral arginine |
790 +--------------+-----------------------------------------------------------+
791 | ASP | negatively charged aspartic acid |
792 +--------------+-----------------------------------------------------------+
793 | ASPH | neutral aspartic acid |
794 +--------------+-----------------------------------------------------------+
795 | CYS | neutral cysteine |
796 +--------------+-----------------------------------------------------------+
797 | CYS2 | cysteine with sulfur bound to another cysteine or a heme |
798 +--------------+-----------------------------------------------------------+
799 | GLU | negatively charged glutamic acid |
800 +--------------+-----------------------------------------------------------+
801 | GLUH | neutral glutamic acid |
802 +--------------+------------------------------+----------------------------+
803 | HISD | neutral histidine with |NDEL| protonated |
804 +--------------+-----------------------------------------------------------+
805 | HISE | neutral histidine with |NEPS| protonated |
806 +--------------+------------------------------+----------------------------+
807 | HISH | positive histidine with both |NDEL| and |NEPS| protonated |
808 +--------------+-----------------------------------------------------------+
809 | HIS1 | histidine bound to a heme |
810 +--------------+-----------------------------------------------------------+
811 | LYSN | neutral lysine |
812 +--------------+-----------------------------------------------------------+
813 | LYS | protonated lysine |
814 +--------------+-----------------------------------------------------------+
816 +--------------+-----------------------------------------------------------+
819 Atom renaming database
820 ~~~~~~~~~~~~~~~~~~~~~~
822 Force fields often use atom names that do not follow IUPAC or PDB
823 convention. The :ref:`arn` database is used to translate the
824 atom names in the coordinate file to the force-field names. Atoms that
825 are not listed keep their names. The file has three columns: the
826 building block name, the old atom name, and the new atom name,
827 respectively. The residue name supports question-mark wildcards that
828 match a single character.
830 An additional general atom renaming file called
831 ``xlateat.dat`` is present in the ``share/top``
832 directory, which translates common non-standard atom names in the
833 coordinate file to IUPAC/PDB convention. Thus, when writing force-field
834 files, you can assume standard atom names and no further atom name
835 translation is required, except for translating from standard atom names
836 to the force-field ones.
841 The hydrogen database is stored in :ref:`hdb` files. It contains information
842 for the :ref:`pdb2gmx <gmx pdb2gmx>` program on how to connect hydrogen atoms to existing
843 atoms. In versions of the database before |Gromacs| 3.3, hydrogen atoms
844 were named after the atom they are connected to: the first letter of the
845 atom name was replaced by an ‘H.’ In the versions from 3.3 onwards, the
846 H atom has to be listed explicitly, because the old behavior was
847 protein-specific and hence could not be generalized to other molecules.
848 If more than one hydrogen atom is connected to the same atom, a number
849 will be added to the end of the hydrogen atom name. For example, adding
850 two hydrogen atoms to ``ND2`` (in asparagine), the hydrogen atoms will
851 be named ``HD21`` and ``HD22``. This is important since atom naming in
852 the :ref:`rtp` file (see :ref:`rtp`) must be the same. The format of the
853 hydrogen database is as follows:
867 On the first line we see the residue name (ALA or ARG) and the number of
868 kinds of hydrogen atoms that may be added to this residue by the
869 hydrogen database. After that follows one line for each addition, on
872 - The number of H atoms added
874 - The method for adding H atoms, which can be any of:
876 #. | *one planar hydrogen*, *e.g.* *rings or peptide bond*
877 | One hydrogen atom (n) is generated, lying in the plane of atoms
878 (i,j,k) on the plane bisecting angle (j-i-k) at a distance of
879 0.1 nm from atom i, such that the angles (n-i-j) and (n-i-k) are
880 :math:`>` 90\ :math:`^{\rm o}`.
882 #. | *one single hydrogen*, *e.g.* *hydroxyl*
883 | One hydrogen atom (n) is generated at a distance of 0.1 nm from
884 atom i, such that angle (n-i-j)=109.5 degrees and dihedral
887 #. | *two planar hydrogens*, *e.g.* *ethylene -C=CH*:math:`_2`, *or amide
889 | Two hydrogens (n1,n2) are generated at a distance of 0.1 nm from
890 atom i, such that angle (n1-i-j)=(n2-i-j)=120 degrees and
891 dihedral (n1-i-j-k)=cis and (n2-i-j-k)=trans, such that names
892 are according to IUPAC standards \ :ref:`129 <refiupac70>`.
894 #. | *two or three tetrahedral hydrogens*, *e.g.* *-CH*:math:`_3`
895 | Three (n1,n2,n3) or two (n1,n2) hydrogens are generated at a
896 distance of 0.1 nm from atom i, such that angle
897 (n1-i-j)=(n2-i-j)=(n3-i-j)=109.47:math:`^{\rm o}`, dihedral
898 (n1-i-j-k)=trans, (n2-i-j-k)=trans+120 and
899 (n3-i-j-k)=trans+240:math:`^{\rm o}`.
901 #. | *one tetrahedral hydrogen*, *e.g.* *C*\ :math:`_3`\* CH*
902 | One hydrogen atom (n:math:`^\prime`) is generated at a distance
903 of 0.1 nm from atom i in tetrahedral conformation such that
905 (n:math:`^\prime`-i-j)=(n:math:`^\prime`-i-k)=(n:math:`^\prime`-i-l)=109.47:math:`^{\rm o}`.
907 #. | *two tetrahedral hydrogens*, *e.g.* *C-CH*\ :math:`_2`\*-C*
908 | Two hydrogen atoms (n1,n2) are generated at a distance of 0.1 nm
909 from atom i in tetrahedral conformation on the plane bisecting
910 angle j-i-k with angle
911 (n1-i-n2)=(n1-i-j)=(n1-i-k)=109.47:math:`^{\rm o}`.
913 #. | *two water hydrogens*
914 | Two hydrogens are generated around atom i according to
915 SPC \ :ref:`80 <refBerendsen81>` water geometry. The symmetry
916 axis will alternate between three coordinate axes in both
919 #. | *three water “hydrogens”*
920 | Two hydrogens are generated around atom i according to
921 SPC \ :ref:`80 <refBerendsen81>` water geometry. The symmetry
922 axis will alternate between three coordinate axes in both
923 directions. In addition, an extra particle is generated on the
924 position of the oxygen with the first letter of the name
925 replaced by ‘M’. This is for use with four-atom water models
926 such as TIP4P \ :ref:`128 <refJorgensen83>`.
928 #. | *four water “hydrogens”*
929 | Same as above, except that two additional particles are
930 generated on the position of the oxygen, with names ‘LP1’ and
931 ‘LP2.’ This is for use with five-atom water models such as
932 TIP5P \ :ref:`130 <refMahoney2000a>`.
934 - The name of the new H atom (or its prefix, *e.g.* ``HD2``
935 for the asparagine example given earlier).
937 - Three or four control atoms (i,j,k,l), where the first always is the
938 atom to which the H atoms are connected. The other two or three
939 depend on the code selected. For water, there is only one control
942 Some more exotic cases can be approximately constructed from the above
943 tools, and with suitable use of energy minimization are good enough for
944 beginning MD simulations. For example secondary amine hydrogen, nitrenyl
945 hydrogen (:math:`\mathrm{C}=\mathrm{NH}`)
946 and even ethynyl hydrogen could be approximately constructed using
947 method 2 above for hydroxyl hydrogen.
954 are stored in ``aminoacids.n.tdb`` and
955 ``aminoacids.c.tdb`` for the N- and C-termini respectively.
956 They contain information for the :ref:`pdb2gmx <gmx pdb2gmx>` program on how
957 to connect new atoms to existing ones, which atoms should be removed or
958 changed, and which bonded interactions should be added. Their format is
959 as follows (from ``gromos43a1.ff/aminoacids.c.tdb``):
968 O O1 OM 15.9994 -0.635
969 OXT O2 OM 15.9994 -0.635
985 The file is organized in blocks, each with a header specifying the name
986 of the block. These blocks correspond to different types of termini that
987 can be added to a molecule. In this example ``[ COO- ]`` is
988 the first block, corresponding to changing the terminal carbon atom into
989 a deprotonated carboxyl group. ``[ None ]`` is the second
990 terminus type, corresponding to a terminus that leaves the molecule as
991 it is. Block names cannot be any of the following:
992 ``replace``, ``add``, ``delete``,
993 ``bonds``, ``angles``,
994 ``dihedrals``, ``impropers``. Doing so would
995 interfere with the parameters of the block, and would probably also be
996 very confusing to human readers.
998 For each block the following options are present:
1001 | Replace an existing atom by one with a different atom type, atom
1002 name, charge, and/or mass. This entry can be used to replace an
1003 atom that is present both in the input coordinates and in the
1004 :ref:`rtp` database, but also to only rename an atom in
1005 the input coordinates such that it matches the name in the force
1006 field. In the latter case, there should also be a corresponding
1007 ``[ add ]`` section present that gives instructions to
1008 add the same atom, such that the position in the sequence and the
1009 bonding is known. Such an atom can be present in the input
1010 coordinates and kept, or not present and constructed by
1011 :ref:`pdb2gmx <gmx pdb2gmx>`. For each atom to be replaced on line
1012 should be entered with the following fields:
1014 - name of the atom to be replaced
1016 - new atom name (optional)
1025 | Add new atoms. For each (group of) added atom(s), a two-line entry
1026 is necessary. The first line contains the same fields as an entry
1027 in the hydrogen database (name of the new atom, number of atoms,
1028 type of addition, control atoms, see :ref:`hdb`), but the
1029 possible types of addition are extended by two more, specifically
1030 for C-terminal additions:
1032 #. | *two carboxyl oxygens, -COO*:math:`^-`
1033 | Two oxygens (n1,n2) are generated according to rule 3, at a
1034 distance of 0.136 nm from atom i and an angle
1035 (n1-i-j)=(n2-i-j)=117 degrees
1037 #. | *carboxyl oxygens and hydrogen, -COOH*
1038 | Two oxygens (n1,n2) are generated according to rule 3, at
1039 distances of 0.123 nm and 0.125 nm from atom i for n1 and n2,
1040 respectively, and angles (n1-i-j)=121 and (n2-i-j)=115 degrees.
1041 One hydrogen (n:math:`^\prime`) is generated around n2 according
1042 to rule 2, where n-i-j and n-i-j-k should be read as
1043 n\ :math:`^\prime`-n2-i and n\ :math:`^\prime`-n2-i-j,
1046 After this line, another line follows that specifies the details of
1047 the added atom(s), in the same way as for replacing atoms, *i.e.*:
1055 - charge group (optional)
1057 Like in the hydrogen database (see :ref:`rtp`), when more than one
1058 atom is connected to an existing one, a number will be appended to
1059 the end of the atom name. **Note** that, like in the hydrogen
1060 database, the atom name is now on the same line as the control atoms,
1061 whereas it was at the beginning of the second line prior to |Gromacs|
1062 version 3.3. When the charge group field is left out, the added atom
1063 will have the same charge group number as the atom that it is bonded
1067 | Delete existing atoms. One atom name per line.
1069 - | ``[ bonds ]``, ``[ angles ]``,
1070 ``[ dihedrals ]`` and ``[ impropers ]``
1071 | Add additional bonded parameters. The format is identical to that
1072 used in the :ref:`rtp` file, see :ref:`rtp`.
1074 Virtual site database
1075 ~~~~~~~~~~~~~~~~~~~~~
1077 Since we cannot rely on the positions of hydrogens in input files, we
1078 need a special input file to decide the geometries and parameters with
1079 which to add virtual site hydrogens. For more complex virtual site
1080 constructs (*e.g.* when entire aromatic side chains are made rigid) we
1081 also need information about the equilibrium bond lengths and angles for
1082 all atoms in the side chain. This information is specified in the
1083 :ref:`vsd` file for each force field. Just as for the termini,
1084 there is one such file for each class of residues in the
1087 The virtual site database is not really a very simple list of
1088 information. The first couple of sections specify which mass centers
1089 (typically called MCH\ :math:`_3`/MNH:math:`_3`) to use for
1090 CH\ :math:`_3`, NH\ :math:`_3`, and NH\ :math:`_2` groups. Depending on
1091 the equilibrium bond lengths and angles between the hydrogens and heavy
1092 atoms we need to apply slightly different constraint distances between
1093 these mass centers. **Note** that we do *not* have to specify the actual
1094 parameters (that is automatic), just the type of mass center to use. To
1095 accomplish this, there are three sections names ``[ CH3 ]``,
1096 ``[ NH3 ]``, and ``[ NH2 ]``. For each of these we expect three columns.
1097 The first column is the atom type bound to the 2/3 hydrogens, the second
1098 column is the next heavy atom type which this is bound, and the third
1099 column the type of mass center to use. As a special case, in the
1100 ``[ NH2 ]`` section it is also possible to specify ``planar`` in the
1101 second column, which will use a different construction without mass
1102 center. There are currently different opinions in some force fields
1103 whether an NH\ :math:`_2` group should be planar or not, but we try hard
1104 to stick to the default equilibrium parameters of the force field.
1106 The second part of the virtual site database contains explicit
1107 equilibrium bond lengths and angles for pairs/triplets of atoms in
1108 aromatic side chains. These entries are currently read by specific
1109 routines in the virtual site generation code, so if you would like to
1110 extend it *e.g.* to nucleic acids you would also need to write new code
1111 there. These sections are named after the short amino acid names
1112 (``[ PHE ]``, ``[ TYR ]``, ``[ TRP ]``, ``[ HID ]``, ``[ HIE ]``,
1113 ``[ HIP ]``), and simply contain 2 or 3 columns with atom names,
1114 followed by a number specifying the bond length (in nm) or angle (in
1115 degrees). **Note** that these are approximations of the equilibrated
1116 geometry for the entire molecule, which might not be identical to the
1117 equilibrium value for a single bond/angle if the molecule is strained.
1124 The primary mechanism used by
1125 :ref:`pdb2gmx <gmx pdb2gmx>` to generate
1126 inter-residue bonds relies on head-to-tail linking of backbone atoms in
1127 different residues to build a macromolecule. In some cases (*e.g.*
1128 disulfide bonds, a heme
1130 polymers), it is necessary to
1131 create inter-residue bonds that do not lie on the backbone. The file
1132 ``specbond.dat`` takes
1133 care of this function. It is necessary that the residues belong to the
1134 same ``[ moleculetype ]``. The ``-merge`` and
1135 ``-chainsep`` functions of :ref:`pdb2gmx <gmx pdb2gmx>` can be
1136 useful when managing special inter-residue bonds between different
1139 The first line of ``specbond.dat`` indicates the number of
1140 entries that are in the file. If you add a new entry, be sure to
1141 increment this number. The remaining lines in the file provide the
1142 specifications for creating bonds. The format of the lines is as
1145 ``resA atomA nbondsA resB atomB nbondsB length newresA
1148 The columns indicate:
1150 #. ``resA`` The name of residue A that participates in the
1153 #. ``atomA`` The name of the atom in residue A that forms
1156 #. ``nbondsA`` The total number of bonds
1159 #. ``resB`` The name of residue B that participates in the
1162 #. ``atomB`` The name of the atom in residue B that forms
1165 #. ``nbondsB`` The total number of bonds
1168 #. ``length`` The reference length for the bond. If
1169 ``atomA`` and ``atomB`` are not within
1170 ``length`` :math:`\pm` 10% in the coordinate file
1171 supplied to :ref:`pdb2gmx <gmx pdb2gmx>`, no bond will be formed.
1173 #. ``newresA`` The new name of residue A, if necessary. Some
1174 force fields use *e.g.* CYS2 for a cysteine in a disulfide or heme
1177 #. ``newresB`` The new name of residue B, likewise.
1187 The topology file is built following the |Gromacs| specification for a
1188 molecular topology. A :ref:`top` file can be generated by
1189 :ref:`pdb2gmx <gmx pdb2gmx>`. All possible entries in the topology file are
1190 listed in :numref:`Tables %s <tab-topfile1>` and
1191 :numref:`%s <tab-topfile2>`. Also tabulated are: all the units of
1192 the parameters, which interactions can be perturbed for free energy
1193 calculations, which bonded interactions are used by
1194 :ref:`grompp <gmx grompp>` for generating exclusions, and which bonded
1195 interactions can be converted to constraints by :ref:`grompp <gmx grompp>`.
1197 .. |VCR| replace:: V\ :math:`^{(cr)}`
1198 .. |WCR| replace:: W\ :math:`^{(cr)}`
1199 .. |CRO| replace:: :math:`^{(cr)}`
1200 .. |TREF| replace:: :numref:`Table %s <tab-topfile2>`
1201 .. |AKJM| replace:: :math:`a~\mathrm{kJ~mol}^{-1}`
1202 .. |KJN6| replace:: :math:`\mathrm{kJ~mol}^{-1}~\mathrm{nm}^{-6}`
1203 .. |BNM| replace:: :math:`b~\mathrm{nm}^{-1}`
1204 .. |C6LJ| replace:: :math:`c_6`
1205 .. |STAR| replace:: :math:`^{(*)}`
1206 .. |NREX| replace:: :math:`n_{ex}^{(nrexcl)}`
1207 .. |QEMU| replace:: :math:`q` (e); :math:`m` (u)
1208 .. |MQM| replace:: :math:`q,m`
1212 .. table:: The :ref:`topology <top>` file.
1214 +------------------------------------------------------------------------------------------------------------+
1216 +===================+===========================+=====+====+=========================================+=======+
1217 | interaction type | directive | # | f. | parameters | F. E. |
1218 | | | at. | tp | | |
1219 +-------------------+---------------------------+-----+----+-----------------------------------------+-------+
1220 | *mandatory* | ``defaults`` | non-bonded function type; |
1221 | | | combination rule\ |CRO|; |
1222 | | | generate pairs (no/yes); |
1223 | | | fudge LJ (); fudge QQ () |
1224 +-------------------+---------------------------+------------------------------------------------------------+
1225 | *mandatory* | ``atomtypes`` | atom type; m (u); q (e); particle type; |
1226 | | | |VCR| ; |WCR| |
1227 +-------------------+---------------------------+------------------------------------------------------------+
1228 | | ``bondtypes`` | (see |TREF|, directive ``bonds``) |
1230 | | ``pairtypes`` | (see |TREF|, directive ``pairs``) |
1232 | | ``angletypes`` | (see |TREF|, directive ``angles``) |
1234 | | ``dihedraltypes``\ |STAR| | (see |TREF|, directive ``dihedrals``) |
1236 | | ``constrainttypes`` | (see |TREF|, directive ``constraints``) |
1237 +-------------------+---------------------------+-----+----+-------------------------------------------------+
1238 | LJ | ``nonbond_params`` | 2 | 1 | |VCR| ; |WCR| |
1240 | Buckingham | ``nonbond_params`` | 2 | 2 | |AKJM| ; |BNM|; |
1241 | | | | | |C6LJ| (|KJN6|) |
1242 +-------------------+---------------------------+-----+----+-------------------------------------------------+
1246 +------------------------------------------------------------------------------------------------------------+
1247 | Molecule definition(s) |
1248 +===================+===========================+============================================================+
1249 | *mandatory* | ``moleculetype`` | molecule name; |NREX| |
1250 +-------------------+---------------------------+-----+----------------------------------------------+-------+
1251 | *mandatory* | ``atoms`` | 1 | atom type; residue number; | type |
1252 | | | | residue name; atom name; | |
1253 | | | | charge group number; |QEMU| | |MQM| |
1254 +-------------------+---------------------------+-----+----------------------------------------------+-------+
1255 | intra-molecular interaction and geometry definitions as described in |TREF| |
1256 +------------------------------------------------------------------------------------------------------------+
1260 +-------------+---------------+------------------------------------+
1262 +=============+===============+====================================+
1263 | *mandatory* | ``system`` | system name |
1264 +-------------+---------------+------------------------------------+
1265 | *mandatory* | ``molecules`` | molecule name; number of molecules |
1266 +-------------+---------------+------------------------------------+
1270 +------------------------------+----------------------------------------------------+
1271 | Inter-molecular interactions | |
1272 +==============================+====================================================+
1273 | *optional* | ``intermolecular_interactions`` |
1274 +------------------------------+----------------------------------------------------+
1275 | one or more bonded interactions as described in |TREF|, with two or more atoms, |
1276 | no interactions that generate exclusions, no constraints, use global atom numbers |
1277 +-----------------------------------------------------------------------------------+
1281 '\# at' is the required number of atom type indices for this directive
1282 'f. tp' is the value used to select this function type
1283 'F. E.' indicates which of the parameters can be interpolated in free energy calculations
1284 |CRO| the combination rule determines the type of LJ parameters, see
1285 |STAR| for ``dihedraltypes`` one can specify 4 atoms or the inner (outer for improper) 2 atoms
1286 |NREX| exclude neighbors :math:`n_{ex}` bonds away for non-bonded interactions
1287 For free energy calculations, type, :math:`q` and :math:`m` or no parameters should be added
1288 for topology 'B' (:math:`\lambda = 1`) on the same line, after the normal parameters.
1290 .. |BZERO| replace:: :math:`b_0`
1291 .. |KB| replace:: :math:`k_b`
1292 .. |KDR| replace:: :math:`k_{dr}`
1293 .. |NM2| replace:: (kJ mol\ :math:`^{-1}`\ nm\ :math:`^{-2}`
1294 .. |NM4| replace:: (kJ mol\ :math:`^{-1}`\ nm\ :math:`^{-4}`
1295 .. |DKJ| replace:: :math:`D` (kJ mol\ :math:`^{-1}`
1296 .. |BETA| replace:: :math:`\beta` (nm\ :math:`^{-1}`
1297 .. |C23| replace:: :math:`C_{i=2,3}` (kJ mol\ :math:`^{-1}\ nm\ :math:`^{-i}`
1298 .. |BMM| replace:: :math:`b_m`
1299 .. |GE0| replace:: :math:`\geq 0`
1300 .. |KO| replace:: :math:`k`
1301 .. |KJM| replace:: kJ mol\ :math:`^{-1}`
1302 .. |LUU| replace:: low, up\ :math:`_1`,\ :math:`_2`
1303 .. |MV| replace:: :math:`V`
1304 .. |MW| replace:: :math:`W`
1305 .. |QIJ| replace:: :math:`q_i`; :math:`q_j`
1306 .. |THE0| replace:: :math:`\theta_0`
1307 .. |KTHE| replace:: :math:`k_\theta`
1308 .. |KJR2| replace:: kJ mol\ :math:`^{-1}`\ rad\ :math:`^{-2}`
1309 .. |RN13| replace:: :math:`r_{13}`
1310 .. |KUB| replace:: :math:`k_{UB}`
1311 .. |C024| replace:: :math:`C_{i=0,1,2,3,4}`
1312 .. |KJRI| replace:: kJ mol\ :math:`^{-1}`\ rad\ :math:`^{-i}`
1313 .. |PHIS| replace:: :math:`\phi_s`
1314 .. |PHI0| replace:: :math:`\phi_0`
1315 .. |KPHI| replace:: :math:`k_\phi`
1316 .. |PHIK| replace:: :math:`\phi,k`
1317 .. |XI0| replace:: :math:`\xi_0`
1318 .. |KXI| replace:: :math:`k_\xi`
1319 .. |C0| replace:: :math:`C_0`
1320 .. |C1| replace:: :math:`C_1`
1321 .. |C2| replace:: :math:`C_2`
1322 .. |C3| replace:: :math:`C_3`
1323 .. |C4| replace:: :math:`C_4`
1324 .. |C5| replace:: :math:`C_5`
1325 .. |A0| replace:: :math:`a_0`
1326 .. |A1| replace:: :math:`a_1`
1327 .. |A2| replace:: :math:`a_2`
1328 .. |A3| replace:: :math:`a_3`
1329 .. |A4| replace:: :math:`a_4`
1330 .. |DOH| replace:: :math:d_{\mbox{\sc oh}}`
1331 .. |DHH| replace:: :math:d_{\mbox{\sc hh}}`
1332 .. |AO| replace:: :math:`a`
1333 .. |BO| replace:: :math:`b`
1334 .. |CO| replace:: :math:`c`
1335 .. |DO| replace:: :math:`d`
1336 .. |KX| replace:: :math:`k_{x}`
1337 .. |KY| replace:: :math:`k_{y}`
1338 .. |KZ| replace:: :math:`k_{z}`
1339 .. |GO| replace:: :math:`g`
1340 .. |RO| replace:: :math:`r`
1341 .. |DPHI| replace:: :math:`\Delta\phi`
1342 .. |DIHR| replace:: :math:`k_{\mathrm{dihr}}`
1343 .. |THET| replace:: :math:`\theta`
1344 .. |NM| replace:: nm\ :math:`^{-1}`
1345 .. |KC| replace:: :math:`k_c`
1346 .. |THEK| replace:: :math:`\theta,k`
1347 .. |R1E| replace:: :math:`r_{1e}`
1348 .. |R2E| replace:: :math:`r_{2e}`
1349 .. |R3E| replace:: :math:`r_{3e}`
1350 .. |KRR| replace:: :math:`k_{rr'}`
1351 .. |KRTH| replace:: :math:`k_{r\theta}`
1352 .. |ALPH| replace:: :math:`\alpha`; |CO| (U nm\ :math:`^{\alpha}`
1353 .. |UM1| replace:: U\ :math:`^{-1}`
1357 .. table:: Details of ``[ moleculetype ]`` directives
1359 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1360 | Name of interaction | Topology file directive | num. | func. | Order of parameters and their units | use in |
1361 | | | atoms [1]_ | type [2]_ | | F.E.? [3]_ |
1362 +====================================+============================+============+===========+=========================================================================+============+
1363 | bond | ``bonds`` [4]_, [5]_ | 2 | 1 | |BZERO| (nm); |KB| |NM2| | all |
1364 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1365 | G96 bond | ``bonds`` [4]_, [5]_ | 2 | 2 | |BZERO| (nm); |KB| |NM4| | all |
1366 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1367 | Morse | ``bonds`` [4]_, [5]_ | 2 | 3 | |BZERO| (nm); |DKJ|; |BETA| | all |
1368 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1369 | cubic bond | ``bonds`` [4]_, [5]_ | 2 | 4 | |BZERO| (nm); |C23| | |
1370 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1371 | connection | ``bonds`` [4]_ | 2 | 5 | | |
1372 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1373 | harmonic potential | ``bonds`` | 2 | 6 | |BZERO| (nm); |KB| |NM2| | all |
1374 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1375 | FENE bond | ``bonds`` [4]_ | 2 | 7 | |BMM| (nm); |KB| |NM2| | |
1376 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1377 | tabulated bond | ``bonds`` [4]_ | 2 | 8 | table number (|GE0|); |KO| |KJM| | |KO| |
1378 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1379 | tabulated bond [6]_ | ``bonds`` | 2 | 9 | table number (|GE0|); |KO| |KJM| | |KO| |
1380 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1381 | restraint potential | ``bonds`` | 2 | 10 | |LUU| (nm); |KDR| (|NM2|) | all |
1382 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1383 | extra LJ or Coulomb | ``pairs`` | 2 | 1 | |MV| [7]_; |MW| [7]_ | all |
1384 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1385 | extra LJ or Coulomb | ``pairs`` | 2 | 2 | fudge QQ (); |QIJ| (e), |MV| [7]_; |MW| [7]_ | |
1386 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1387 | extra LJ or Coulomb | ``pairs_nb`` | 2 | 1 | |QIJ| (e); |MV| [7]_; |MW| [7]_ | |
1388 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1389 | angle | ``angles`` [5]_ | 3 | 1 | |THE0| (deg); |KTHE| (|KJR2|) | all |
1390 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1391 | G96 angle | ``angles`` [5]_ | 3 | 2 | |THE0| (deg); |KTHE| (|KJM|) | all |
1392 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1393 | cross bond-bond | ``angles`` | 3 | 3 | |R1E|, |R2E| (nm); |KRR| (|NM2|) | |
1394 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1395 | cross bond-angle | ``angles`` | 3 | 4 | |R1E|, |R2E|, |R3E| (nm); |KRTH| (|NM2|) | |
1396 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1397 | Urey-Bradley | ``angles`` [5]_ | 3 | 5 | |THE0| (deg); |KTHE| (|KJR2|); |RN13| (nm); |KUB| (|NM2|) | all |
1398 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1399 | quartic angle | ``angles`` [5]_ | 3 | 6 | |THE0| (deg); |C024| (|KJRI|) | |
1400 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1401 | tabulated angle | ``angles`` | 3 | 8 | table number (|GE0|); |KO| (|KJM|) | |KO| |
1402 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1403 | | restricted | | | | | |
1404 | | bending potential | ``angles`` | 3 | 10 | |THE0| (deg); |KTHE| (|KJM|) | |
1405 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1406 | proper dihedral | ``dihedrals`` | 4 | 1 | |PHIS| (deg); |KPHI| (|KJM|); multiplicity | |PHIK| |
1407 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1408 | improper dihedral | ``dihedrals`` | 4 | 2 | |XI0| (deg); |KXI| (|KJR2|) | all |
1409 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1410 | Ryckaert-Bellemans dihedral | ``dihedrals`` | 4 | 3 | |C0|, |C1|, |C2|, |C3|, |C4|, |C5| (|KJM|) | all |
1411 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1412 | periodic improper dihedral | ``dihedrals`` | 4 | 4 | |PHIS| (deg); |KPHI| (|KJM|); multiplicity | |PHIK| |
1413 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1414 | Fourier dihedral | ``dihedrals`` | 4 | 5 | |C1|, |C2|, |C3|, |C4|, |C5| (|KJM|) | all |
1415 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1416 | tabulated dihedral | ``dihedrals`` | 4 | 8 | table number (|GE0|); |KO| (|KJM|) | |KO| |
1417 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1418 | proper dihedral (multiple) | ``dihedrals`` | 4 | 9 | |PHIS| (deg); |KPHI| (|KJM|); multiplicity | |PHIK| |
1419 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1420 | restricted dihedral | ``dihedrals`` | 4 | 10 | |PHI0| (deg); |KPHI| (|KJM|) | |
1421 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1422 | combined bending-torsion potential | ``dihedrals`` | 4 | 11 | |A0|, |A1|, |A2|, |A3|, |A4| (|KJM|) | |
1423 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1424 | exclusions | ``exclusions`` | 1 | | one or more atom indices | |
1425 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1426 | constraint | ``constraints`` [4]_ | 2 | 1 | |BZERO| (nm) | all |
1427 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1428 | constraint [6]_ | ``constraints`` | 2 | 2 | |BZERO| (nm) | all |
1429 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1430 | SETTLE | ``settles`` | 1 | 1 | |DOH|, |DHH| (nm) | |
1431 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1432 | 2-body virtual site | ``virtual_sites2`` | 3 | 1 | |AO| () | |
1433 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1434 | 3-body virtual site | ``virtual_sites3`` | 4 | 1 | |AO|, |BO| () | |
1435 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1436 | 3-body virtual site (fd) | ``virtual_sites3`` | 4 | 2 | |AO| (); |DO| (nm) | |
1437 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1438 | 3-body virtual site (fad) | ``virtual_sites3`` | 4 | 3 | |THET| (deg); |DO| (nm) | |
1439 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1440 | 3-body virtual site (out) | ``virtual_sites3`` | 4 | 4 | |AO|, |BO| (); |CO| (|NM|) | |
1441 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1442 | 4-body virtual site (fdn) | ``virtual_sites4`` | 5 | 2 | |AO|, |BO| (); |CO| (nm) | |
1443 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1444 | N-body virtual site (COG) | ``virtual_sitesn`` | 1 | 1 | one or more constructing atom indices | |
1445 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1446 | N-body virtual site (COM) | ``virtual_sitesn`` | 1 | 2 | one or more constructing atom indices | |
1447 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1448 | N-body virtual site (COW) | ``virtual_sitesn`` | 1 | 3 | | one or more pairs consisting of | |
1449 | | | | | | constructing atom index and weight | |
1450 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1451 | position restraint | ``position_restraints`` | 1 | 1 | |KX|, |KY|, |KZ| (|NM2|) | all |
1452 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1453 | flat-bottomed position restraint | ``position_restraints`` | 1 | 2 | |GO|, |RO| (nm), |KO| (|NM2|) | |
1454 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1455 | distance restraint | ``distance_restraints`` | 2 | 1 | type; label; |LUU| (nm); weight () | |
1456 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1457 | dihedral restraint | ``dihedral_restraints`` | 4 | 1 | |PHI0| (deg); |DPHI| (deg); |DIHR| (|KJR2|) | all |
1458 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1459 | orientation restraint | ``orientation_restraints`` | 2 | 1 | exp.; label; |ALPH|; obs. (U); weight (|UM1|) | |
1460 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1461 | angle restraint | ``angle_restraints`` | 4 | 1 | |THE0| (deg); |KC| (|KJM|); multiplicity | |THEK| |
1462 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1463 | angle restraint (z) | ``angle_restraints_z`` | 2 | 1 | |THE0| (deg); |KC| (|KJM|); multiplicity | |THEK| |
1464 +------------------------------------+----------------------------+------------+-----------+-------------------------------------------------------------------------+------------+
1467 The required number of atom indices for this directive
1470 The index to use to select this function type
1473 Indicates which of the parameters can be interpolated in free energy calculations
1476 This interaction type will be used by :ref:`grompp <gmx grompp>` for generating exclusions
1479 This interaction type can be converted to constraints by :ref:`grompp <gmx grompp>`
1482 The combination rule determines the type of LJ parameters, see
1485 No connection, and so no exclusions, are generated for this interaction
1487 Description of the file layout:
1489 - Semicolon (;) and newline characters surround comments
1491 - On a line ending with :math:`\backslash` the newline character is
1494 - Directives are surrounded by ``[`` and ``]``
1496 - The topology hierarchy (which must be followed) consists of three
1499 - the parameter level, which defines certain force-field
1500 specifications (see :numref:`Table %s <tab-topfile1>`)
1502 - the molecule level, which should contain one or more molecule
1503 definitions (see :numref:`Table %s <tab-topfile2>`)
1505 - the system level, containing only system-specific information
1506 (``[ system ]`` and ``[ molecules ]``)
1508 - Items should be separated by spaces or tabs, not commas
1510 - Atoms in molecules should be numbered consecutively starting at 1
1512 - Atoms in the same charge group must be listed consecutively
1514 - The file is parsed only once, which implies that no forward
1515 references can be treated: items must be defined before they can be
1518 - Exclusions can be generated from the bonds or overridden manually
1520 - The bonded force types can be generated from the atom types or
1523 - It is possible to apply multiple bonded interactions of the same type
1526 - Descriptive comment lines and empty lines are highly recommended
1528 - Starting with |Gromacs| version 3.1.3, all directives at the parameter
1529 level can be used multiple times and there are no restrictions on the
1530 order, except that an atom type needs to be defined before it can be
1531 used in other parameter definitions
1533 - If parameters for a certain interaction are defined multiple times
1534 for the same combination of atom types the last definition is used;
1535 starting with |Gromacs| version 3.1.3 :ref:`grompp <gmx grompp>` generates
1536 a warning for parameter redefinitions with different values
1538 - Using one of the ``[ atoms ]``,
1539 ``[ bonds ]``, ``[ pairs ]``,
1540 ``[ angles ]``, etc. without having used
1541 ``[ moleculetype ]`` before is meaningless and generates
1544 - Using ``[ molecules ]`` without having used
1545 ``[ system ]`` before is meaningless and generates a
1548 - After ``[ system ]`` the only allowed directive is
1551 - Using an unknown string in ``[ ]`` causes all the data
1552 until the next directive to be ignored and generates a warning
1554 Here is an example of a topology file, ``urea.top``:
1559 ; Example topology file
1561 ; The force-field files to be included
1562 #include "amber99.ff/forcefield.itp"
1569 1 C 1 URE C 1 0.880229 12.01000 ; amber C type
1570 2 O 1 URE O 2 -0.613359 16.00000 ; amber O type
1571 3 N 1 URE N1 3 -0.923545 14.01000 ; amber N type
1572 4 H 1 URE H11 4 0.395055 1.00800 ; amber H type
1573 5 H 1 URE H12 5 0.395055 1.00800 ; amber H type
1574 6 N 1 URE N2 6 -0.923545 14.01000 ; amber N type
1575 7 H 1 URE H21 7 0.395055 1.00800 ; amber H type
1576 8 H 1 URE H22 8 0.395055 1.00800 ; amber H type
1588 ; ai aj ak al funct definition
1603 [ position_restraints ]
1604 ; you wouldn't normally use this for a molecule like Urea,
1605 ; but we include it here for didactic purposes
1607 1 1 1000 1000 1000 ; Restrain to a point
1608 2 1 1000 0 1000 ; Restrain to a line (Y-axis)
1609 3 1 1000 0 0 ; Restrain to a plane (Y-Z-plane)
1611 [ dihedral_restraints ]
1612 ; ai aj ak al type phi dphi fc
1616 ; Include TIP3P water topology
1617 #include "amber99/tip3p.itp"
1627 Here follows the explanatory text.
1629 **#include “amber99.ff/forcefield.itp” :** this includes
1630 the information for the force field you are using, including bonded and
1631 non-bonded parameters. This example uses the AMBER99 force field, but
1632 your simulation may use a different force field. :ref:`grompp <gmx grompp>`
1633 will automatically go and find this file and copy-and-paste its content.
1634 That content can be seen in
1635 ``share/top/amber99.ff/forcefield.itp}``, and it
1644 ; nbfunc comb-rule gen-pairs fudgeLJ fudgeQQ
1647 #include "ffnonbonded.itp"
1648 #include "ffbonded.itp"
1650 The two ``#define`` statements set up the conditions so that
1651 future parts of the topology can know that the AMBER 99 force field is
1656 - ``nbfunc`` is the non-bonded function type. Use 1 (Lennard-Jones) or 2
1659 - ``comb-rule`` is the number of the combination rule (see :ref:`nbpar`).
1661 - ``gen-pairs`` is for pair generation. The default is
1662 ‘no’, *i.e.* get 1-4 parameters from the pairtypes list. When
1663 parameters are not present in the list, stop with a fatal error.
1664 Setting ‘yes’ generates 1-4 parameters that are not present in the
1665 pair list from normal Lennard-Jones parameters using
1668 - ``fudgeLJ`` is the factor by which to multiply
1669 Lennard-Jones 1-4 interactions, default 1
1671 - ``fudgeQQ`` is the factor by which to multiply
1672 electrostatic 1-4 interactions, default 1
1674 - :math:`N` is the power for the repulsion term in a 6-\ :math:`N`
1675 potential (with nonbonded-type Lennard-Jones only), starting with
1676 |Gromacs| version 4.5, :ref:`grompp <gmx mdrun>` also reads and applies
1677 :math:`N`, for values not equal to 12 tabulated interaction functions
1678 are used (in older version you would have to use user tabulated
1681 **Note** that ``gen-pairs``, ``fudgeLJ``,
1682 ``fudgeQQ``, and :math:`N` are optional.
1683 ``fudgeLJ`` is only used when generate pairs is set to
1684 ‘yes’, and ``fudgeQQ`` is always used. However, if you want
1685 to specify :math:`N` you need to give a value for the other parameters
1688 Then some other ``#include`` statements add in the large
1689 amount of data needed to describe the rest of the force field. We will
1690 skip these and return to ``urea.top``. There we will see
1692 **[ moleculetype ] :** defines the name of your molecule
1693 in this :ref:`top` and nrexcl = 3 stands for excluding
1694 non-bonded interactions between atoms that are no further than 3 bonds
1697 **[ atoms ] :** defines the molecule, where
1698 ``nr`` and ``type`` are fixed, the rest is user
1699 defined. So ``atom`` can be named as you like,
1700 ``cgnr`` made larger or smaller (if possible, the total
1701 charge of a charge group should be zero), and charges can be changed
1704 **[ bonds ] :** no comment.
1706 **[ pairs ] :** LJ and Coulomb 1-4 interactions
1708 **[ angles ] :** no comment
1710 **[ dihedrals ] :** in this case there are 9 proper
1711 dihedrals (funct = 1), 3 improper (funct = 4) and no Ryckaert-Bellemans
1712 type dihedrals. If you want to include Ryckaert-Bellemans type dihedrals
1713 in a topology, do the following (in case of *e.g.* decane):
1718 ; ai aj ak al funct c0 c1 c2
1722 In the original implementation of the potential for
1723 alkanes \ :ref:`131 <refRyckaert78>` no 1-4 interactions were used, which means that in
1724 order to implement that particular force field you need to remove the
1725 1-4 interactions from the ``[ pairs ]`` section of your
1726 topology. In most modern force fields, like OPLS/AA or Amber the rules
1727 are different, and the Ryckaert-Bellemans potential is used as a cosine
1728 series in combination with 1-4 interactions.
1730 **[ position_restraints ] :** harmonically restrain the selected particles to reference
1731 positions (:ref:`positionrestraint`). The reference positions are read
1732 from a separate coordinate file by :ref:`grompp <gmx grompp>`.
1734 **[ dihedral_restraints ] :** restrain selected dihedrals to a reference value. The
1735 implementation of dihedral restraints is described in section
1736 :ref:`dihedralrestraint` of the manual. The parameters specified in
1737 the ``[dihedral_restraints]`` directive are as follows:
1739 - ``type`` has only one possible value which is 1
1741 - ``phi`` is the value of :math:`\phi_0` in :eq:`eqn. %s <eqndphi>` and
1742 :eq:`eqn. %s <eqndihre>` of the manual.
1744 - ``dphi`` is the value of :math:`\Delta\phi` in :eq:`eqn. %s <eqndihre>` of the
1747 - ``fc`` is the force constant :math:`k_{dihr}` in :eq:`eqn. %s <eqndihre>` of the
1750 **#include “tip3p.itp” :** includes a topology file that was already
1751 constructed (see section :ref:`molitp`).
1753 **[ system ] :** title of your system, user-defined
1755 **[ molecules ] :** this defines the total number of (sub)molecules in your system
1756 that are defined in this :ref:`top`. In this example file, it stands for 1
1757 urea molecule dissolved in 1000 water molecules. The molecule type ``SOL``
1758 is defined in the ``tip3p.itp`` file. Each name here must correspond to a
1759 name given with ``[ moleculetype ]`` earlier in the topology. The order of the blocks of
1760 molecule types and the numbers of such molecules must match the
1761 coordinate file that accompanies the topology when supplied to :ref:`grompp <gmx grompp>`.
1762 The blocks of molecules do not need to be contiguous, but some tools
1763 (e.g. :ref:`genion <gmx genion>`) may act only on the first or last such block of a
1764 particular molecule type. Also, these blocks have nothing to do with the
1765 definition of groups (see sec. :ref:`groupconcept` and
1766 sec. :ref:`usinggroups`).
1773 If you construct a topology file you will use frequently (like the water
1774 molecule, ``tip3p.itp``, which is already constructed for
1775 you) it is good to make a ``molecule.itp`` file. This only
1776 lists the information of one particular molecule and allows you to
1777 re-use the ``[ moleculetype ]`` in multiple systems without
1778 re-invoking :ref:`pdb2gmx <gmx pdb2gmx>` or manually copying and pasting. An
1779 example ``urea.itp`` follows:
1788 1 C 1 URE C 1 0.880229 12.01000 ; amber C type
1790 8 H 1 URE H22 8 0.395055 1.00800 ; amber H type
1797 ; ai aj ak al funct definition
1806 Using :ref:`itp` files results in a very short
1812 ; Example topology file
1814 ; The force field files to be included
1815 #include "amber99.ff/forcefield.itp"
1819 ; Include TIP3P water topology
1820 #include "amber99/tip3p.itp"
1833 A very powerful feature in |Gromacs| is the use of ``#ifdef``
1834 statements in your :ref:`top` file. By making use of this
1835 statement, and associated ``#define`` statements like were
1836 seen in ``amber99.ff/forcefield.itp`` earlier,
1837 different parameters for one molecule can be used in the same
1838 :ref:`top` file. An example is given for TFE, where there is
1839 an option to use different charges on the atoms: charges derived by De
1840 Loof et al. :ref:`132 <refLoof92>` or by Van Buuren and
1841 Berendsen \ :ref:`133 <refBuuren93a>`. In fact, you can use much of the
1842 functionality of the C preprocessor, ``cpp``, because
1843 :ref:`grompp <gmx grompp>` contains similar pre-processing functions to scan
1844 the file. The way to make use of the ``#ifdef`` option is as
1847 - either use the option ``define = -DDeLoof`` in the
1848 :ref:`mdp` file (containing :ref:`grompp <gmx grompp>` input
1849 parameters), or use the line ``#define DeLoof`` early in
1850 your :ref:`top` or :ref:`itp` file; and
1852 - put the ``#ifdef`` statements in your
1853 :ref:`top`, as shown below:
1863 ; nr type resnr residu atom cgnr charge mass
1865 ; Use Charges from DeLoof
1870 5 CH2 1 TFE CH2 1 0.25
1871 6 OA 1 TFE OA 1 -0.65
1872 7 HO 1 TFE HO 1 0.41
1874 ; Use Charges from VanBuuren
1879 5 CH2 1 TFE CH2 1 0.26
1880 6 OA 1 TFE OA 1 -0.55
1886 6 7 1 1.000000e-01 3.138000e+05
1887 1 2 1 1.360000e-01 4.184000e+05
1888 1 3 1 1.360000e-01 4.184000e+05
1889 1 4 1 1.360000e-01 4.184000e+05
1890 1 5 1 1.530000e-01 3.347000e+05
1891 5 6 1 1.430000e-01 3.347000e+05
1894 This mechanism is used by :ref:`pdb2gmx <gmx pdb2gmx>` to implement optional position
1895 restraints (:ref:`positionrestraint`) by ``#include``-ing an :ref:`itp` file
1896 whose contents will be meaningful only if a particular ``#define`` is set
1897 (and spelled correctly!)
1899 Topologies for free energy calculations
1900 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
1902 Free energy differences between two systems, A and B, can be calculated
1903 as described in sec. :ref:`fecalc`. Systems A and B are described by
1904 topologies consisting of the same number of molecules with the same
1905 number of atoms. Masses and non-bonded interactions can be perturbed by
1906 adding B parameters under the ``[ atoms ]`` directive. Bonded interactions can be
1907 perturbed by adding B parameters to the bonded types or the bonded
1908 interactions. The parameters that can be perturbed are listed in
1909 :numref:`Tables %s <tab-topfile1>` and :numref:`%s <tab-topfile2>`.
1910 The :math:`\lambda`-dependence of the
1911 interactions is described in section sec. :ref:`feia`. The bonded
1912 parameters that are used (on the line of the bonded interaction
1913 definition, or the ones looked up on atom types in the bonded type
1914 lists) is explained in :numref:`Table %s <tab-topfe>`. In most cases, things should
1915 work intuitively. When the A and B atom types in a bonded interaction
1916 are not all identical and parameters are not present for the B-state,
1917 either on the line or in the bonded types, :ref:`grompp <gmx grompp>` uses the A-state
1918 parameters and issues a warning. For free energy calculations, all or no
1919 parameters for topology B (:math:`\lambda = 1`) should be added on the
1920 same line, after the normal parameters, in the same order as the normal
1921 parameters. From |Gromacs| 4.6 onward, if :math:`\lambda` is treated as a
1922 vector, then the ``bonded-lambdas`` component controls all bonded terms that
1923 are not explicitly labeled as restraints. Restrain terms are controlled
1924 by the ``restraint-lambdas`` component.
1926 .. |NOT| replace:: :math:`-`
1930 .. table:: The bonded parameters that are used for free energy topologies,
1931 on the line of the bonded interaction definition or looked up
1932 in the bond types section based on atom types. A and B indicate the
1933 parameters used for state A and B respectively, + and |NOT| indicate
1934 the (non-)presence of parameters in the topology, x indicates that
1935 the presence has no influence.
1937 +--------------------+---------------+---------------------------------+---------+
1938 | B-state atom types | parameters | parameters in bonded types | |
1939 + + +-----------------+---------------+ +
1940 | all identical to | on line | A atom types | B atom types | message |
1941 + +-------+-------+-------+---------+-------+-------+ +
1942 | A-state atom types | A | B | A | B | A | B | |
1943 +====================+=======+=======+=======+=========+=======+=======+=========+
1944 | | +AB | |NOT| | x | x | | | |
1945 | | +A | +B | x | x | | | |
1946 | yes | |NOT| | |NOT| | |NOT| | |NOT| | | | error |
1947 | | |NOT| | |NOT| | +AB | |NOT| | | | |
1948 | | |NOT| | |NOT| | +A | +B | | | |
1949 +--------------------+-------+-------+-------+---------+-------+-------+---------+
1950 | | +AB | |NOT| | x | x | x | x | warning |
1951 | | +A | +B | x | x | x | x | |
1952 | | |NOT| | |NOT| | |NOT| | |NOT| | x | x | error |
1953 | no | |NOT| | |NOT| | +AB | |NOT| | |NOT| | |NOT| | warning |
1954 | | |NOT| | |NOT| | +A | +B | |NOT| | |NOT| | warning |
1955 | | |NOT| | |NOT| | +A | x | +B | |NOT| | |
1956 | | |NOT| | |NOT| | +A | x | + | +B | |
1957 +--------------------+-------+-------+-------+---------+-------+-------+---------+
1961 Below is an example of a topology which changes from 200 propanols to
1962 200 pentanes using the GROMOS-96 force field.
1967 ; Include force field parameters
1968 #include "gromos43a1.ff/forcefield.itp"
1975 ; nr type resnr residue atom cgnr charge mass typeB chargeB massB
1976 1 H 1 PROP PH 1 0.398 1.008 CH3 0.0 15.035
1977 2 OA 1 PROP PO 1 -0.548 15.9994 CH2 0.0 14.027
1978 3 CH2 1 PROP PC1 1 0.150 14.027 CH2 0.0 14.027
1979 4 CH2 1 PROP PC2 2 0.000 14.027
1980 5 CH3 1 PROP PC3 2 0.000 15.035
1983 ; ai aj funct par_A par_B
1995 ; ai aj ak funct par_A par_B
2001 ; ai aj ak al funct par_A par_B
2002 1 2 3 4 1 gd_12 gd_17
2003 2 3 4 5 1 gd_17 gd_17
2013 Atoms that are not perturbed, ``PC2`` and
2014 ``PC3``, do not need B-state parameter specifications, since
2015 the B parameters will be copied from the A parameters. Bonded
2016 interactions between atoms that are not perturbed do not need B
2017 parameter specifications, as is the case for the last bond in the
2018 example topology. Topologies using the OPLS/AA force field need no
2019 bonded parameters at all, since both the A and B parameters are
2020 determined by the atom types. Non-bonded interactions involving one or
2021 two perturbed atoms use the free-energy perturbation functional forms.
2022 Non-bonded interactions between two non-perturbed atoms use the normal
2023 functional forms. This means that when, for instance, only the charge of
2024 a particle is perturbed, its Lennard-Jones interactions will also be
2025 affected when lambda is not equal to zero or one.
2027 **Note** that this topology uses the GROMOS-96 force field, in which the
2028 bonded interactions are not determined by the atom types. The bonded
2029 interaction strings are converted by the C-preprocessor. The force-field
2030 parameter files contain lines like:
2034 #define gb_26 0.1530 7.1500e+06
2036 #define gd_17 0.000 5.86 3
2038 .. _constraintforce:
2043 | The constraint force between two atoms in one molecule can be
2044 calculated with the free energy perturbation code by adding a
2045 constraint between the two atoms, with a different length in the A and
2046 B topology. When the B length is 1 nm longer than the A length and
2047 lambda is kept constant at zero, the derivative of the Hamiltonian
2048 with respect to lambda is the constraint force. For constraints
2049 between molecules, the pull code can be used, see sec. :ref:`pull`.
2050 Below is an example for calculating the constraint force at 0.7 nm
2051 between two methanes in water, by combining the two methanes into one
2052 “molecule.” **Note** that the definition of a “molecule” in |Gromacs|
2053 does not necessarily correspond to the chemical definition of a
2054 molecule. In |Gromacs|, a “molecule” can be defined as any group of
2055 atoms that one wishes to consider simultaneously. The added constraint
2056 is of function type 2, which means that it is not used for generating
2057 exclusions (see sec. :ref:`excl`). Note that the constraint free energy
2058 term is included in the derivative term, and is specifically included
2059 in the ``bonded-lambdas`` component. However, the free energy for changing
2060 constraints is *not* included in the potential energy differences used
2061 for BAR and MBAR, as this requires reevaluating the energy at each of
2062 the constraint components. This functionality is planned for later
2067 ; Include force-field parameters
2068 #include "gromos43a1.ff/forcefield.itp"
2075 ; nr type resnr residu atom cgnr charge mass
2076 1 CH4 1 CH4 C1 1 0 16.043
2077 2 CH4 1 CH4 C2 2 0 16.043
2079 ; ai aj funct length_A length_B
2082 #include "gromos43a1.ff/spc.itp"
2096 Files with the :ref:`gro` file extension contain a molecular
2097 structure in GROMOS-87 format. A sample piece is included below:
2101 MD of 2 waters, reformat step, PA aug-91
2103 1WATER OW1 1 0.126 1.624 1.679 0.1227 -0.0580 0.0434
2104 1WATER HW2 2 0.190 1.661 1.747 0.8085 0.3191 -0.7791
2105 1WATER HW3 3 0.177 1.568 1.613 -0.9045 -2.6469 1.3180
2106 2WATER OW1 4 1.275 0.053 0.622 0.2519 0.3140 -0.1734
2107 2WATER HW2 5 1.337 0.002 0.680 -1.0641 -1.1349 0.0257
2108 2WATER HW3 6 1.326 0.120 0.568 1.9427 -0.8216 -0.0244
2109 1.82060 1.82060 1.82060
2111 This format is fixed, *i.e.* all columns are in a fixed position. If you
2112 want to read such a file in your own program without using the |Gromacs|
2113 libraries you can use the following formats:
2116 ``“%5i%5s%5s%5i%8.3f%8.3f%8.3f%8.4f%8.4f%8.4f”``
2118 Or to be more precise, with title *etc.* it looks like this:
2124 for (i=0; (i<natoms); i++) {
2125 "%5d%-5s%5s%5d%8.3f%8.3f%8.3f%8.4f%8.4f%8.4f\n",
2126 residuenr,residuename,atomname,atomnr,x,y,z,vx,vy,vz
2128 "%10.5f%10.5f%10.5f%10.5f%10.5f%10.5f%10.5f%10.5f%10.5f\n",
2129 box[X][X],box[Y][Y],box[Z][Z],
2130 box[X][Y],box[X][Z],box[Y][X],box[Y][Z],box[Z][X],box[Z][Y]
2133 ``(i5,2a5,i5,3f8.3,3f8.4)``
2135 So ``confin.gro`` is the |Gromacs| coordinate file and is
2136 almost the same as the GROMOS-87 file (for GROMOS users: when used with
2137 ``ntx=7``). The only difference is the box for which |Gromacs|
2138 uses a tensor, not a vector.
2142 Force field organization
2143 ------------------------
2150 Many force fields are available by default. Force fields are detected by
2151 the presence of ``<name>.ff`` directories in the
2152 ``$GMXLIB/share/gromacs/top`` sub-directory and/or the
2153 working directory. The information regarding the location of the force
2154 field files is printed by :ref:`pdb2gmx <gmx pdb2gmx>` so you can easily keep
2155 track of which version of a force field is being called, in case you
2156 have made modifications in one location or another. The force fields
2157 included with |Gromacs| are:
2159 - AMBER03 protein, nucleic AMBER94 (Duan et al., J. Comp. Chem. 24,
2162 - AMBER94 force field (Cornell et al., JACS 117, 5179-5197, 1995)
2164 - AMBER96 protein, nucleic AMBER94 (Kollman et al., Acc. Chem. Res. 29,
2167 - AMBER99 protein, nucleic AMBER94 (Wang et al., J. Comp. Chem. 21,
2170 - AMBER99SB protein, nucleic AMBER94 (Hornak et al., Proteins 65,
2173 - AMBER99SB-ILDN protein, nucleic AMBER94 (Lindorff-Larsen et al.,
2174 Proteins 78, 1950-58, 2010)
2176 - AMBERGS force field (Garcia & Sanbonmatsu, PNAS 99, 2782-2787, 2002)
2178 - CHARMM27 all-atom force field (CHARM22 plus CMAP for proteins)
2180 - GROMOS96 43a1 force field
2182 - GROMOS96 43a2 force field (improved alkane dihedrals)
2184 - GROMOS96 45a3 force field (Schuler JCC 2001 22 1205)
2186 - GROMOS96 53a5 force field (JCC 2004 vol 25 pag 1656)
2188 - GROMOS96 53a6 force field (JCC 2004 vol 25 pag 1656)
2190 - GROMOS96 54a7 force field (Eur. Biophys. J. (2011), 40,, 843-856,
2191 DOI: 10.1007/s00249-011-0700-9)
2193 - OPLS-AA/L all-atom force field (2001 aminoacid dihedrals)
2195 A force field is included at the beginning of a topology file with an
2196 ``#include`` statement followed by
2197 ``<name>.ff/forcefield.itp``. This statement includes the
2198 force-field file, which, in turn, may include other force-field files.
2199 All the force fields are organized in the same way. An example of the
2200 ``amber99.ff/forcefield.itp`` was shown in
2203 For each force field, there several files which are only used by
2204 :ref:`pdb2gmx <gmx pdb2gmx>`. These are: residue databases
2205 (:ref:`rtp`) the hydrogen
2206 database (:ref:`hdb`), two
2207 termini databases (``.n.tdb`` and ``.c.tdb``,
2208 see ) and the atom type database
2210 contains only the masses. Other optional files are described in sec. :ref:`pdb2gmxfiles`.
2212 Changing force-field parameters
2213 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
2215 If one wants to change the parameters of few bonded interactions in a
2216 molecule, this is most easily accomplished by typing the parameters
2217 behind the definition of the bonded interaction directly in the
2218 :ref:`top` file under the ``[ moleculetype ]``
2219 section (see :ref:`topfile` for the format and units).
2220 If one wants to change the parameters for all instances of a
2221 certain interaction one can change them in the force-field file or add a
2222 new ``[ ???types ]`` section after including the force
2223 field. When parameters for a certain interaction are defined multiple
2224 times, the last definition is used. As of |Gromacs| version 3.1.3, a
2225 warning is generated when parameters are redefined with a different
2226 value. Changing the Lennard-Jones parameters of an atom type is not
2227 recommended, because in the GROMOS force fields the Lennard-Jones
2228 parameters for several combinations of atom types are not generated
2229 according to the standard combination rules. Such combinations (and
2230 possibly others that do follow the combination rules) are defined in the
2231 ``[ nonbond_params ]`` section, and changing the
2232 Lennard-Jones parameters of an atom type has no effect on these
2238 As of |Gromacs| version 3.1.3, atom types can be added in an extra
2239 ``[ atomtypes ]`` section after the inclusion of the
2240 normal force field. After the definition of the new atom type(s),
2241 additional non-bonded and pair parameters can be defined. In pre-3.1.3
2242 versions of |Gromacs|, the new atom types needed to be added in the
2243 ``[ atomtypes ]`` section of the force-field files, because
2244 all non-bonded parameters above the last ``[ atomtypes ]``
2245 section would be overwritten using the standard combination rules.