correct volume for supercritical CO2 system

Dear Lammps-Users,

I am trying to simulate 1000 CO2 molecule to get the correct density
at super critical condition.
I am using ff reax for hydrocarbon of chenoweth.
These are the steps I ran:
1. minimize the initial system.
2. equilibrate the system with NVT for 10ps
3. raise the temperature to critical temperature (304K) with NVT for 20ps.
4. equilibrate the system at Tcrit with NVT for 20ps
5. Do NPT for figuring out the density at critical condition (T = 304,
P = 72 atm) until the volume is stable.

However, the volume I got is very large roughly about 1,000,000 A^3
which give the density about 0.073 g/cc.
The value of density from NIST of CO2 at supercritical condition is
about 0.4674 g/cc

I am doing this to check if I could get a correct density for
supercritical CO2 at higher pressure and temperature.
I attached my input and the last log file when I did my NPT.

Erik (148 KB)

in.min (2.83 KB)

log.npt.tpcrit.5 (217 KB)

ffield.reax (7.02 KB)

Two suggestions:

a) check if all your runs are long enough that the system came to a stable
equilibrium. Particularly the NPT. Since your system is small - see if
the volume and pressure are stable.

b) Try setting the volume at 0.46 g/cc^3 and running NVT
and monitoring the pressure, to see what its stable value
is. Compare that to your value at the much lower density.


Steps 1-4 do little for you because your initial minimization does not
include changes to the box volume.
If your initial volume is too far the actual critical value, all you
have achieved from 1-4 is to get the CO2 molecules
bouncing around.


In addition, running at or around the critical point with npt is in general a very bad idea. The critical point is accompanied by large density fluctuations. This cannot be represented by a 1000 molecule system. In addition, there is very little to gain by running with npt over nvt for such a system. I would just set the density and at best run a larger system, if you can afford it and run at a state with sufficient distance to the critical point.


Hi All,

Thank you for the suggestion. I only run this to check of I could get
a correct density at given T and P of supercritical fluid.
The real simulation I want to run is that we know the T and P but we
not sure the density because its a mixture of some molecules.
Is there any suggestion how I could get a correct density if we only
know the T and P of interest ?


Hi Mr. Plimpton,

I did run several times before I sent the email to the mailist.
a. I run enough time to let the volume and pressure stabilize and only
fluctuate a little.
b. That box is initially already set as the correct density, but it
become a lot bigger when I run NPT.

I also did run a NPT at atm condition (300K, 1 atm). The density that
I got when I ran this condition is almost correct to the theoretical
I am guessing that its not working for supercritical fluid value.

The problem is, I would like to run a mixture at a High Pressure which
CO2 will be in supercritical phase.
Because of lack of theoretical value at certain P and T, I couldnt get
the value of density of the mixture.
Thats why I thought I could run this to correct the density of my system.


Why using reaxff for a phase coexistence problem where reactions don't
matter? If there is no further justification to employing reaxff,
ditch it, as it just going to make slower any progress to be made.
As to your question you could:

1-) Look at this paper
J. G. Harris; K. H. Yung; "Carbon Dioxide's Liquid-Vapor Coexistence
Curve and Critical Properties As Predicted by a Simple Molecular
Model", J. Phys. Chem. 99 12021-12024 (1995).
2-) Use configurational-bias-monte-carlo in addition to/instead of MD.
The paper above is actually one of the references within the Towhee
code freely available on the web. There might even be examples within
the Towhee distribution on how do accomplish this task. Yes, more to
be learned but you could use equilibrated configs from Towhee as input
for Lammps (will work best if the forcefield to be used is the same in
both codes)

Working with the two methods in parallel will certainly help you
validate the consistence of your results.
That's how I would address your problem.


Dear Mr. Campana,

I will look at the paper ant that method.
I am using ReaxFF because I was trying to do reaction latter with it.
My plan before is to get the system to correct density first before I
started the reaction process.
So you think I won't be able to do this in lammps ?


Dear Mr. Campana,

I will look at the paper ant that method.
I am using ReaxFF because I was trying to do reaction latter with it.
My plan before is to get the system to correct density first before I
started the reaction process.
So you think I won’t be able to do this in lammps ?

You certainly can use lammps to do it. But for just the phase coexistence CBMC is another way to reach the same goal.
I just wanted to point out alternative ways for you to cross validate your results.

Again, I would try running NVT at 2 densities (volumes).
The expt correct one, and the larger volume that you are getting.
Monitor the pressure in both cases. If they are not that different
and the pressure is fluctuating a lot (which it will for a small system),
then you are asking NPT to do something difficult. You might also
play with your Pdamp constant when you run NPT.