Crate dynamics

Crate dynamics 

Source
Expand description

See the https://github.com/David-OConnor/dynamics/blob/main/README.md for a general overview. The textual information here is informal, and aimed at code maintenance; not library use.

This module contains a traditional molecular dynamics approach.

Good article A summary paper

Amber Force Fields reference Small molucules using GAFF2 Amber RM 2025

To download .dat files (GAFF2), download Amber source (Option 2) here. Files are in dat -> leap -> parm

Base units: Å, ps (10^-12), Dalton (AMU), native charge units (derive from other base units; not a traditional named unit).

We are using f64, and CPU-only for now, unless we confirm f32 will work. Maybe a mixed approach: Coordinates, velocities, and forces in 32-bit; sensitive global reductions (energy, virial, integration) in 64-bit.

We use Verlet integration. todo: Velocity verlet? Other techniques that improve and build upon it?

Amber: ff19SB for proteins, gaff2 for ligands. (Based on recs from https://ambermd.org/AmberModels.php).

We use the term “Non-bonded” interactions to refer to Coulomb, and Lennard Interactions, the latter of which is an approximation for Van der Waals force.

§A broad list of components of this simulation:

  • Water: Rigid OPC water molecules that have mutual non-bonded interactions with dynamic atoms and water
  • Thermostat/barostat, with a way to specify temp, pressure, water density
  • OPC water model
  • Cell wrapping
  • Verlet integration (Water and non-water)
  • Amber parameters for mass, partial charge, VdW (via LJ), dihedral/improper, angle, bond len
  • Optimizations for Coulomb: Ewald/SPME.
  • Optimizations for LJ: Dist cutoff for now.
  • Amber 1-2, 1-3 exclusions, and 1-4 scaling of covalently-bonded atoms.
  • Rayon parallelization of non-bonded forces
  • WIP SIMD and CUDA parallelization of non-bonded forces, depending on hardware availability. todo
  • A thermostat and barostat
  • An energy-measuring system.

A timing test, using bond-stretching forces between two atoms only. Measure the period of oscillation for these atom combinations, e.g. using custom Mol2 files. c6-c6: 35fs (correct). os-os: 47fs nc-nc: 34fs hw-hw: 9fs Our measurements, 2025-08-04 c6-c6: 35fs os-os: 31fs nc-nc: 34fs (Correct) hw-hw: 6fs


We use traditional MD non-bonded terms to maintain geometry: Bond length, valence angle between 3 bonded atoms, dihedral angle between 4 bonded atoms (linear), and improper dihedral angle between each hub and 3 spokes. (E.g. at ring intersections). We also apply Coulomb force between atom-centered partial charges, and Lennard Jones potentials to simulate Van der Waals forces. These use spring-like forces to retain most geometry, while allowing for flexibility.

We use the OPC water model. (See water_opc.rs). For both maintaining the geometry of each water molecule, and for maintaining Hydrogen atom positions, we do not apply typical non-bonded interactions: We use SHAKE + RATTLE algorithms for these. In the case of water, it’s required for OPC compliance. For H, it allows us to maintain integrator stability with a greater timestep, e.g. 2fs instead of 1fs.

On f32 vs f64 floating point precision: f32 may be good enough for most things, and typical MD packages use mixed precision. Long-range electrostatics are a good candidate for using f64. Or, very long runs.

Note on performance: It appears that non-bonded forces dominate computation time. This is my observation, and it’s confirmed by an LLM. Both LJ and Coulomb take up most of the time; bonded forces are comparatively insignificant. Building neighbor lists are also significant. These are the areas we focus on for parallel computation (Thread pools, SIMD, CUDA)

Re-exports§

pub use integrate::Integrator;

Modules§

integrate
Contains integration code, including the primary time step.
minimize_energy
params
For Amber and other parameters.
snapshot
Related to storing snapshots (also known as trajectories) of MD runs.

Structs§

AtomDynamics
A trimmed-down atom for use with molecular dynamics. Contains parameters for single-atom, but we use ParametersIndex for multi-atom parameters.
Dihedral
An amino acid in a protein structure, including all dihedral angles required to determine the conformation. Includes backbone and side chain dihedral angles. Doesn’t store coordinates, but coordinates can be generated using forward kinematics from the angles.
ForcesOnWaterMol
Per-water, per-site force accumulator. Used transiently when applying nonbonded forces. This is the force on each atom in the molecule.
MdConfig
MdOverrides
These are primarily used for debugging and testing, but may be used for specific scenarios as well, e.g. if wishing to speed up computations for real-time use by removing long range forces.
MdState
MolDynamics
Packages information required to perform dynamics on a Molecule. This is used to initialize the simulation with atoms and related; one or more of these is passed at init.
ParamError
Represents problems loading parameters. For example, if an atom is missing a force field type or partial charge, or has a force field type that hasn’t been loaded.

Enums§

ComputationDevice
FfMolType
This is used to assign the correct force field parameters to a molecule.
HydrogenConstraint
We use this variant in the configuration API. Deferrs to HydrogenConstraintInner for holding constraints.
SimBoxInit
Note: The shortest edge should be > 2(r_cutoff + r_skin), to prevent atoms from interacting with their own image in the real-space component.

Functions§

find_planar_posit
Given two positions, find the third planar position. Uses the central_posit_0 distance for the computed distances.
find_tetra_posit_final
Given 3 satellite atoms, find the 4th, in tetrahedral config.
find_tetra_posits
Given two satelite positions, find the third and fourth tetrahedral position. Uses the central, posit_0 distance for the computed distances.
load_snapshots
merge_params
Add items from one parameter set to the other. If there are duplicates, the second set’s overrides the baseline.
populate_hydrogens_dihedrals
Adds hydrogens to a molecule, and populdates residue dihedral angles. This is useful in particular for mmCIF files from RCSB PDB, as they don’t have these. Uses Amber (or similar)-provided parameters as a guide.
save_snapshots