Potential energy per water molecule

Dear LAMMPS users ,

I do NVT simulation on water and I use " compute pe " for calculating potential .
In my system I have also one fix SHAKE and three fix wall but As I don’t use pair_modify energy yes so the potential energy of these fixes will not contribute .
The potential energy per water molecules that I calculate ( at around 10K ) is around -79 kcal/mol while Average Configurational energy per water molecules is -46.5 kj/mol .

I use :
units real ,
TIP4P/2005 water model ,
Three fix wall/lj126 ,
time_step 0.1 .

Can you please let me know what is the reason of this problem ?
Thanks for your attention .

Regards,
Saeed.

there is not enough info here to know what

your are computing, or why you think one number is wrong and

one is right.

Steve

Dear Steve ,
Thanks for your reply .
I have a system of water molecules and 3 lj/126 wall that effect on water’s oxygen .
When I do NVT simulation on these molecules I get potential energy around -79 kcal/mol per water molecule .
Average Configurational energy per water molecules is -46.5 kj/mol and I believe we can’t get potential energy for water molecule lower than this amount !

I use compute pe command for potential energy of entire system .

Model

Dipole moment e

Dielectric constant

self-diffusion, 10-5cm2/s

Average configurational energy, kJ mol-1

Density maximum,
°C

Expansion coefficient,
10-4 °C-1

SSD

2.35 [511]

72 [511]

2.13 [511]

-40.2 [511] -13 [511] -

SPC

2.27 [181]

65 [185]

3.85 [182]

-41.0 [185] -45 [983[]](http://www1.lsbu.ac.uk/water/ref10.html#r983) 7.3 [704] **

SPC/E

2.35 [3]

71 [3]

2.49 [182]

-41.5 [3] -38 [183] 5.14 [994]

SPC/Fw

2.39 [994]

79.63 [994]

2.32 [994]

PPC

2.52 [3]

77 [3]

2.6 [3]

-43.2 [3] +4 [184] -

TIP3P

2.35 [180]

82 [3]

5.19 [182] -41.1 [180] -91 [983[]](http://www1.lsbu.ac.uk/water/ref2.html#r180) 9.2 [180]

TIP3P/Fw

2.57 [994]

193 [994]

3.53 [994] - - 7.81 [994]

iAMOEBA

2.78 [2031]

80.7 [2031]

2.54 [2031] - 4 [2031] 2.5 [2031]

uAMOEBA

2.80 [2401]

76.3 [2401]

2.41 [2401] - - 3.38 [2401]

QCT **
1.85 [1251] -
1.5 [1251]

-42.7 [1251]
+10 [1251] 3.5 [1251]

TIP4P

2.18 [3,[180]](http://www1.lsbu.ac.uk/water/ref2.html#r180)

53 a [3]

3.29 [182] -41.8 [180] -25 [180] 4.4 [180]

TIP4P-Ew
2.32 [649] 62.9 [649] 2.4 [649] -46.5 [649] +1[649] 3.1[649]

TIP4P-FQ

2.64[197]

79 [197]

1.93 [197] -41.4 [201] +7 [197] -

TIP4P/2005

2.305 [984]

60 [984]

2.08 [984] - +5 [984] 2.8 [984]

TIP4P/2005f

2.319 [1765]

55.3 [1765]

1.93 [1765 ] - +7 [1765 ] -

TIP4P/ε

2.4345 [2444]

78.3 [2444]

2.10 [2444] - +4 [2444] -

OPC

2.48 [2168]

78.4 [2168]

2.3 [2168] - -1 [2168] 2.7 [2168]

OPC3

2.43 [2722]

78.4 [2722]

2.3 [2722] - -13 [2722] 4.3 [2722]

SWFLEX-AI

2.69 [201]

116 [201] 3.66 [201] -41.7 [201] - -

COS/G3 **

2.57 [704]

88 [704]

2.6 [704] -41.1 [704] -78[1939] 7.0 [704]

COS/D2

2.55[ 1617[]](http://www1.lsbu.ac.uk/water/ref8.html#r704)

78.9 [1617[]](http://www1.lsbu.ac.uk/water/ref8.html#r704)

2.2 [1617[]](http://www1.lsbu.ac.uk/water/ref8.html#r704) -41.8 [1617[]](http://www1.lsbu.ac.uk/water/ref8.html#r704) - 4.9 [1617[]](http://www1.lsbu.ac.uk/water/ref8.html#r704)

GCPM

2.723 [859[]](http://www1.lsbu.ac.uk/water/ref9.html#r859)

84.3 [859[]](http://www1.lsbu.ac.uk/water/ref9.html#r859)

2.26 [859[]](http://www1.lsbu.ac.uk/water/ref9.html#r859) -44.8 [859] -13 [859[]](http://www1.lsbu.ac.uk/water/ref9.html#r859) -

SWM4-NDP

2.461 [933[]](http://www1.lsbu.ac.uk/water/ref9.html#r859)

79 [933[]](http://www1.lsbu.ac.uk/water/ref9.html#r859)

2.33 [933[]](http://www1.lsbu.ac.uk/water/ref9.html#r859) -41.5 [933[]](http://www1.lsbu.ac.uk/water/ref9.html#r859)
<-53 [1999]

Is average config Eng in the table, the same as what compute pe is calculating?

That table is for bulk water. Do you get similar values with LAMMPS if you

compute buik water - no walls? Are you at the same density as in the table?

Steve

There are a lot of problems with your question, but the main one is: you can not compare results at 10 K with results at 25 C. You should take the time to read some of the papers cited in the table that you posted, instead of looking for a quick fix on this mailing list.

Hello,

I have the same kind of question concerning the energy per H2O in liquid at normal conditions of temperature and pressure. When using the TIP4P/2005 potential, i get -9.6 kcal/mol/H2O, which is good when compared to exp. value in the same conditions, of -9.92 kcal/mol/H2O. But when using REAX/C with ffield.reax.ZnOH (which seems a good choice because this ffield belongs to the “water branch” family of REAX parameter sets), i get the stunning value of -255 kcal/mol/H2O (final density in NPT is 0.9 g/cm3). I checked my units are “real”. This value is given by etotal/(N_H2O).
One can have access to 14 contributions of this energy, using compute pair reax/c. I may be wrong, but only some of these contributions may be useful for the calculation of Energy/H2O with Reax/c? If so, could someone tell me which contributions should I add up?
I think my input script is right. It works well with TIP4P/2005. Both are in real units. The only difference with TIP4P/2005 is the potential part. No shake with Reax/c, and :
timestep 0.2
pair_style reax/c NULL
pair_coeff * * ffield.reax.ZnOH O H
fix myQeqReax all qeq/reax 1 0.0 10.0 1.0e-6 reax/c

ffield.reax.ZnOH is available in lammps/example
Here type 1 is Oxygen and type 2 is Hydrogen.

Regards,
Xavier

There are a variety of ReaxFF containing O and H. Try the ones being used for aquated and battery simulations.

Jim

Hello,

I have the same kind of question concerning the energy per H2O in liquid at normal conditions of temperature and pressure. When using the TIP4P/2005 potential, i get -9.6 kcal/mol/H2O, which is good when compared to exp. value in the same conditions, of -9.92 kcal/mol/H2O. But when using REAX/C with ffield.reax.ZnOH (which seems a good choice because this ffield belongs to the "water branch" family of REAX parameter sets), i get the stunning value of -255 kcal/mol/H2O (final density in NPT is 0.9 g/cm3). I checked my units are "real". This value is given by etotal/(N_H2O).
One can have access to 14 contributions of this energy, using compute pair reax/c. I may be wrong, but only some of these contributions may be useful for the calculation of Energy/H2O with Reax/c? If so, could someone tell me which contributions should I add up?

xavier,

absolute potential energy values have no meaning, only energy
differences. thus the *only* somewhat meaningful value for
etotal_bulk/n_h2o would happen when the value of etotal for an
isolated water is 0kcal/mol. this should indeed be the case, when you
use fix shake and fixed boundary conditions. this should already break
down for flexible water potentials.

the reaxff formulation is very complex with many inter- and
intra-molecular interactions.

axel.

Indeed, when substracting the energy of an isolated "Reax" water molecule, this yields the right order of magnitude, with -9.1 kcal/mol/H2O.
Thank you Axel!

Xavier

-----Message d'origine-----