Crate dynamics

Crate dynamics 

Source
Expand description

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:

  • Atoms are divided into three categories: – Dynamic: Atoms that move – Static: Atoms that don’t move, but have mutual non-bonded interactions with dynamic atoms and water – 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/PME/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+barostat for the whole system. (Is water and dyn separate here?) todo

  • 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;
pub use params::ProtFFTypeChargeMap;
pub use params::ProtFfMap;

Modules§

integrate
Contains integration code, including the primary time step.
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
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.
PMEIndex
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§

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