Skip to main content
U.S. flag

An official website of the United States government

Official websites use .gov
A .gov website belongs to an official government organization in the United States.

Secure .gov websites use HTTPS
A lock ( ) or https:// means you’ve safely connected to the .gov website. Share sensitive information only on official, secure websites.

LAMMPS MD: Equation of State (pressure vs. density) - Linear-Force Shifted Potential at 2.5σ

The main purpose of the following data set is to present equation of state (density-pressure-temperature) data for a version of the Lennard-Jones fluid that was obtained using the LAMMPS Molecular Dynamics (MD) simulation suite. The secondary purpose of this data set is to provide sample LAMMPS input and initial configuration files that an end user may use in LAMMPS to obtain the same equation of state data, to either verify output from an installation of LAMMPS or educate the user about the basic features of LAMMPS. We first present output from LAMMPS, and compare it to Monte-Carlo derived results, and then provide a tutorial on how to reproduce that same data using LAMMPS.

Simulation details

Simulation Package LAMMPS - see installation tutorial
Method Canonical Ensemble (fixed N, V, T) Molecular Dynamics using default LAMMPS "nvt" ensemble option
Number of LJ particles 10000
Simulation Cell Cubic cell, volume set with constant N=10000 to achieve desired reduced density
Truncation Linear-force Shifted at 2.5σ  (σ = Lennard-Jones length)
Time Step 0.005τ   (τ = Lennard-Jones dimensionless time unit)
Thermostat Temperature imposed via Nosé-Hoover chained thermostats (three, LAMMPS default) with Tdamp=0.500τ
Simulation length 1.1x106 time steps (intended as 105 equilibration steps followed by 106 production steps)

Additional simulation details are in the sample input file, see below.

Results

1. Results for select temperatures (tabular)

The following table contains pressure and intermolecular potential energy per LJ particle (both appropriately reduced) as a function of density for three reduced temperatures, kT/ε=0.75, 1.00, and 1.50. The ensemble averages were calculated the LAMMPS runs described above and uncertainty, stated by the standard error, were estimated using block analysis (see below). For kT/ε=0.75, which is below the critical temperature kTc/ε=0.931, the unstable state densities are clearly identified.
  kT/ε = 0.75 kT/ε = 1.00 kT/ε = 1.50
ρσ3 3 +/- U/ε +/- 3 +/- U/ε +/- 3 +/- U/ε +/-
0.001 7.457E-04 9.501E-08 -8.279E-03 1.260E-04 9.967E-04 4.782E-08 -6.704E-03 -7.084E-05 1.498E-03 1.017E-07 -5.333E-03 3.559E-05
0.005 3.641E-03 3.170E-07 -4.200E-02 1.346E-04 4.919E-03 4.316E-07 -3.362E-02 9.019E-05 7.454E-03 4.140E-07 -2.706E-02 6.825E-05
0.010 7.068E-03 8.143E-07 -8.455E-02 1.355E-04 9.677E-03 5.384E-07 -6.727E-02 1.125E-04 1.482E-02 8.775E-07 -5.394E-02 8.261E-05
0.050 2.672E-02 1.519E-05 -4.535E-01 6.865E-04 4.230E-02 4.402E-06 -3.344E-01 1.353E-04 7.078E-02 6.713E-06 -2.664E-01 1.032E-04
0.100 unstable 7.113E-02 1.505E-05 -6.599E-01 3.522E-04 1.349E-01 1.551E-05 -5.241E-01 1.260E-04
0.200 unstable 1.015E-01 3.248E-05 -1.261E+00 5.928E-04 2.542E-01 2.655E-05 -1.015E+00 1.845E-04
0.300 unstable 1.156E-01 4.205E-05 -1.768E+00 6.053E-04 3.834E-01 6.212E-05 -1.484E+00 1.531E-04
0.400 unstable 1.314E-01 1.049E-04 -2.207E+00 2.404E-04 5.613E-01 1.139E-04 -1.945E+00 1.057E-04
0.500 unstable 1.801E-01 1.348E-04 -2.643E+00 9.893E-05 8.608E-01 1.970E-04 -2.410E+00 6.261E-05
0.600 unstable 3.649E-01 1.580E-04 -3.121E+00 6.989E-05 1.422E+00 1.237E-04 -2.877E+00 3.458E-05
0.700 5.198E-02 2.775E-04 -3.801E+00 7.699E-05 9.132E-01 2.570E-04 -3.624E+00 6.092E-05 2.495E+00 4.709E-04 -3.316E+00 1.007E-04
0.800 9.281E-01 4.624E-04 -4.319E+00 7.030E-05 2.206E+00 3.675E-04 -4.087E+00 6.115E-05 4.466E+00 5.145E-04 -3.674E+00 9.027E-05
0.900 3.010E+00 2.906E-04 -4.732E+00 4.946E-05 4.782E+00 5.599E-04 -4.423E+00 9.073E-05 7.886E+00 7.175E-04 -3.876E+00 1.352E-04
1.000 4.254E+00 2.191E-01 -5.440E+00 3.997E-02 9.410E+00 7.370E-04 -4.524E+00 1.291E-04 1.350E+01 1.483E-03 -3.829E+00 2.519E-04

 

2. Complete Pressure-Density Phase Diagram

The following graphics show pressure-density equations of state for select temperatures (0.65 < kT/ε < 1.20) in the form of a phase diagram. In all four graphics, the dashed line indicates the vapor-liquid coexistence boundary as calculated from Grand Canonical-Transition Matrix Monte Carlo Simulations (GC-TMMC, see results elsewhere in the NIST SRSW). Solid symbols indicate pressure-density-temperature results from LAMMPS simulations, with the standard error shown by error bars. Solid lines are equation of state data also from GC-TMMC. There is outstanding agreement between the LAMMPS MD results and GC-TMMC results, except for the isotherm just above the critical temperature (i.e., the red line and points, at kT/ε= 0.95) where system-size effects are expected to lead to disagreement between different techniques. (The LAMMPS results used N=10000 whereas the GC-TMMC simulations used 0 < N < 435.)

Lennard-Jones Fluid (LFS 2.5$$\sigma$$) Phase Diagram

Lennard-Jones Fluid (LFS 2.5$$\sigma$$) Phase Diagram - Vapor Side

Lennard-Jones Fluid (LFS 2.5$$\sigma$$) Phase Diagram

Lennard-Jones Fluid (LFS 2.5$$\sigma$$) Phase Diagram - Liquid Side

Reproducing output

The LAMMPS MD results in the preceding table and graphics may be reproduced using example LAMMPS runs described as follows.

1. LAMMPS Installation

The data shown above were generated using a generic installation of LAMMPS, the executables of which may be reproduced as described here. This tutorial assumes some level of familiarity with POSIX-compliant operating systems (e.g., Unix, Linux, or Mac OSX) at the command-prompt level and access to a system with the GCC compiler, OpenMPI parallelization suite, Python, and the git version control system. This tutorial uses ">" to indicate the shell command prompt and $LAMMPS_DIR to identify the directory where LAMMPS is downloaded and later compiled. First, LAMMPS can be obtained using the instructions at http://lammps.sandia.gov/download.html#git via:

   >git clone https://github.com/lammps/lammps.git

Second, LAMMPS executables may be compiled via:

   >cd $LAMMPS_DIR/src
   >git checkout r15061
   >make purge
   >make package-update
   >make yes-user-misc
   >make -j8 mpi

The installation sequence 1) switches to the "r15061" commit of LAMMPS (to ensure that a user is using the same code used to generate the data shown above), 2) removes any existing installation of LAMMPS, 3) updates any out of date packages, 4) installs the "USER-MISC" package to enable the force-shift version of Lennard-Jones potential ("lj/fs" in LAMMPS terminology), and 5) builds an MPI-enabled executable using eight processor cores for parallel compilation. The end result is an executable named "lmp_mpi" located in the $LAMMPS_DIR/src directory.

Finally, we used the following operating system, compiler, and MPI libraries to build the LAMMPS executable:
   Linux OS: CentOS 7; 3.10.0-327.13.1.el7.x86_64 kernel
   GCC: v4.8.5 (2015-06-23, listed as Red Hat 4.8.5-4)
   OpenMPI: v1.10.0 (2015-08-24)
   LAMMPS: Git Checkout Tag r15061 (2016-05-14)

2. Example LAMMPS Scripts

Example LAMMPS scripts and initial configurations that will yield the data shown above may be obtained from a git repository at: https://github.com/dwsideriusNIST/LAMMPS_Examples, e.g. via:

   >git clone git [at] github.com:dwsideriusNIST/LAMMPS_Examples.git

This tutorial assumes that the git repository resides in $EXAMPLES_DIR. The relevant inner directories are:

LJ_initial_cfgs

Initial configurations of N=10000 Lennard-Jones particles at the specified density
i.e., in.nvt.dens_0.2000 is an initial configuration at density ρσ3=0.2000.

run_scripts Relevant LAMMPS Scripts
analysis Python scripts to post-process LAMMPS output
LJ_example Example shell script that will run a LAMMPS simulation at ρσ3=0.4000

3. Example LAMMPS Run

A short example LAMMPS run is available in $EXAMPLES_DIR/LJ_example, which may be run to confirm that LAMMPS is operating properly. This simulation is executed via:

   >cd $EXAMPLES_DIR/LJ_example
   >sh example.sh

In the "example.sh" script, one can see that this script runs a LJ simulation at ρσ3=0.4000 and kT/ε=1.50 for 200000 time steps (all LAMMPS settings are as shown previously) then analyzes the MD trajectory to report the average temperature, pressure, and potential energy. The final output, from the block analysis script, will report the following ensemble averages and uncertainty:

Thermodynamic Ensemble Averages from: ave.dens_0.4000.out
   Discarded Timesteps: 100000
   Number of Blocks: 5
   Block size (timesteps): 20000
   Stated uncertainty is the standard error of the thermodynamic property


Temperature      1.50006 +/- 0.00018
PotentialEnergy -1.94506 +/- 0.00030
Pressure         0.56070 +/- 0.00024

We note that the average pressure and energy differ from those in the table above because this example run was only 200000 timesteps.

4. Production LAMMPS Runs

Production LAMMPS runs to generate LJ equation of state data may be run using the scripts provided in $EXAMPLES_DIR/run_scripts and the initial configurations in $EXAMPLES_DIR/LJ_initial_cfgs. Syntax for execution of a run is:
 
   >mpirun -np $NP lmp_mpi -in LJ.NVT.startfromrestart -var rho $FOO -var temp $BAR
 
where $NP is the number of processors to utilize, $FOO is the reduced density, and $BAR is the reduced temperature. We note that the initial configuration file (e.g., in.nvt.dens_0.4000) must be in the run directory. Block Analysis of the results may be done via the script provided in $EXAMPLES_DIR/analysis:
 
   >$EXAMPLES_DIR/analysis/block_analysis.py -f $FILENAME -b $BLOCKS -m $STEPS_SKIP
 
where $FILENAME is the LAMMPS output file containing MD trajectory data, $BLOCKS is the number of blocks to use for uncertainty analysis, and $STEPS_SKIP is the number of timesteps from the start of the simulation to discard (i.e., "equilibration" timesteps).
 
Lastly, an example production run at kT/ε=0.75, ρσ3=0.8000, using 12 processors may done via:
 
   >cd $EXAMPLES_DIR
   >mkdir test ; cd test
   >cp ../run_scripts/LJ.NVT.startfromrestart ./ ; cp ../LJ_initial_cfgs/in.nvt.dens_0.8000 ./
   >mpirun -np 12 lmp_mpi -in LJ.NVT.startfromrestart -var rho 0.8000 -var temp 0.75
 
Based on the input file, this simulation will run for 1100000 time steps, and produce a MD trajectory file named ave.dens_0.8000.out. The trajectory may be analyzed (discarding the first 100000 timesteps and using 10 blocks to estimate the uncertainty in the remaining 1000000 timesteps) via:
 
   >../analysis/block_analysis.py -f ave.dens_0.8000.out -b 10 -m 100000
 
The final output is:
 

Thermodynamic Ensemble Averages from: ave.dens_0.8000.out
  Discarded Timesteps:    100000
  Number of Blocks:       10
  Block size (timesteps): 100000
  Stated uncertainty is the standard error of the thermodynamic property

Temperature      0.74998  +/-   0.00007
InternalEnergy  -4.31895  +/-   0.00007
Pressure         0.92811  +/-   0.00046

 
The results shown are identical to those for kT/ε=0.75 and ρσ3=0.8000 in the previous table. Finally, we note that these results will be generated exactly as shown provided that the number of parallel processors is set to 12. Choosing a different number of processors will affect the resultant ensemble average and uncertainty due to differences in the random number generator in LAMMPS as well as round-off truncations.
Created August 15, 2016, Updated June 2, 2021