OpenMM

What it simulates

OpenMM runs classical molecular dynamics from Python. You give it a structure, a force field and an integrator. It minimises the energy, then integrates Newton's equations with a Langevin thermostat. Atoms interact through the AMBER or CHARMM force fields, and PME handles long-range electrostatics in periodic boxes.

The recipe wraps this in one helper script. Input is a PDB file with hydrogens and, for a periodic run, a CRYST1 box line. Output is a final PDB, a DCD trajectory and a CSV log of energies.

The key assumption is a fixed-charge force field with rigid bonds to hydrogen. Bonds do not break, and the PDB must already contain every hydrogen.

When to use it

• Short runs of a solvated structure when you want Python control of every step.
• Custom forces and integrators, and GPU speed in a native install.
• Checking a structure with two force fields.

Alternatives: GROMACS is faster for long runs of large biomolecular systems once you have a .tpr. LAMMPS handles non-biomolecular systems and custom pair potentials.

Install it yourself

OpenMM 8.6.1 from conda-forge. A dry run resolved openmm 8.6.1 (Python 3.14 build, 270 MB in total):

conda install -c conda-forge openmm=8.6.1
python -m openmm.testInstallation

The second command lists the platforms OpenMM can use. The minimal Python script follows the steps in the OpenMM user guide, and the helper makes the same calls:

import openmm as mm, openmm.app as app, openmm.unit as u
pdb = app.PDBFile("water.pdb")
ff = app.ForceField("amber14-all.xml", "amber14/tip3p.xml")
system = ff.createSystem(pdb.topology, nonbondedMethod=app.PME,
                         nonbondedCutoff=0.9*u.nanometer, constraints=app.HBonds)
integ = mm.LangevinMiddleIntegrator(300*u.kelvin, 1/u.picosecond, 2*u.femtosecond)
sim = app.Simulation(pdb.topology, system, integ)
sim.context.setPositions(pdb.positions)
sim.minimizeEnergy()
sim.reporters.append(app.DCDReporter("traj.dcd", 500))
sim.step(5000)

Reading the output

The helper run_md.py reads the PDB, builds the system (PME with a 0.9 nm cutoff if the PDB has a box, no cutoff if it does not), constrains bonds to hydrogen, and uses a Langevin middle integrator with friction 1/ps. It runs on the CPU platform with two threads. It minimises for at most 200 iterations, assigns velocities at the target temperature using the seed, then runs. It writes a DCD frame and a CSV row every report_every steps.

The example box is 20 A on a side (CRYST1 20.000) with 774 atoms, which is 258 water molecules. At 2 fs per step, 5000 steps is 10 ps.

• md_log.csv (955 bytes) has 10 rows, at steps 500 to 5000. Columns: step, time in ps, potential energy and kinetic energy in kJ/mole, temperature in K, speed in ns/day.
• At 1 ps the potential energy is -10515 kJ/mole and the temperature 255.7 K. At 10 ps they are -10238 kJ/mole and 313.5 K. Frames in between sit between 292 and 307 K.
• trajectory.dcd (93956 bytes) holds 10 frames. final.pdb (62846 bytes) holds the last coordinates, wrapped into the box.
• openmm_stdout.txt repeats step, potential energy and temperature every 500 steps and ends with done.

The temperature of 774 atoms fluctuates by about 10 K from frame to frame. Judge the average over many frames.

Your own runs

dr() { docker run --rm --user "$(id -u):$(id -g)" -e HOME=/tmp -v "$PWD/in:/work/in:ro" -v "$PWD/out:/work/out" -w /work/out "$@"; }
dr biosimhub/openmm:8.6.1 python /opt/driver/run_md.py /work/in/structure.pdb --forcefield charmm36-tip3p --steps 5000

Each run below used 5000 steps and was averaged over its last five frames.

1. Default (amber14-tip3p, 300 K, seed 42): mean 301.4 K, potential energy -10353.9 kJ/mole.
2. --forcefield charmm36-tip3p: 299.4 K and -10373.7 kJ/mole.
3. --temperature 350: 348.7 K and -9735.2 kJ/mole. The potential energy rose by 619 kJ/mole.
4. --seed 7: 298.9 K and -10309.3 kJ/mole.

The seed alone moved the potential energy by 45 kJ/mole, which is more than the 20 kJ/mole gap between the two force fields. Compare force fields only over longer runs with several seeds.

Pitfalls

• OpenMM needs a structure with hydrogens. Add every hydrogen before you upload.
• Without a CRYST1 line the helper runs with no cutoff and no periodic box. A bulk liquid needs the periodic box, so put it in CRYST1.
• Energies are in kJ/mole and the box is in nm inside OpenMM. PDB files use angstroms.
• The first frame is not at equilibrium. The default log reads 255.7 K at 1 ps. Discard early frames.
• Time step above 2 fs needs the bond-to-hydrogen constraints that the helper turns on. Check that the energy stays stable at larger steps.

Further reading

• OpenMM user guide, running simulations: https://docs.openmm.org/latest/userguide/application/02_running_sims.html
• Getting started and installation: https://docs.openmm.org/latest/userguide/application/01_getting_started.html
• Cookbook: https://openmm.github.io/openmm-cookbook/latest/index.html
• Source: https://github.com/openmm/openmm
• Eastman P et al. (2017) PLoS Comput Biol 13:e1005659. https://doi.org/10.1371/journal.pcbi.1005659