A fast and accurate molecular dynamics engine built in Zig that simulates phospholipid bilayers using the Martini force field - the gold standard for coarse-grained biological membrane simulations. This tool helps researchers understand how cell membranes behave and interact, which is crucial for drug delivery, membrane fusion, and understanding cellular processes. Largely a learning project to understand Zig. Formerly created around PAIR.camp 2024, this work prepares for a port to CUDA by building out a more accurate simulation framework that is run on CPU.
Think of this as a "digital microscope" that lets you watch how lipid molecules move and interact in real-time. Instead of simulating every single atom (which would be impossibly slow), we group atoms together into "beads" using the Martini force field - like zooming out on a map. This gives us the big picture of how membranes behave while still being scientifically accurate. The goal of this project though in this stage is largely educational, although in future perhaps I may build it out such that it may be useful generally either in computational efficiency or some other manner.
- Realistic lipid molecules: Each phospholipid has a water-loving head and two oil-loving tails, just like real cell membranes
- Fast simulation: Optimized to run thousands of steps per second on your computer
- Accurate physics: Uses the Martini force field with real molecular forces and biological temperatures
- Easy visualization: Outputs files you can view with standard molecular visualization software
We simulate three main types of forces:
- Attraction/Repulsion: Molecules push each other away when too close, pull together when at the right distance
- Bonds: Keep the head and tails of each lipid molecule connected
- Temperature: Keeps everything moving at body temperature (37°C)
Each phospholipid molecule is made up of:
- Head bead (4.7 Å): Represents the phosphate and choline groups - these love water
- Tail beads (4.7 Å each): Represent the fatty acid chains - these avoid water
- Bonds: Keep everything connected with realistic spring-like forces
- Zig compiler (version 0.13.0 or newer)
# Build the project
zig build
# Run a lipid bilayer simulation
zig build run
# Run tests to make sure everything works
zig build testWhen you run the simulation, you'll get:
- Console output: Shows progress, timing, and final statistics
- XYZ files: One file per time step showing where each molecule is
- Visualization: You can open these files in molecular viewers to see the bilayer
A lipid bilayer looks like a sandwich:
- Outer layers: Water-loving heads facing the water
- Middle layer: Oil-loving tails pointing toward each other
- Thickness: About 28 Å (very thin - about 1/1000th the width of a human hair!)
The simulation tracks:
- Bilayer thickness: How thick the membrane is
- Temperature: Whether it's at the right temperature
- Bond lengths: Whether molecules are staying connected properly
- Particle positions: Where each molecule is at each time step
You can visualize the results using:
# VMD (Visual Molecular Dynamics) - professional software
vmd -xyz frame_000000.xyz
# Ovito - free, user-friendly viewer
ovito frame_000000.xyz
# Python script
python3 visualize.py
# Note that visualize.py is largely deprecated, and it is reccomended to utilise series based visualization.- Q (Charged): Lipid heads - phosphate and choline groups (water-loving)
- P (Polar): Glycerol backbone and ester groups (water-loving)
- C (Apolar): Lipid tails - saturated hydrocarbon chains (oil-loving)
- N (Nonpolar): Cholesterol and unsaturated hydrocarbons
- O (Water): Water molecules (P4 type)
The current version runs about 3,500 steps per second on a typical computer. This means:
- A 1-second simulation shows 3,500 snapshots of molecular motion
- You can simulate realistic membrane behavior in minutes
- The simulation is stable and doesn't break down over time
The simulation produces realistic results:
- Bilayer thickness matches experimental measurements
- Molecules maintain proper connections
- Temperature stays near body temperature
- Structure remains stable over time
We use the Lennard-Jones potential to model molecular interactions:
Force = 24ε/r × [2(σ/r)¹² - (σ/r)⁶]
Where:
- ε (epsilon) = interaction strength (different for heads vs tails)
- σ (sigma) = particle size
- r = distance between particles
I used the Velocity Verlet algorithm, which:
- Updates positions using current velocities and forces
- Calculates new forces based on new positions
- Updates velocities using new forces
- Repeats thousands of times per second
// Martini particle dimensions (4:1 mapping)
Q0_DIAMETER: 4.7 Å // Charged particles (phosphate groups)
P1_DIAMETER: 4.7 Å // Polar particles (glycerol backbone)
C1_DIAMETER: 4.7 Å // Apolar particles (hydrocarbon tails)
BOND_LENGTH: 4.7 Å // Standard Martini bond length
// Masses (atomic mass units)
Q0_MASS: 72.0 amu // ~72 atoms per charged bead
P1_MASS: 72.0 amu // ~72 atoms per polar bead
C1_MASS: 72.0 amu // ~72 atoms per apolar bead
// Interaction strengths (Martini 2.2)
Q0_Q0_EPSILON: 5.0 // Strong electrostatic interaction
P_P_EPSILON: 4.0 // Strong polar interaction
C_C_EPSILON: 1.0 // Weak apolar interaction- Neighbor lists: Make it even faster for larger systems, this will be extended with CUDA in future.
- More lipid types: Different kinds of phospholipids, with an aim to simulate more complex lipid raft interactions.
- Proteins: Add membrane proteins to the simulation, both embedded between layers and into single layers.
- Fusion events: Study how membranes merge together, potentially including SNARE protein behaviour through (3).
- Better analysis: More tools to understand the results with visualisation available "in library".
This simulator could help with:
- Drug delivery: Understanding how drugs interact with cell membranes
- Membrane fusion: Studying how cells merge (important for fertilization, viral entry)
- Lipid rafts: Understanding specialized membrane domains
- Membrane proteins: How proteins interact with lipids
The project is organized into:
src/
├── root.zig # Main simulation engine
├── main.zig # Program that runs the simulation
└── visualize.py # Python script for basic visualization
Key components:
- MARTINI_CONSTANTS: Martini force field parameters for lipids
- Vec3: Fast 3D math operations
- Particle: Individual molecular beads
- Phospholipid: Complete lipid molecules
- SimulationBox: The main simulation container
We welcome contributions! Here are some areas where help would be great:
- Performance improvements: Make it run even faster
- New features: Add different types of molecules or interactions
- Documentation: Help make it easier for others to use
- Bug fixes: Find and fix any issues
- Testing: Add more tests to ensure accuracy
If you run into problems:
- Check that you have the right version of Zig
- Make sure all tests pass:
zig build test - Look at the console output for error messages
- Check that the XYZ files are being created properly
This project is open source under the MIT License. Feel free to use it for research, education, or commercial applications.
This simulator follows established practices in computational biophysics and membrane science, implementing the widely-validated Martini force field. It builds on decades of research into how to accurately model biological membranes at the molecular level. Many existing approaches have been created to model phospholipid interaction and this work aims to build on those areas in the long term by creating a more readily accessible out of the box environment for visualisation of these phenomena, starting from building out the backend CPU logic.