The vast majority of the CVs implemented in PLUMED are calculated from a list of atom positions. Within PLUMED atoms are specified using their numerical indices in the molecular dynamics input file.
In PLUMED lists of atoms can be either provided directly inside the definition of each collective variable, or predefined as a GROUP that can be reused multiple times. Lists of atoms can be written as:
Some collective variable must accept a fixed number of atoms, for example a DISTANCE is calculated using two atoms only, an ANGLE is calcuated using either 3 or 4 atoms and TORSION is calculated using 4 atoms.
Additional material and examples can be also found in the tutorial Belfast tutorial: Analyzing CVs.
In addition, for certain colvars, pdb files can be read in using the following keywords and used to select ATOMS:
MOLINFO | This command is used to provide information on the molecules that are present in your system. |
The information on the molecules in your system can either be provided in the form of a pdb file or as a set of lists of atoms that describe the various chains in your system using MOLINFO. If a pdb file is used plumed the MOLINFO command will endeavor to recognize the various chains and residues that make up the molecules in your system using the chainIDs and resnumbers from the pdb file. You can then use this information in commands where this has been implemented to specify atom lists. One place where this is particularly useful is when using the commands ALPHARMSD, ANTIBETARMSD and PARABETARMSD.
MOLINFO also introduces special groups that can be used in atom selection. These special groups always begin with a @ symbol. The following special groups are currently available in PLUMED:
Symbol | Topology type | Despription |
@phi-# | protein | The torsional angle defined by the C, CA, N and C atoms of the protein backbone in the #th residue. See http://en.wikipedia.org/wiki/Ramachandran_plot |
@psi-# | protein | The torsional angle defined by the N, C, CA and N atoms of the protein backbone in the #th residue. See http://en.wikipedia.org/wiki/Ramachandran_plot |
@omega-# | protein | The torsional angle defined by the CA, N, C and CA atoms of the protein backbone in the #th residue. See http://en.wikipedia.org/wiki/Ramachandran_plot |
@chi1-# | protein | The first torsional angle of the sidechain of the #th residue. Be aware that this angle is not defined for GLY or ALA residues. See http://en.wikipedia.org/wiki/Ramachandran_plot |
The following example shows how to use MOLINFO with TORSION to calculate the torsion angles phi and psi for the first and fourth residue of the protein:
MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb t1: TORSION ATOMS=@phi-3 t2: TORSION ATOMS=@psi-4 PRINT ARG=t1,t2 FILE=colvar STRIDE=10
PLUMED is designed so that for the majority of the CVs implemented the periodic boundary conditions are treated in the same manner as they would be treated in the host code. In some codes this can be problematic when the colvars you are using involve some property of a molecule. These codes allow the atoms in the molecules to become separated by periodic boundaries, a fact which PLUMED could only deal with were the topology passed from the MD code to PLUMED. Making this work would involve a lot laborious programming and goes against our original aim of having a general patch that can be implemented in a wide variety of MD codes. Consequentially, we have implemented a more pragmatic solution to this probem - the user specifies in input any molecules (or parts of molecules) that must be kept in tact throughout the simulation run. In PLUMED 1 this was done using the ALIGN_ATOMS keyword. In PLUMED 2 the same effect can be achieved using the WHOLEMOLECULES command.
The following input computes the end-to-end distance for a polymer of 100 atoms and keeps it at a value around 5.
WHOLEMOLECULES ENTITY0=1-100 e2e: DISTANCE ATOMS=1,100 NOPBC RESTRAINT ARG=e2e KAPPA=1 AT=5
Notice that NOPBC is used to be sure in DISTANCE that if the end-to-end distance is larger than half the simulation box the distance is compute properly. Also notice that, since many MD codes break molecules across cell boundary, it might be necessary to use the WHOLEMOLECULES keyword (also notice that it should be before distance).
Notice that most expressions are invariant with respect to a change in the order of the atoms, but some of them depend on that order. E.g., with WHOLEMOLECULES it could be useful to specify atom lists in a reversed order.
# to see the effect, one could dump the atoms as they were before molecule reconstruction: # DUMPATOMS FILE=dump-broken.xyz ATOMS=1-20 WHOLEMOLECULES STRIDE=1 ENTITY0=1-20 DUMPATOMS FILE=dump.xyz ATOMS=1-20
Notice that since PLUMED 2.1 it is also possible to shift coordinates stored within PLUMED so as to align them to a template structure, using the FIT_TO_TEMPLATE keyword.
Finally, you can bring a set of atom as close as possible to another set of atoms using the WRAPAROUND keyword.
Sometimes, when calculating a colvar, you may not want to use the positions of a number of atoms directly. Instead you may wish to use the position of a virtual atom whose position is generated based on the positions of a collection of other atoms. For example you might want to use the center of mass of a group of atoms. Plumed has a number of routines for calculating the positions of these virtual atoms from lists of atoms:
CENTER | Calculate the center for a group of atoms, with arbitrary weights. |
COM | Calculate the center of mass for a group of atoms. |
GHOST | Calculate the absolute position of a ghost atom with fixed coordinatesin the local reference frame formed by three atoms. The computed ghost atom is stored as a virtual atom that can be accessed inan atom list through the the label for the GHOST action that creates it. |
To specify to a colvar that you want to use the position of a virtual atom to calculate a colvar rather than one of the atoms in your system you simply use the label for your virtual atom in place of the usual numerical index. Virtual atoms and normal atoms can be mixed together in the input to colvars as shown below:
COM ATOMS=1,10 LABEL=com1 DISTANCE ATOMS=11,com1
If you don't want to calculate CVs from the virtual atom. That is to say you just want to monitor the position of a virtual atom (or any set of atoms) over the course of your trajectory you can do this using DUMPATOMS.