Getting Academics Coding in LAMMPS

Unofficial Thoughts on the Good, the Bad, and the Way Forward

Shern Ren Tee

AIBN and CTCMS, University of Queensland s.tee@uq.edu.au

26 Apr 2023

About Me

I did (PhD) weird DNA simulations, (2019-2020) weird protein simulations, and (2021-now) weird electrode-electrolyte simulations. Spot the pattern! Find out more at my (very poorly maintained) website.

The People I Work With

Prof Debra Bernhardt and the Bernhardt Group

The Turing Way

Research Software Engineers Au-NZ

Introduction: LAMMPS as a Community MD Code

Molecular Dynamics Today: Bigger…

(Alder and Wainwright 1959) Molecular dynamics in 1959 …

(Turoňová et al. 2020) … and today! (4.1 million atoms, 2.5 microseconds!)

Molecular Dynamics Today: … and weirder

A digital model of Birmingham

Comparison between real data and model output

A virtual pandemic with realistic geography and temporal mobility, modelled using molecular dynamics! (Alexiadis et al. 2021)

A Case Study: ELECTRODE for Conductive Molecular Dynamics

Capacitor modelling with dynamically-responsive charges using the USER-CONP2 package (Tee 2021)

Calculating the capacitance of molecular nanotubes against theory using the ELECTRODE package (Ahrens-Iwers et al. 2022)

Conductive electrode molecular dynamics is implemented in other packages (GROMACS, OpenMM, MetalWalls, …) but so far only the LAMMPS packages are being widely used outside their developing group!

Method Development Cycle in Molecular Dynamics

  • Implementation
  • Benchmarking
  • Documentation
  • Distribution

LAMMPS is a great community code for this cycle!

A Brief LAMMPS History (Thompson et al. 2022)

  • Mid-1990s:
    • Collab (two US DoE labs and three companies)
    • free Fortran code with license agreement
    • 10 years … 100 downloads
  • 2004:
    • C++ GPL code, 50k lines of code
    • downloaded more times in first month than in previous ten years
  • Today:
    • 1 million lines of code, several hundred contributors
    • 20-30k downloads per year

Lessons: ASAP and AMAP

Lessons Learned from LAMMPS (Plimpton 2019)

ASAP: Make your code as simple as possible to understand and extend

AMAP: Enable it to be used in as many as possible ways

Licensing: LGPL can be a good compromise

See also: (Plimpton and Gale 2013)

ASAP and AMAP in LAMMPS

  • Drop a .cpp and .h into lammps/src and compile, and you get an add-on!
  • LAMMPS can be called as a library:
    • Interfaces: C-style, Python, SWIG
    • Internally-implemented fix external
    • Caller can start multiple LAMMPS instances and coordinate messages between them

AMAP: LAMMPS Inputs As Workflows

Other packages: “turn these dials”

LAMMPS: “put these dials on your system and turn them”

(For example, a script won’t run any steps without run and won’t move any atoms without an integrator like fix nvt!)

Resembles “fluent interface” philosophy

in.peptide
# Solvated 5-mer peptide

pair_style      lj/charmm/coul/long 8.0 10.0 10.0
bond_style      harmonic
angle_style     charmm
dihedral_style  charmm
improper_style  harmonic
kspace_style    pppm 0.0001

read_data   data.peptide

neighbor    2.0 bin
timestep    2.0

fix     1 all nvt temp 275.0 275.0 100.0 tchain 1
fix     2 all shake 0.0001 10 100 b 4 6 8 10 12 14 18 a 31

dump    1 all atom 10 atom.peptide

run     300

ASAP 1: Simple, Clear Coding

LAMMPS uses “C with classes”:

  • Very little overloading, templating, STL (slowly increasing)
  • Low-level structs and kernels are C-style (e.g. n-dim arrays)
  • Detailed comments wherever needed

ASAP 2-4: Straightforward extensibility

What do you think this code does?

fix_jumble.cpp
// ...
int FixJumble::setmask()
{
  int mask = 0;
  mask |= FixConst::PRE_REVERSE;
  // enum'ed in parent class Fix
  return mask;
}

void FixJumble::pre_reverse(int, int)
{
  double f** = atom->f; // forces
  for (int i = 0; i < nlocal; ++i)
  {
    f[i][0] += rand(); // x-direction
    f[i][1] += rand(); // y
    f[i][2] += rand(); // z
  }
}
// ...

ASAP 2: Robust core, clear “sockets”

Timestep loop:

  • timestep initialization
  • fix->initial()
  • (sometimes) rebuild neighbor list
  • or send ghost atoms
  • force initialization
  • fix->pre_force()
  • compute forces
  • fix->pre_reverse()
  • receive ghost forces
  • fix->post_force()
  • finalize and output
fix_jumble.cpp
// ...
int FixJumble::setmask()
{
  int mask = 0;
  mask |= FixConst::PRE_REVERSE;
  // enum'ed in parent class Fix
  return mask;
}

void FixJumble::pre_reverse(int, int)
{
  double f** = atom->f; // forces
  for (int i = 0; i < nlocal; ++i)
  {
    f[i][0] += rand(); // x-direction
    f[i][1] += rand(); // y
    f[i][2] += rand(); // z
  }
}
// ...

ASAP 3: Good parenting, easy children

Parent “styles” define interfaces:

  • pair for MD potentials
  • compute for diagnostics (temperature, pressure)
  • fix for doing anything

Children styles inherit interface from parents; just override specific parent virtual functions as needed

fix_jumble.cpp
// ...
int FixJumble::setmask()
{
  int mask = 0;
  mask |= FixConst::PRE_REVERSE;
  // enum'ed in parent class Fix
  return mask;
}

void FixJumble::pre_reverse(int, int)
{
  double f** = atom->f; // forces
  for (int i = 0; i < nlocal; ++i)
  {
    f[i][0] += rand(); // x-direction
    f[i][1] += rand(); // y
    f[i][2] += rand(); // z
  }
}
// ...

ASAP 3: Parenting Pointers

Pointers class makes shared data “quasi-static”:

pointers.h
class Pointers {
 protected:
  LAMMPS *lmp;     // ...
  Atom *&atom;     // ...

 public:
  Pointers(LAMMPS *ptr) :
   lmp(ptr),       // ...
   atom(ptr->atom),// ...
}

Now anything inheriting Pointers auto-gets references to the necessary data contained in lmp!

fix_jumble.cpp
// ...
int FixJumble::setmask()
{
  int mask = 0;
  mask |= FixConst::PRE_REVERSE;
  // enum'ed in parent class Fix
  return mask;
}

void FixJumble::pre_reverse(int, int)
{
  double f** = atom->f; // forces
  for (int i = 0; i < nlocal; ++i)
  {
    f[i][0] += rand(); // x-direction
    f[i][1] += rand(); // y
    f[i][2] += rand(); // z
  }
}
// ...

ASAP 4: Style Factories

Factory pattern creates a derived FixJumble and returns a Fix pointer to LAMMPS:

modify.cpp
// ...
#define FIX_CLASS

#define FixStyle(_key, _Class)
   (*fix_map)[#_key] 
     = &style_creator<Fix, _Class>;

// C/Make'd list of fixes
#include "style_fix.h" // IWYU pragma: keep
#undef FixStyle
#undef FIX_CLASS
// ...
fix_jumble.h
// ...
#ifdef FIX_CLASS
FixStyle(jumble,FixJumble);
#else
// ...
fix_jumble.cpp
// ...
int FixJumble::setmask()
{
  int mask = 0;
  mask |= FixConst::PRE_REVERSE;
  // enum'ed in parent class Fix
  return mask;
}

void FixJumble::pre_reverse(int, int)
{
  double f** = atom->f; // forces
  for (int i = 0; i < nlocal; ++i)
  {
    f[i][0] += rand(); // x-direction
    f[i][1] += rand(); // y
    f[i][2] += rand(); // z
  }
}
// ...

AMAP: Making LAMMPS a Library

  • LAMMPS_NS namespace prevents code collisions
  • Global (static) variables replaced with Pointer class
  • MPI_COMM_WORLD replaced with communicator used to initialize LAMMPS
main.cpp
LAMMPS *lammps = new LAMMPS(argc, argv, lammps_comm);
lammps->input->file();
delete lammps;

Challenges and the Future

Successes …

User almost completely modified code on their own, made a suitable test, and was just missing some technical details on atom indexing in LAMMPS:

… and Failures

User had an idea for a complicated compute to be used in evaluating a potential, and was asking basic and strange questions a year after his original post:

Deep Thought: Programming as Theory-Building

Programming isn’t about writing the code; it’s about understanding the problem and expressing that understanding through code. (emphasis added)

… Why does having a theory of the program matter? Because this enables rapid and effective modification of the program to respond to changing requirements without piling up technical debt or hacks. (C J Silverio), (Peter Naur)

Diagram of “neighbor-listing”, a fundamental molecular dynamics operation and a frequent confuser of LAMMPS coders, from (Thompson et al. 2022)

The Future: LAMMPS as a Coding Community?

If LAMMPS is a program and not just code, then its continued survival and growth relies on introducing people to the theory of the program of molecular dynamics in LAMMPS!

This illustration is created by Scriberia with The Turing Way community. Used under a CC-BY 4.0 licence. DOI: 10.5281/zenodo.3332807

Acknowledgements

Prof Debra Bernhardt and the Bernhardt Group

The Turing Way

Research Software Engineers Au-NZ

References

Ahrens-Iwers, Ludwig J. V., Mathijs Janssen, Shern R. Tee, and Robert H. Meißner. 2022. ELECTRODE: An electrochemistry package for atomistic simulations.” The Journal of Chemical Physics 157 (8). https://doi.org/10.1063/5.0099239.
Alder, B. J., and T. E. Wainwright. 1959. Studies in Molecular Dynamics. I. General Method.” The Journal of Chemical Physics 31 (2): 459–66. https://doi.org/10.1063/1.1730376.
Alexiadis, A., A. Albano, A. Rahmat, M. Yildiz, A. Kefal, M. Ozbulut, N. Bakirci, D. A. Garzón-Alvarado, C. A. Duque-Daza, and J. H. Eslava-Schmalbach. 2021. “Simulation of Pandemics in Real Cities: Enhanced and Accurate Digital Laboratories.” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 477 (2245): 20200653. https://doi.org/10.1098/rspa.2020.0653.
Plimpton, Steven J. 2019. “Lessons Learned from Developing LAMMPS.” Sandia National Lab.(SNL-NM), Albuquerque, NM (United States); https://www.osti.gov/servlets/purl/1639472.
Plimpton, Steven J., and Julian D. Gale. 2013. “Developing Community Codes for Materials Modeling.” Current Opinion in Solid State and Materials Science 17 (6): 271–76. https://doi.org/https://doi.org/10.1016/j.cossms.2013.09.005.
Tee, Shern Ren. 2021. Theory and Practice in Constant Potential Molecular Dynamics Simulations.” In Multiscale Modeling of Electrochemical Reactions and Processes. AIP Publishing LLC. https://doi.org/10.1063/9780735422377_004.
Thompson, Aidan P., H. Metin Aktulga, Richard Berger, Dan S. Bolintineanu, W. Michael Brown, Paul S. Crozier, Pieter J. in ’t Veld, et al. 2022. “LAMMPS - a Flexible Simulation Tool for Particle-Based Materials Modeling at the Atomic, Meso, and Continuum Scales.” Computer Physics Communications 271: 108171. https://doi.org/https://doi.org/10.1016/j.cpc.2021.108171.
Turoňová, Beata, Mateusz Sikora, Christoph Schürmann, Wim J. H. Hagen, Sonja Welsch, Florian E. C. Blanc, Sören von Bülow, et al. 2020. “In Situ Structural Analysis of SARS-CoV-2 Spike Reveals Flexibility Mediated by Three Hinges.” Science 370 (6513): 203–8. https://doi.org/10.1126/science.abd5223.