Suggested run parameters for a crystal simulation

As has been stated already, crystal simulations are in principle no different than other condensed-phase, periodic molecular simulations. However, in solution-phase simulations, the protein is surrounded by a large box of water and the system can be treated more or less isotropically. In a crystal simulation, the system is almost certainly not isotropic in terms of its stiffness, electrostatics, and other critical features. This distinction prompts the use of certain parameters in crystal simulations that may be slightly more expensive than typical solution-phase simulations. A good comparison can be made to simulations of membranes, which also require consideration of the anisotropic nature of the system.

The simplest change to make is to ensure that the long-ranged van-der Waals forces are treated accurately. Protein:water interfaces are found throughout the crystal lattice, so long-range corrections to the virial (and energy) based on the assumption that the system composition around each particle is homogeneous beyond the Lennard-Jones cutoff are not as reliable as they are in a solution-phase simulation. We suggest maintaining this approximation (which is engaged by default in sander and pmemd) at some level, but pushing back the cutoff for explicit calculation of Lennard-Jones interactions to at least 10 or as much as 12 Angstroms. The pmemd program offers the option of setting different cutoffs for Lennard-Jones and electrostatic terms, which can mitigate the cost of the longer cutoff. Setting es_cutoff=9.0 and vdw_cutoff=11.0 will deliver considerably more precise Lennard-Jones forces than simply specifying cut=9.0, but only increase the cost of the simulation by 5-10%, a cost which can often be recoverd by optimizing electrostatic parameters in the &ewald namelist (though such optimizations should only be attempted, as stated in the manual, by people who really know what they are doing).

The next critical change to make is to set ntp=2, which will specify constant pressure dynamics with anisotropic volume rescaling depending on all three elements of the virial trace. However, in the current implementations of sander and pmemd this is only possible with orthorhombic simulation cells (when all unit cell angles are 90 degrees). If the unit cell is mono- or tri-clinic, ntp=1 should be used. Performing constant-pressure dynamics requires a thermostat, and for conditions in which the simulation cell dimensions vary at each step the Langevin thermostat (ntt=3) is the best option, although users must take care to specify new random seeds for the thermostat each time the simulation is checkpointed and restarted.

View a list of references related to simulations of crystal structures.

Return to the beginning of the tutorial.