Argon Simulations with EM3, a New Modular Molecular Mechanics Program

Andrew Rohskopf

REVIEW 4 major objections 5 minor 1 references

Argon Simulations with EM3, a New Modular Molecular Mechanics Program

T0 review · 4 major / 5 minor · reviewed 2026-08-14 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read This paper reports a new modular molecular dynamics program, EM3, and argues that its Lennard-Jones argon simulations are validated by agreement with Monte Carlo and literature values to within about 10 percent, by energy-conserving NVE…

desk verdict A genuinely useful open-source educational MD code with credible validation, undermined by wrong equations and missing details that should be fixed before it can serve as a teaching reference. read the letter →

arxiv 1908.00601 v1 pith:2GENPYNX submitted 2019-08-01 physics.comp-ph

Andrew Rohskopf This is my paper · ORCID

classification physics.comp-ph

keywords moleculardynamicsLennard-Jonespotentialliquidargoncodevalidationobject-orientedC++AndersenthermostatEinsteinrelationself-diffusioncoefficient

verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This report describes a newly written object-oriented molecular dynamics program, EM3, and argues that the program is correct. The argument is a validation campaign on Lennard-Jones argon: solid-state runs conserve total energy in the NVE stage, hold temperature at the set point during Andersen-thermostatted equilibration, and reproduce the known cohesive energy, while three liquid states at reduced densities and temperatures $(n=0.5, T=5)$, $(n=0.9, T=2)$, and $(n=0.8, T=4)$ match Monte Carlo and published equation-of-state values within about 10% in potential energy and within a few percent in pressure. A self-diffusion coefficient for liquid argon at 158 K is then computed from the mean-square displacement, showing that the modular code can be extended beyond the equilibrium properties used for validation. If the validation is sound, EM3 is a working open-source MD program whose modular C++ structure makes it useful as a learning tool and as a base for future development.

What carries the argument

The load-bearing machinery is the modular class structure of EM3: an Input class, a Neighbor class using the minimum-image convention, a Potential class for the Lennard-Jones pair interaction, an Update class running velocity Verlet, a Compute class for kinetic energy, temperature, pressure, and mean-square displacement, plus memory, output, and timer utilities. The validation is carried by the standard identities these classes implement: the truncated Lennard-Jones potential $U=\sum_{j>i}4(r_{ij}^{-12}-r_{ij}^{-6})$, the force and virial-pressure formulas derived from it, the Einstein relation $\partial\langle r^2\rangle/\partial t = 2dD$ connecting MSD to diffusion, and energy conservation as a check that integration and neighbor-list bookkeeping are correct. The modular boundaries are what make the validation meaningful for future use: swapping the Potential or Compute class is the intended route to extending the code.

What would settle it

Re-run one of the validated liquid cases, such as $n=0.9$, $T=2$, from a fresh 500-atom fcc configuration with a longer production window, and require that the NVE total energy per atom show no drift above the statistical noise and that the production-window averages of potential energy and pressure stay within the reported error band; a visible energy drift or an average outside that band would falsify the validation claim.

Watch

Extended reading notes

Core claim

The central claim is that EM3 correctly implements a standard Lennard-Jones molecular dynamics pipeline and is thereby validated as a simulation tool. In the paper's own terms, the agreement with literature and Monte Carlo values within 10 percent, together with energy-conserving NVE dynamics and a functioning Andersen thermostat, is sufficient to validate the code. On that basis EM3 gives the correct potential energy and pressure for the three liquid-argon states, the correct cohesive energy for the 0 K fcc solid, qualitatively correct solid-state MSD behavior at 10, 30, and 50 K, and a statistically clean Einstein-relation diffusion coefficient ($D = 0.2641\,\sigma^2/\mathrm{ps}$) for argon at 158 K and $8\times10^{27}\,\mathrm{m^{-3}}$.

Load-bearing premise

The validation rests on assuming that the truncated Lennard-Jones model with a $3\sigma$ cutoff and no long-range correction is comparable to the Monte Carlo and equation-of-state references, so the roughly 10% agreement reflects a correct code rather than a cancellation of errors between two different models.

Editorial extensions

If this is right

If the validation holds, EM3 can be trusted to compute equilibrium properties of Lennard-Jones argon at states not explicitly tested, within the same cutoff and system-size conventions.
The reported energy conservation in NVE runs implies that the velocity Verlet integration and neighbor-list/minimum-image bookkeeping are error-free at the precision shown, so new potentials added through the same classes inherit a tested integrator.
The diffusion-coefficient calculation shows that the code can produce transport properties, not just static averages, by computing MSD during an NVE production run.
The modular, open-source structure means the program can serve as a pedagogical bridge: a user can trace how input, neighbor lists, potential, integration, and computed quantities connect in a few hundred lines of C++.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

Not in the paper: adding the analytic long-range correction to the three liquid averages would show whether the residual error shrinks, locating any remaining discrepancy in the cutoff truncation rather than in the forcing.
Not in the paper: the strongest test of the modularity claim would be to implement a second potential class, such as a shifted-force or embedded-atom-style potential, and repeat the same validation without touching the integrator.
Not in the paper: the same 158 K trajectories could also yield a velocity autocorrelation function and therefore an independent estimate of the diffusion coefficient, checking whether the Einstein-relation slope is robust.

Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Desk Editor's Note

read the letter

EM3 is a real, open-source educational MD program with the code and example inputs actually shipped, which puts it ahead of many code-description papers. The validation core is credible: for three liquid argon state points the MD energies and pressures agree with the Allen–Tildesley MC runs and the Johnson et al. EOS within about 7%, and the energy conservation and thermostat plots are consistent with honest simulation behavior. That is genuine evidence the force loop and integrator work. The diffusion coefficient at 158 K is a sensible demonstration of the code's extensibility.

The soft spots are real but mostly in the prose, not the simulation core. Equations 6 and 7 are not the velocity Verlet; they are missing position updates and dimensionally wrong, and they contradict the C++ in Figure 3. Equation 11 for the total pressure is dimensionally inconsistent (PV on one side, a stress on the other). The offset-tag description contradicts the example INPUT file: the text says setting offset to zero subtracts the cutoff energy, but Figure 4 shows offset: 1 for the 3 K run while the text says that run had the offset subtracted. For a paper whose stated goal is teaching, these errors are not cosmetic. Also missing: error bars on all reported averages, and any description of the MC reference runs' cutoff and system size, so 'within 10%' is a loose statement.

On the long-range correction point: the stress-test worry about 'cancellation of biases' does not survive contact with the actual numbers. Adding the standard tail corrections to the MD results would move them closer to Johnson et al., not further away; the discrepancy is simply the omitted tail, and 3σ cutoff makes it modest. The more defensible criticism is that the paper never declares what model the MC references actually solved.

Bottom line: this is a legitimate code paper, not a physics advance. It deserves a serious referee, and with the equations fixed and some uncertainty quantification it could be a decent teaching resource. I would not desk-reject it, but I would send it back for revision.

Referee Report

4 major / 5 minor

Summary. The paper presents EM3, a new modular C++ molecular mechanics program, and validates it on Lennard-Jones argon by comparing solid-state energy conservation, cohesive energy, liquid-state potential energy, and pressure against Monte Carlo simulations and the Johnson et al. equation of state, followed by a self-diffusion coefficient calculation. The central claim, stated in Section 2.3, is that agreement with literature and MC values within 10%, together with energy-conserving dynamics and a working thermostat, is sufficient to validate the code. The paper also emphasizes the educational value of the modular, open-source implementation and includes the source code and examples on GitHub.

Significance. If the central validation claim holds, EM3 could serve as a useful educational tool for newcomers to molecular simulation and object-oriented programming in computational physics. The open-source availability of the code, the reproducible input files, and the comparison against independent external benchmarks (Allen-Tildesley MC and Johnson et al.) are genuine strengths. However, the validation has important caveats: the truncated Lennard-Jones potential is compared with references whose cutoff treatment is not specified, and several equations in the manuscript are internally inconsistent or dimensionally wrong. These issues weaken confidence in the reported agreement, so the significance of the work is currently conditional on a careful revision of both the validation protocol and the presentation of the methods.

major comments (4)

[Section 1.5, Eqs. (6)-(7)] The equations describing the velocity Verlet integrator do not match the code in Figure 3. Eq. (6) should include the old position on the right-hand side: x_i^α(t+δt) = x_i^α(t) + δt v_i^α(t+δt/2). More seriously, Eq. (7) is dimensionally inconsistent: it equates a velocity to a product of a position and a velocity, whereas the correct final velocity update is v_i^α(t+δt) = v_i^α(t+δt/2) + (δt/2) a_i^α(t+δt). Since the velocity Verlet scheme is the core algorithm of the proposed program, this derivation error is load-bearing and should be corrected.
[Section 1.6, Eq. (11)] The pressure equation P V = N k_B T + (1/3)(P_xx+P_yy+P_zz) is dimensionally inconsistent because the stress components from Eq. (4) have units of pressure (energy per volume), not energy. The correct relation is P = (N k_B T)/V + (1/3)(P_xx+P_yy+P_zz), or the virial term should be expressed as the energy-like sum (1/V) multiplied by the stress components. Pressure is one of the primary validated quantities in Tables 3 and 4, so this equation needs to be fixed.
[Section 2.3, Tables 3-4] The 'agreement within 10%' is the central evidence for code validation, but the comparison is between EM3's truncated Lennard-Jones model (r_c = 3σ, no long-range corrections) and Johnson et al.'s full Lennard-Jones equation of state, while the MC runs' cutoff and system size are not stated. The paper itself acknowledges in Section 2.3 that long-range corrections could have been added but declines to do so. At the n=0.9, T=2 state point, the uniform-density tail corrections to energy and pressure are comparable to or larger than the reported 4.8% and 5.1% discrepancies, so the apparent agreement may result from a cancellation of truncation effects rather than from correct force and pressure computation. The authors should either add long-range corrections to the EM3 averages or explicitly demonstrate that the MC reference uses the same truncation and system size as EM3, and report those details. Without this, the validation remains conditional.
[Section 1.5, Eq. (8)] The thermostat described by Eq. (8) is a deterministic velocity rescaling applied every timestep. This is not the Andersen thermostat, and it does not generate a canonical NVT distribution because it fixes the instantaneous kinetic energy and suppresses its fluctuations. The statement in Section 1.5 that 'The resulting dynamics samples phase space consistent with a canonical NVT ensemble' is therefore incorrect. The authors should either implement the stochastic Andersen reassignment or rename the thermostat a Berendsen-type rescaling and qualify the ensemble statement accordingly.

minor comments (5)

[Section 2.4, near Figure 18] The text refers to 'Figure X' when discussing the temperature versus time plot; this should be 'Figure 18'.
[Table 4 caption] The caption reads 'Table 4 – Average pressure (LJ units) comparisons from this report, MC calculations literature MC calculations.' It should be completed, e.g., '... comparisons from this report, MC calculations, and literature values.'
[Section 1.2] The description of the INPUT file contains a garbled sentence: 'The maximum number of neighbors are declared on the next line by declares are the number of timesteps...'. This passage should be rewritten for clarity.
[Section 2.4, Eq. (12)] The definition of r(0) in Eq. (12) as 'the atomic position vector at the beginning of the simulation' is ambiguous given that a 100 ps equilibration precedes the NVE sampling. Clarify that the MSD is evaluated with respect to the initial configuration at the start of the production run, or describe the time-origin convention actually used.
[Section 2.3] The MC simulation details are insufficiently reported: the version of the Allen-Tildesley code, the number of particles, the cutoff, and the averaging protocol are all unspecified. These details are essential to interpret the comparison in Tables 3 and 4.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: all validation targets are external benchmarks, and no fitted parameter is renamed as a prediction.

full rationale

The paper's central claim is that EM3 is a validated modular MD program for Lennard-Jones argon, supported by comparisons against external Monte Carlo results and literature values. The potential-energy and pressure comparisons in Tables 3 and 4 are benchmarked against the Allen-Tildesley MC code and Johnson et al.'s equation of state, which are independent references. The cutoff energy offset is adjusted in the comparisons, but the MD averages are not fitted to any of the reference values; the reported agreement is a direct output of the simulation. The diffusion coefficient is obtained from the Einstein relation applied to the MSD, and the linear regression is standard data analysis rather than a fit to an external target. The cohesive-energy comparison likewise references experimental and ab initio values, and no parameter is tuned to reproduce them. The note about long-range corrections is a correctness or estimation caveat, not a circular step: omitting tail corrections may bias the comparison, but it does not make the derivation depend on the conclusion. The paper does cite prior work by others and one utility code by Tadano, but none of these citations is load-bearing in a circular way, and there is no self-citation chain that defines the target result into existence. Because the validations are against independent external sources and the central quantities are computed from the dynamics without fitting to the benchmark values, no circularity is present.

Assumptions & free parameters 0 free parameters · 4 assumptions · 0 invented entities

The central claim depends on standard domain assumptions: the LJ model for argon, the adequacy of the 3 sigma cutoff without long-range corrections, and the correctness of external reference data. No free parameters are fitted to the target data, and no new entities are introduced.

assumptions (4)

domain assumption Lennard-Jones potential accurately models argon interactions in the simulated state points.
The entire validation assumes the LJ form with epsilon = 1.65e-21 J and sigma = 3.4e-10 m is a valid model for argon. This is standard for argon MD, but it is an external input from literature, not derived in the paper.
domain assumption A cutoff of rc = 3 sigma without long-range corrections is sufficient to reproduce reference values.
Section 2.3 states long-range corrections could be added but were not, assuming 10% agreement is sufficient. The reference MC and Johnson et al. values may be computed with different truncation or corrections.
domain assumption The Allen-Tildesley MC code and Johnson et al. equation of state are correct and directly comparable.
The paper uses these as ground truth for validation without independent verification of their settings or accuracy for the chosen state points.
domain assumption Velocity rescaling every timestep produces proper canonical sampling.
Section 1.5 describes rescaling velocities by sqrt(T0/T) each timestep. This is a simple velocity rescale thermostat, not the standard Andersen stochastic collision scheme, but the paper assumes it samples the NVT ensemble correctly.

how reviews work

0 comments

Cite this review

Pith. "Pith review of Argon Simulations with EM3, a New Modular Molecular Mechanics Program." pith.science (2026). https://pith.science/paper/2GENPYNX

@misc{pith2026190800601,
  author       = {Pith},
  title        = {Pith review of: Argon Simulations with EM3, a New Modular Molecular Mechanics Program},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/2GENPYNX}},
  note         = {Machine review of arXiv:1908.00601}
}

read the original abstract

Argon molecular dynamics (MD) simulations are performed with a newly developed MD program, Easy M(1)odular M(2)olecular M(3)echanics (EM3). The program was developed in an object-oriented fashion containing classes for each critical part of a functioning MD program. An organizational scheme for a general molecular mechanics program is therefore presented, along with the framework of the EM3 program. With the modular nature and open-source availability, EM3 can serve as a learning tool for newcomers to molecular simulations and code organization via object-oriented programming. Validations of the code are presented in comparison with Monte Carlo (MC) simulations of liquid argon at different densities and temperatures. A calculation of the self-diffusion coefficient for liquid argon is also performed, exhibiting the extendibility of EM3. This report comes packaged with the EM3 source code and examples, located at https://github.com/rohskopf/em3.

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

1 extracted references · 1 canonical work pages

[1]

Fast parallel algorithms for short-range molecular dynamics

1 Plimpton, S. Fast parallel algorithms for short-range molecular dynamics. Journal of computational physics 117, 1-19 (1995). 2 Lindahl, E., Hess, B. & Van Der Spoel, D. GROMACS 3.0: a package for molecular simulation and trajectory analysis. Molecular modeling annual 7, 306- 317 (2001). 3 Johnston, M. A., Galván, I. F. & Villà‐Freixa, J. Framework‐based...

work page 1995

Pith tools

Reviewed August 14, 2026 · model on record in the stance chip above.

[1] [1]

Fast parallel algorithms for short-range molecular dynamics

1 Plimpton, S. Fast parallel algorithms for short-range molecular dynamics. Journal of computational physics 117, 1-19 (1995). 2 Lindahl, E., Hess, B. & Van Der Spoel, D. GROMACS 3.0: a package for molecular simulation and trajectory analysis. Molecular modeling annual 7, 306- 317 (2001). 3 Johnston, M. A., Galván, I. F. & Villà‐Freixa, J. Framework‐based...

work page 1995