MDDEM

Molecular Dynamics / Discrete Element Method in Rust

Note: The initial working example of MDDEM was hand-written, but all code has since been touched and expanded by Claude Code. This project explores alternative coding patterns to LAMMPS, prioritizing ergonomics — composable plugins, dependency injections, typed configs, and a Rust-native API.

What is MDDEM?

MDDEM (pronounced like "Madem" without the 'a') is a particle simulation engine written in Rust. It supports both Discrete Element Method (DEM) for granular materials and Molecular Dynamics (MD) for continuous-potential systems like Lennard-Jones fluids.

The design is built around composability. A dependency-injection scheduler (inspired by Bevy) and plugin system let you assemble simulations from independent, reusable pieces. Physics models, integrators, neighbor lists, and output formats are all plugins. Systems declare what resources they need as function arguments; the scheduler injects them automatically and handles execution order.

Configuration follows a two-tier approach. Tier 1 is declarative TOML config for standard simulations — named fields, typed values, validated at startup. Tier 2 is the Rust API for complex simulations — main.rs composes plugins, and custom systems are real functions with full type safety and IDE autocomplete. Both tiers can be mixed: the hopper example uses TOML config with custom Rust systems for runtime wall control.

Why is MDDEM?

At first, I wanted to learn about LAMMPS communication more through rewriting it into Rust. I also have had many pain points with editing LAMMPS code, and working with LAMMPS scripts, and wanted to see if a scheduler with dependency injection would work for something like this (big fan of bevy).

Now that Claude code is good enough to debug MPI communication neighbor list problems (it still struggles a lot with these, don't we all), I have expanded the scope to hopefully be a nice starting place for anyone wanting to try and vibe code a MD or DEM simulation in rust. I would not trust this code for anything you want to publish. I also would not contribute to this code in a serious manual fashion. View it as a playground to test out AI agents for whatever work you're doing. That being said, it's producing reasonable physics results, at about 90% the performance of LAMMPS.

If you're unsure about why this is even a repo (I kinda agree with you), I would look at the hopper example first. It's very simple, but shows how easy it is to add your own code into the mix.

Installation

MDDEM is not on crates.io yet. Add it as a git dependency:

[dependencies]
mddem = { git = "https://github.com/SueHeir/MDDEM" }

To build with MPI support (default), you need an MPI implementation installed (e.g., OpenMPI, MPICH). To build without MPI:

[dependencies]
mddem = { git = "https://github.com/SueHeir/MDDEM", default-features = false }

To clone and work on MDDEM directly:

git clone https://github.com/SueHeir/MDDEM.git
cd MDDEM
cargo build --release

Quick Start

DEM: Granular Gas

main.rs

use mddem::prelude::*;

fn main() {
    let mut app = App::new();
    app.add_plugins(CorePlugins)
        .add_plugins(GranularDefaultPlugins);
    app.start();
}

config.toml

[domain]
x_high = 0.025
y_high = 0.025
z_high = 0.025
periodic_x = true
periodic_y = true
periodic_z = true

[neighbor]
skin_fraction = 1.1
bin_size = 0.005

[[dem.materials]]
name = "glass"
youngs_mod = 8.7e9
poisson_ratio = 0.3
restitution = 0.95
friction = 0.4

[[particles.insert]]
material = "glass"
count = 500
radius = 0.001
density = 2500.0
velocity = 0.5

[run]
steps = 10000
thermo = 100

MD: Lennard-Jones Fluid

main.rs

use mddem::prelude::*;

fn main() {
    let mut app = App::new();
    app.add_plugins(CorePlugins).add_plugins(LJDefaultPlugins);
    app.start();
}

config.toml

[domain]
x_high = 10.06
y_high = 10.06
z_high = 10.06
periodic_x = true
periodic_y = true
periodic_z = true

[neighbor]
skin_fraction = 1.12
bin_size = 1.0

[lattice]
style = "fcc"
density = 0.85
temperature = 0.85
mass = 1.0
skin = 1.25

[lj]
epsilon = 1.0
sigma = 1.0
cutoff = 2.5

[thermostat]
temperature = 0.85
coupling = 1.0

[run]
steps = 100000
thermo = 1000

Running

# Single process
cargo run --release -- config.toml

# With MPI (build first, then launch)
cargo build --release
mpiexec -n 4 ./target/release/my_simulation config.toml

# Print the compiled schedule (Graphviz DOT)
cargo run --release -- config.toml --schedule

What's Included

DEM (Granular)

• Hertz elastic normal contact with viscoelastic damping (LAMMPS hertz/material equivalent)
• Mindlin tangential spring-history with Coulomb friction cap
• Velocity Verlet for translational + rotational degrees of freedom
• Quaternion-based orientation tracking
• Automatic timestep (5% of Rayleigh wave period)
• Configurable gravity body force
• General plane wall contacts, toggleable at runtime
• Named material types with per-pair mixing

MD (Molecular)

• Lennard-Jones 12-6 with cutoff, virial accumulator, and tail corrections
• Nose-Hoover NVT thermostat (symmetric Liouville splitting)
• Langevin thermostat (stochastic friction + random force, group-aware)
• FCC lattice initialization with Maxwell-Boltzmann velocities
• Radial distribution function, mean square displacement, virial pressure

Infrastructure

• Optional 3D MPI domain decomposition with multi-hop ghost forwarding
• Single-process mode with ghost atoms for periodic boundaries
• Bin-based neighbor lists with CSR storage and forward-only stencil
• Brute force and sweep-and-prune neighbor lists also available
• Named atom groups ([[group]]) with type and region filters
• Fixes: AddForce, SetForce, Freeze, MoveLinear (group-based)
• TOML config with serde validation and deny_unknown_fields on all config structs
• Dump files (CSV/binary), restart files (bincode/JSON), VTP visualization
• Generic restart serialization — any registered AtomData extension is automatically saved/restored
• Multi-stage runs with per-stage output control
• #[derive(AtomData)] proc macro for zero-boilerplate per-atom extension structs

Examples

DEM examples include validate.py scripts for physics checks (Haff's law cooling, hopper settling). The lj_argon example validates against known liquid Argon properties (RDF, MSD, pressure) and generates diagnostic plots. Run ./validate.sh to execute all tests and validations.

Performance

Single-core LJ fluid benchmark comparing MDDEM to LAMMPS (29 Sep 2024 release). Identical physics: LJ 12-6 with cutoff 2.5 sigma, FCC lattice at rho*=0.8442, Nose-Hoover NVT at T*=1.44, neighbor rebuild every 20 steps, 200 timesteps. Compiled with --release on Apple M1 Pro. RDF/MSD disabled in both codes for fair comparison.

LAMMPS is ~1.1x faster at scale, with consistent O(N) scaling in both codes. The force loop (~63% of MDDEM runtime) uses LAMMPS-style precomputed constants with a single reciprocal per pair. The neighbor list build (~29% of runtime) uses CSR bins with a forward stencil, sorted position caches, sorted neighbor indices for sequential cache access, and unsafe bounds-check elimination. The DI scheduler caches downcast pointers at system entry so Res<T>/ResMut<T> access is a direct dereference with no dynamic dispatch in hot loops.

MPI benchmark (4 processes, 2x2x1 decomposition) on the same hardware:

* The 108-atom MPI case is excluded — with 2x2x1 decomposition, subdomains are smaller than the ghost cutoff, producing invalid physics (NaN energies, zero neighbors). Probably a bug.

MDDEM MPI achieves 2.4-3.5x speedup over single-core at scale, with the ratio to LAMMPS narrowing to 1.1x at 32k+ atoms. Spatial sorting of atoms by bin is enabled in both single-core and MPI modes, improving cache locality for force and neighbor list computations. Communication uses per-dimension exchange with non-blocking sends, multi-hop ghost forwarding when needed, and lightweight ghost position updates between neighbor rebuilds.

Roadmap

Planned features, organized by implementation wave:

Completed

1. Groups — Named atom subsets ([[group]]) for selective operations
2. Custom thermo computes — ThermoCompute trait + configurable thermo columns
3. Langevin thermostat — Stochastic friction + random force
4. AddForce / SetForce / Freeze / MoveLinear — Group-based atom manipulation fixes

Wave 3: Medium Complexity

5. Shrink-wrap boundaries — Auto-expanding domain per axis
6. Energy minimization — CG + FIRE, runs as a [[run]] stage
7. Fix ave/time — Running time averages to file

Wave 4: Box Manipulation

8. NPT barostat — Nose-Hoover pressure control with box rescaling
9. Fix deform — Box size change over time

Wave 5: New Physics

10. Short-range Coulomb — k_e*q_i*q_j/r^2 with cutoff (no PPPM)
11. Bonded particle model — Bond lists between DEM spheres (parallel bonds, breakage criteria) for modeling cemented/rock-like materials

Every feature ships with an example and validation against analytical solutions or LAMMPS.

Out of Scope

These are specialized features that won't be in core. Users can write plugins for them:

• EAM / Tersoff / many-body potentials
• ReaxFF / ML potentials (SNAP, ACE)
• Angles / dihedrals (bonds are planned — see roadmap)
• SHAKE / RATTLE constraints
• Units system (TOML config uses explicit units in docs)
• rRESPA multi-timescale integration
• PPPM / Ewald long-range electrostatics
• Rigid body dynamics
• GCMC, NEB, spin dynamics
• Triclinic (non-orthogonal) boxes
• Multiple dump formats (XYZ, DCD, NetCDF) — CSV + VTP + binary covers most needs

Architecture

A simulation is composed by adding plugin groups to an App:

use mddem::prelude::*;

fn main() {
    let mut app = App::new();
    app.add_plugins(CorePlugins).add_plugins(GranularDefaultPlugins);
    app.start();
}

CorePlugins bundles config loading, communication, domain decomposition, neighbor lists, Velocity Verlet, and output. GranularDefaultPlugins adds DEM atom data, insertion, contact forces, gravity, and walls. LJDefaultPlugins adds FCC lattice, LJ forces, thermostat, and measurements. Individual plugins can be added separately for custom configurations.

Testing

cargo test --workspace                    # All tests (with MPI)
cargo test --workspace --no-default-features  # All tests (without MPI)
cargo test -p dem_granular                # Single crate

./validate.sh                             # Full validation: tests + examples + physics checks
./validate.sh --long                      # Production-length runs with full physics validation
./validate.sh --dem                       # DEM examples only
./validate.sh --md                        # MD examples only

Unit tests cover: single-process communication, domain decomposition, Velocity Verlet integration, neighbor lists (all three algorithms), Hertz/Mindlin contact forces, rotational dynamics, material mixing, LJ force and virial, Nose-Hoover thermostat, FCC lattice insertion, and RDF/MSD measurement.

Physics validation scripts (validate.py per example) check simulation output against analytical solutions: Haff's cooling law for granular gas, settling behavior for hopper, and RDF/MSD/pressure against known liquid Argon properties for LJ fluid.

Contributing

Vibe-coded PRs are accepted, provided the contributor is qualified in the relevant domain and has personally reviewed the code being submitted. Report issues at github.com/SueHeir/MDDEM/issues.

License

MIT OR Apache-2.0