Simulations

Introduction

Assorted topics relevant to programming particle-based simulation codes and to using these codes for the modeling of a wide range of systems, notably complex fluids.

Molecular dynamics simulations

For many projects, the CSML uses LAMMPS for molecular dynamics, although (depending on the application) we also employ NAMD and GROMACS.

For educational purposes, Moldy is strongly recommended, not in the least because it offers an excellent manual. For class purposes, we maintain a separate list of executables for various operating systems.

Many questions regarding LAMMPS can be resolved by consulting the manual, but we address some common problems below.

LAMMPS Special Usage Notes

Temperature Normalization

By default LAMMPS normalizes the temperature by an amount n_dof - d, where n_dof is the system's total number of degrees of freedom and d its dimensionality. Subtracting d accounts for the center-of-mass motion of the system. This leads to an incorrect reported value if the system has a proper frame of reference, e.g., when using a Langevin thermostat in which all particles interact with a stationary background solvent. In this case it is necessary to ensure n_dof is used instead of n_dof - d. To do this, use compute_modify as follows

compute myTemp all temp
compute_modify myTemp extra 0
thermo_modify temp myTemp

As a note, the above only affects the reported temperature. The dynamics are computed correctly regardless.

Compute RDF using 'rerun' command

The 'rerun' command in LAMMPS performs a post-processing simulation by reading the atom information line-by-line from the dump file(s) created from a previous simulation. The command syntax is as follows:

 rerun file1 file2 ... keyword args ... 

A detailed description of the syntax can be found on the LAMMPS website.

Besides the fact that the atoms' positions (and possibly velocities, etc.) are pre-determined from the dump file(s), we use the rerun command as if we are running a normal simulation (with some differences and limitations, explained below). When the rerun command is called, it invokes the read_dump command to read in lines from the dumpfile(s) line-by-line, each time invoking the run command to output computed energy, forces, and any thermo output or diagnostic info the user has defined. Thus, in the input file for this pseudo simulation, we must define a system, units, dimensions, box, etc, and these will typically be identical to the original simulation.

Commands from the original simulation that will not be included are ones such as dump commands and time integration fixes (e.g. fix nve; rerun only looks at single moments in time and cannot perform time integration). Fixes that constrain forces on atoms (such as fix langevin) can be invoked in general, but it does not make sense to do this for computing the RDF (even though the langevin thermostat may be employed in the original simulation).

As an example, let us consider computing the RDF for a typical Lennard-Jones fluid past the interaction cutoff, and let us assume that we have already generated a dumpfile containing information on the atom positions over some set of timesteps. Then we will run a second simulation that reads in the particle positions from the dumpfile(s) over some subset of the original recorded timesteps (see arguments for the rerun command), and will compute and output the RDF with the following commands:

compute    rdfID groupID rdf N #computes rdf with N bins
fix        fixID groupID ave/time Nevery Nrepeat Nfreq c_[rdfID] file rdf.dat mode vector # see note below
rerun      dump.dat dump x y z # 'dump.dat' is the dumpfile to be read 

Notes:

• Since the compute rdf command will only compute the RDF up to the interaction cutoff distance, we must change this parameter in the pair_style and pair_coeff commands so that we can obtain the RDF over the desired domain (i.e. if we want to compute the RDF up to a cutoff of 4.0, we would set the 'cutoff' arguments in those commands to 4.0).
• While the rerun command creates a set of atoms at every snapshot of the dumpfile that it reads, the compute rdf command expects a set of atoms to be present at the start of the rerun simulation (remember, the compute command comes before the rerun command) and will produce an error if no atoms are present. To avoid this, one can use the create atoms command (or read in the data file via the read data command) used for creating atoms in the original simulation at the beginning of the rerun.
• When using rerun for 2D simulations, set dimension to 2 in the input file. For the rerun command, only read in two coordinates from dump file, for example "rerun dump.dat dump x y".