Much of science can be explained by the movement and interaction of molecules. Molecular dynamics (MD) is a computational technique used to explore these phenomena, from noble gases to biological macromolecules. Molly.jl is a pure Julia package for MD, and for the simulation of physical systems more broadly.
At the minute the package is a proof of concept for MD in Julia. It is not production ready. It can simulate a system of atoms with arbitrary interactions as defined by the user. Implemented features include:
- Interface to allow definition of new forces, simulators, thermostats, neighbor finders, loggers etc.
- Read in pre-computed Gromacs topology and coordinate files with the OPLS-AA forcefield and run MD on proteins with given parameters. In theory it can do this for any regular protein, but in practice this is untested.
- Non-bonded interactions - Lennard-Jones Van der Waals/repulsion force, electrostatic Coulomb potential, gravitational potential, soft sphere potential, Mie potential.
- Bonded interactions - covalent bonds, bond angles, torsion angles.
- Andersen thermostat.
- Velocity Verlet and velocity-free Verlet integration.
- Explicit solvent.
- Periodic boundary conditions in a cubic box.
- Neighbor list to speed up calculation of non-bonded forces.
- Automatic multithreading.
- GPU acceleration on CUDA-enabled devices.
- Run with Float64 or Float32.
- Some analysis functions, e.g. RDF.
- Physical agent-based modelling.
- Visualise simulations as animations.
- Differentiable molecular simulation on an experimental branch - see the relevant docs.
Features not yet implemented include:
- Protein force fields other than OPLS-AA.
- Water models.
- Energy minimisation.
- Other temperature or pressure coupling methods.
- Cell-based neighbor list.
- Protein preparation - solvent box, add hydrogens etc.
- Trajectory/topology file format readers/writers.
- Quantum mechanical modelling.
- High test coverage.
Julia is required, with Julia v1.5 or later required to get the latest version of Molly.
Install Molly from the Julia REPL.
Enter the package mode by pressing ]
and run add Molly
.
Some examples are given here, see the documentation for more on how to use the package.
Simulation of a Lennard-Jones gas:
using Molly
n_atoms = 100
box_size = 2.0 # nm
temp = 298 # K
mass = 10.0 # Relative atomic mass
atoms = [Atom(mass=mass, σ=0.3, ϵ=0.2) for i in 1:n_atoms]
coords = [box_size .* rand(SVector{3}) for i in 1:n_atoms]
velocities = [velocity(mass, temp) for i in 1:n_atoms]
general_inters = (LennardJones(),)
s = Simulation(
simulator=VelocityVerlet(),
atoms=atoms,
general_inters=general_inters,
coords=coords,
velocities=velocities,
temperature=temp,
box_size=box_size,
thermostat=AndersenThermostat(1.0),
loggers=Dict("temp" => TemperatureLogger(100)),
timestep=0.002, # ps
n_steps=10_000
)
simulate!(s)
Simulation of a protein:
using Molly
timestep = 0.0002 # ps
temp = 298 # K
n_steps = 5_000
atoms, specific_inter_lists, general_inters, nb_matrix, coords, box_size = readinputs(
joinpath(dirname(pathof(Molly)), "..", "data", "5XER", "gmx_top_ff.top"),
joinpath(dirname(pathof(Molly)), "..", "data", "5XER", "gmx_coords.gro"))
s = Simulation(
simulator=VelocityVerlet(),
atoms=atoms,
specific_inter_lists=specific_inter_lists,
general_inters=general_inters,
coords=coords,
velocities=[velocity(a.mass, temp) for a in atoms],
temperature=temp,
box_size=box_size,
neighbor_finder=DistanceNeighborFinder(nb_matrix, 10),
thermostat=AndersenThermostat(1.0),
loggers=Dict("temp" => TemperatureLogger(10),
"writer" => StructureWriter(10, "traj_5XER_1ps.pdb")),
timestep=timestep,
n_steps=n_steps
)
simulate!(s)
The above 1 ps simulation looks something like this when you view it in VMD:
Contributions are very welcome - see the roadmap issue for more.