Third-party backends#

LAMMPS#

Atooms provides a simulation backend for LAMMPS, an efficient and feature-rich molecular dynamics simulation package. The backend accepts a string variable containing regular LAMMPS commands and initial configuration to start the simulation. The latter can be provided in any of the following forms:

  • a System object

  • a Trajectory object

  • the path to an xyz trajectory

In the last two cases, the last configuration will be used to start the simulation.

Here we we use the first configuration of an existing trajectory for a Lennard-Jones fluid

import os
from atooms.core.utils import download
import atooms.trajectory as trj
from atooms.backends import lammps

# You can change it so that it points to the LAMMPS executable
lammps.lammps_command = 'lmp'

download('https://framagit.org/atooms/atooms/raw/master/data/lj_N1000_rho1.0.xyz', "/tmp")
system = trj.TrajectoryXYZ('/tmp/lj_N1000_rho1.0.xyz')[0]
cmd = """
pair_style      lj/cut 2.5
pair_coeff      1 1 1.0 1.0  2.5
neighbor        0.3 bin
neigh_modify    check yes
timestep        0.002
"""
backend = lammps.LAMMPS(system, cmd)

We now wrap the backend in a simulation instance. This way we can rely on atooms to write thermodynamic data and configurations to disk during the simulation: we just add the write_config() and write_thermo() callbacks to the simulation. You can add your own functions as callbacks to perform arbitrary manipulations on the system during the simulation. Keep in mind that calling these functions causes some overhead, so avoid calling them at too short intervals.

from atooms.simulation import Simulation
from atooms.system import Thermostat
from atooms.simulation.observers import store, write_config

# We create the simulation instance and set the output path
sim = Simulation(backend, output_path='/tmp/lammps.xyz')
# Write configurations every 1000 steps in xyz format
sim.add(write_config, 1000, trajectory_class=trj.TrajectoryXYZ)
# Store thermodynamic properties every 500 steps
sim.add(store, 100, ['steps', 'potential energy per particle', 'temperature'])

We add a thermostat to keep the system temperature at T=2.0 and run the simulations for 10000 steps.

backend.system.thermostat = Thermostat(temperature=2.0, relaxation_time=0.1)
sim.run(4000)

Note that we use atooms Thermostat object here: the backend will take care of adding appropriate commands to the LAMMPS script.

We have a quick look at the kinetic temperature as function of time to make sure the thermostat is working

plt.plot(sim.data['steps'], sim.data['temperature'])
plt.xlabel('Steps')
plt.ylabel('Temperature')
_images/lammps.png

We can then use the postprocessing package to compute the radial distribution function or any other correlation function from the trajectory.

RUMD#

There is native support for an efficient MD molecular dynamics code running entirely on GPU called RUMD, developed by the Glass and Time group in Roskilde. It is optimized for small and medium-size systems.

Here we pick the last frame of the trajectory, change the density of the system to unity and write this new configuration to a trajectory format suitable for RUMD

from atooms.trajectory import Trajectory

with Trajectory('/tmp/lj_N1000_rho1.0.xyz') as trajectory:
    system = trajectory[-1]
    system.density = 1.0
    print('New density:', round(len(system.particle) / system.cell.volume, 2))

from atooms.trajectory import TrajectoryRUMD
with TrajectoryRUMD('rescaled.xyz.gz', 'w') as trajectory:
    trajectory.write(system)

New density: 1.0

Now we run a short molecular dynamics simulation with the RUMD backend, using a Lennard-Jones potential:

import rumd
from atooms.backends.rumd import RUMD
from atooms.simulation import Simulation

potential = rumd.Pot_LJ_12_6(cutoff_method=rumd.ShiftedPotential)
potential.SetParams(i=0, j=0, Epsilon=1.0, Sigma=1.0, Rcut=2.5)
backend = RUMD('rescaled.xyz.gz', [potential], integrator='nve'
sim = Simulation(backend)
sim.run(1000)

A repository of interaction models for simple liquids and glasses is available in the atooms-models component package. It generates RUMD potentials automatically from standardized json file or Python dictionaries.