# High negative pressure during simulation

Respected peers,
I am trying to simulate the lignin molecule with a mixture of ozone and water with reaxff potential (I used the reaxff potential available on the LAMMPS example.). I am getting some questionable values for pressure. While running for NVE simulation on the following system: lignin + 500 ozone + 3500 water for 1 ns and 0K starting temperature, the system recorded a very high negative pressure of negative tens of thousand. Later this value saturated to around -2000 as the simulation reaches the end of 1 ns (the temperature at the end of 1 ns is around 473 K). According to LAMMPS documentation, the pressure is supposed to be atm units (I used ‘real’ units). Considering the negative pressure, I incremented water molecules to 5000 but the pressure trend remained the same (I simulated this particular system for 10 ns). I am pretty new with LAMMPS and I will be grateful if someone will help me understand these results
Below is my input file (for lignin + 500 ozone + 3500 water system).

units real
boundary p p p

atom_style charge

pair_style reaxff lmp_control
pair_coeff * * ffield.reax.cho C H O

group lig id 1:201
group ozn id 202:1701
group water id >= 1702
group lig+ozn id 1:1701

compute reax all pair reaxff

#compute pressure all pressure thermo_temp

variable eb equal c_reax[1]
variable ea equal c_reax[2]
variable elp equal c_reax[3]
variable emol equal c_reax[4]
variable ev equal c_reax[5]
variable epen equal c_reax[6]
variable ecoa equal c_reax[7]
variable ehb equal c_reax[8]
variable et equal c_reax[9]
variable eco equal c_reax[10]
variable ew equal c_reax[11]
variable ep equal c_reax[12]
variable efi equal c_reax[13]
variable eqeq equal c_reax[14]

thermo 1

neighbor 2.5 bin
neigh_modify every 10 delay 0 check yes

#velocity all create 300.0 4928459 mom yes rot yes dist gaussian

#fix 1 all nve temp 300.0 300.0 0.1
fix 1 all nve
fix 2 all qeq/reaxff 1 0.0 10.0 1.0e-6 reaxff
fix 4 all reaxff/bonds 20 bonds.reaxff.lignin+ozone+water_nve300_1ns
fix 5 lig reaxff/bonds 20 bonds.reaxff.lignin+ozone+water(lignin)_nve300_1ns
fix 6 ozn reaxff/bonds 20 bonds.reaxff.lignin+ozone+water(ozone)_nve300_1ns
fix 7 lig+ozn reaxff/bonds 20 bonds.reaxff.lignin+ozone+water(lignin+ozone)_nve300_1ns
#fix 8 all press/berendsen iso 0.0 0.0 1000.0
#fix 5 all bond/break 5 1 1.9

minimize 1.0e-6 1.0e-9 10000 100000

log log.lignin+ozone+water.nve300_1ns

variable nqeq equal f_2

thermo_style custom step temp epair etotal press &
v_eb v_ea v_elp v_emol v_ev v_epen v_ecoa &
v_ehb v_et v_eco v_ew v_ep v_efi v_eqeq v_nqeq

timestep 0.1

dump 1 all atom 20 dump.reaxff.lignin+ozone+water_nve300_1ns
dump 2 lig atom 20 dump.reaxff.lignin+ozone+water(lig)_nve300_1ns
dump 3 ozn atom 20 dump.reaxff.lignin+ozone+water(ozn)_nve300_1ns
dump 4 lig+ozn atom 20 dump.reaxff.lignin+ozone+water(lig+ozn)_nve300_1ns

restart 10000 restart.lignin+ozone+water_nve300_1ns

run 10000

Please let me know if the above-mentioned information is enough to pinpoint the issue. My apologies as being a new user, I couldn’t upload my log file.

a) you are using a combustion parameter set for a molecular system. That is a very bad idea.
b) a negative pressure is an indication that your system “wants” to shrink. Does your system have a proper density? You may need to shrink the box. But also you need to equilibrate your system. Also, condensed systems are not very compressible; so small changes in density will have large pressure changes.

If you are also new to MD, then it is not a good idea to start with such a complex system and an advanced force field as ReaxFF. You are likely to make many simple “procedural” mistakes that you can avoid by practicing your MD skils with reproducing results from simpler systems with simpler force fields.

Hi,