Run simulation

The mdrun program executes molecular dynamics simulations using the DROPPS/CGPS OpenMM-based backend. It consumes a runtime file generated by dps grompp and produces:

  • trajectory files (XTC)

  • energy & log files

  • checkpoint files

  • final wrapped configuration (PDB)

This is analogous to gmx mdrun in GROMACS.

Implemented in mdrun.py.

Overview

dps mdrun:

  1. Loads all runtime data (system, topology, integrator parameters, positions) from a .tpr pickle generated by dps grompp.

  2. Sets up the OpenMM integrator (Langevin).

  3. Automatically selects the best available compute platform: CUDA → OpenCL → CPU.

  4. Applies optional: - energy minimization - velocity initialization - center-of-mass motion removal - pressure coupling (MC barostat)

  5. Configures reporters for: - screen log - file log - trajectory (XTC) - checkpoint

  6. Runs warming phase (if enabled).

  7. Runs production dynamics.

  8. Writes the final wrapped PDB file to disk.

Usage

Run a simulation:

dps mdrun -s run.tpr -o output/run

Resume from a checkpoint:

dps mdrun -s run.tpr -cpi output/run.chk -o output/run2

Arguments

Required

-s, --run-input FILE

Runtime file produced by dps grompp (Python pickle format).

-o, --output-prefix PREFIX

Prefix for all output files:

  • PREFIX.log

  • PREFIX.xtc

  • PREFIX.chk

  • PREFIX.pdb

Optional

-cpi, --checkpoint FILE

Checkpoint file from a previous run. If provided:

  • no minimization is performed

  • no warming is performed

  • velocities & box come from checkpoint

Simulation Workflow

### 1. Load runtime data

with open(args.run_input, "rb") as f:
    loaded = pickle.load(f)

Loaded fields include:

  • parameters (MDP dictionary)

  • mdsystem (OpenMM System)

  • mdtopology

  • positions

  • pdb_raw (for final PDB writing)

### 2. Set up Langevin integrator

temperature = parameters["production_temperature"] * kelvin
friction    = parameters["friction"] / picosecond
timestep    = parameters["dt"] * picosecond

### 3. Automatic platform selection

Preferred order:

['CUDA', 'OpenCL', 'CPU']

Platform-specific features:

  • CUDA: DeviceIndex=0, Precision=mixed

  • OpenCL: similar

  • CPU: no properties

Aborts if no platform is available.

### 4. Center-of-mass motion removal

If comm_mode = Linear:

mdsystem.addForce(openmm.CMMotionRemover(...))

Otherwise a warning is printed.

### 5. Pressure coupling (optional)

If pcoulp = True:

mdsystem.addForce(openmm.MonteCarloBarostat(pressure, temperature, tau_P))

### 6. Energy minimization

Performed only if:

  • minimize = True in parameters

  • no checkpoint is used

Done using:

simulation.minimizeEnergy(force_tol, max_steps)

Maximum force and RMS force are reported.

### 7. Reporters

Supported reporters:

  • screen log (printed to STDOUT)

  • file log (written to PREFIX.log)

  • XTC trajectory (PREFIX.xtc)

  • checkpoint files (PREFIX.chk)

Configuration is driven by parameters:

  • nst_screenlog

  • screenlog_grps

  • nst_filelog

  • filelog_grps

  • nst_xout

  • nst_cp

### 8. Velocity initialization

If gen_vel = True:

simulation.context.setVelocitiesToTemperature(initial_T)

### 9. Warming schedule (optional)

If:

  • initial temperature < production temperature

  • no checkpoint is used

Then temperature is raised gradually:

for T in range(T_initial, T_final):
    simulation.integrator.setTemperature(T)
    simulation.step(step_per_kelvin)

Temperature update every 10 K printed.

### 10. Production dynamics

Simulation length:

parameters["nsteps"]

Temperature is set to the final production temperature.

If a checkpoint is supplied:

simulation.loadCheckpoint(args.checkpoint)

### 11. Final wrapping and PDB writing

Final coordinates:

positions = simulation.context.getState(getPositions=True)
                 .getPositions(asNumpy=True).value_in_unit(nanometer)

Box vectors:

box_vec = [Lx, Ly, Lz]

Coordinates wrapped:

wrapped = positions % box_vec

Final PDB written with preserved metadata from pdb_raw:

write_pdb(PREFIX.pdb, box_vec, ...)