Free energy calculations for atomic solids through the Einstein crystal/molecule methodology using GROMACS and LAMMPS

In this work the free energy of solid phases is computed for the Lennard-Jones potential and for a model of NaCl. The free energy is evaluated through the Einstein crystal/molecule methodologies using the Molecular Dynamics programs: GROMACS and LAMMPS. The obtained results are compared with the results obtained from Monte Carlo. Good agreement between the different programs and methodologies was found. The procedure to perform the free energy calculations for the solid phase in the Molecular Dynamic programs is described. Since these programs allow to study any continuous intermolecular potential (when given in a tabulated form) this work shows that for isotropic potentials (describing for instance atomic solids or colloidal particles) free energy calculations can be performed on a routinely basis using GROMACS and/or LAMMPS.

💡 Research Summary

This paper presents a significant methodological advance by demonstrating that the absolute free energy of solids can be reliably calculated using two widely adopted, open-source molecular dynamics (MD) software packages: GROMACS and LAMMPS. The study bridges a gap between specialized, custom-built Monte Carlo (MC) codes and general-purpose MD simulators for a computationally challenging task.

The core of the work is the adaptation of the established Einstein crystal and Einstein molecule methodologies for free energy calculation into the workflow of GROMACS and LAMMPS. These methods compute the Helmholtz free energy of a solid by defining a reversible thermodynamic integration path from a reference system—an ideal Einstein crystal where particles are tethered to their lattice sites with harmonic springs—to the real solid of interest. The authors detail the practical implementation steps within the MD frameworks, leveraging built-in features such as position restraints/constraints and thermostats capable of correctly handling harmonic oscillators (e.g., the Bussi thermostat).

To validate their approach, the authors performed rigorous tests on two types of model systems with isotropic potentials. First, atomic solids interacting via a Lennard-Jones potential, tested with different truncation schemes (shifted and non-shifted) and system sizes (256 and 1372 particles). Second, an ionic solid modeled by the Joung-Cheatham potential for NaCl, which includes both Lennard-Jones and Coulombic interactions. For each case, the free energy results obtained from GROMACS and LAMMPS simulations were compared against benchmark results from their own MC code.

The results, summarized in a comprehensive table, show excellent agreement across all platforms and methodologies. The computed free energies from MD simulations are within 0.03 NkBT of the MC results, which is within the typical statistical uncertainty of such calculations. This confirms that the proposed MD implementation is both accurate and robust. The paper also provides practical recommendations for reliable calculations, suggesting the use of relatively large system sizes (over 1000 particles), extensive sampling, and a large potential cut-off distance.

The broader impact of this work is substantial. By enabling free energy calculations for solids within GROMACS and LAMMPS, it makes this powerful technique accessible to a much wider community of researchers who are already familiar with these tools. Since both packages support tabulated numerical potentials, the methodology is applicable to any continuous isotropic potential, paving the way for routine free energy studies of atomic solids, colloidal particles, and similar systems. The authors conclude by suggesting future work to explore the extension of this approach to molecular solids.