equation of state with ReaxFF

Dear all,

I’m testing ReaxFF for fcc Al equation of state (eos). i.e. trying to calculate energy as a function of lattice constant or volume. Strangely, I get only linear increase in energy with increasing the lattice constant, starting from the smaller lattice constant than the equilibrium one. I tried two ways, first just varied lattice constant and calculated energies. Next, used fix_deform with remapping coordinates. Both gave the same results.

As a simple and quick check I tried the same for Si with Si.sw and it gave correct eos. Then I just changed potential to ReaxFF and again got linear increase only. I doubt there is something else that I can’t notice than the problem with the potential. All files attached.

Thanks for any help.


in.Si_bulk (615 Bytes)

al_control (512 Bytes)

Si.sw (653 Bytes)

ffield.reax.Al (21.7 KB)

Aidan or Ray can probably comment on this.


Good idea checking with Si.sw first, so we can eliminate some basic stuff.
If the energy is ioncreasing linearly, then the ReaxFF model you are using
must be predicting a smaller lattice constants than the one you started
with. For sufficiently small lattice constant, you will see the pressure
become positive, due to the ReaxFF repulsives cores. Always with ReaxFF,
you should check whether you are using a good parameterization, something
that is outside the scope of LAMMPS.


Hi Manana,

Although potentials are always intentionally fit to experimental
equilibrium values, they usually end up with some deviations. You
should start with the equilibrium lattice constant that the specific
ReaxFF parameterization you are using is fitted to, instead of the
experimental constant. Also first try a broader sampling range for
the lattice constants.

p.s. Apart from units should be "real" for reax/c, your script looks good.


Thanks all for your replies. Yes, I did try smaller values too. I
started with a<a0 and went to a>a0 with different steps, tried
different ranges. a0 is the lattice constant predicted with ReaxFF.
What makes me shocked is that even for very small a, I mean less than
the lattice constant, energy is just decreasing, no turning point at
all in energy curve (which line in my case).

I also did expansion of the minimized structure with fix deform but
it's the same, pure linear increase.

I used two sets of parameters, one included in LAMMPS and the other
from the author. No help. I think if parametrization is not good it
could only give quantitative inaccuracy but It's very hard for me to
believe that ReaxFF produces eos qualitatively that much wrong,
especially that they show in paper that it fits DFT quite well, for
the same Al.

Best wishes,

P.S. I just checked pressure vs volume and it seems reasonable. Could
there be any problem in computing energy itself?

The pressure should be numerically equal to the limit of deltaU/deltaV, as
deltaV goes to zero. Possibly, what you think is the energy, is not the
energy. Send me a sample script and logfile output.


Hi Aidan,

Here they are, thanks. Log file recorded only last run because I was
using 'jump in file', but there is also Al.txt outputing all the
needed variables.

This is the version where I varied lattice constants and got nice P(V) curve.

Best wishes,

Al.txt (1.12 KB)

al_control (512 Bytes)

in.Al (767 Bytes)

ffield.reax.Al (21.7 KB)

Al.log (1.61 KB)

Thanks for sending this. The nonbond energy was not being tallied, because
you had energy_update_freq=0 in your control file. Normally, this would
have no effect on nonbond energy, but in the Tabulated version of the
nonbond energy, there was an if-test on energy_update_freq. I have removed
this, and it will be in the next release. In the meantime, I suggest
setting energy_update_freq=1.


Thank you very much! Now everything’s fine.