[lammps-users] large change in volume during simulation of water using npt

Dear lammps developers and users
First of all thanks to Axel Kohimeyer for his suggestion.
I have simulated water at constant temperature and pressure using nve with Berendsen thermostat and berendsen barostat . The system is equilibrating , but the volume of the simulation box changes to 6-7 times the initial volume of the box.Input file is attaching with this. If any correction is required ,please help me.
#Input file for SPC/E water

boundary p p p
restart 50000 restart

#Simulation Units andSimulation Type

units metal
atom_style full

#Reads configuration from file

read_data water500.txt
#read_restart restart.*

#Mass relative mass of atoms

mass 2 15.99943
mass 1 1.00794
group oxygen type 2
group hydrogen type 1
set group oxygen charge -0.8476
set group hydrogen charge 0.4238

#Define potential type

pair_style lj/charmm/coul/long 10.0 12.0 12.0
kspace_style pppm 5.0e-5
kspace_modify mesh 36 36 36 order 5

#Lennard-Jones Parameters

pair_coeff 2 2 0.0067371579 3.166
pair_coeff 1 1 0.0000000000 0.000
pair_coeff 1 2 0.0000000000 0.000
bond_style harmonic
bond_coeff 1 100.0 1.0
angle_style harmonic
angle_coeff 1 300.0 109.47
group water type 1 2
velocity all create 800 234324 dist gaussian mom yes rot yes
neighbor 2.0 bin
neigh_modify delay 0 page 200000 one 4000
#minimize 0 0 1000 10000
thermo_style custom step pe ke etotal temp press

#Applies or “fixes” NVT MD

fix 1 water shake 0.0001 20 500 b 1 a 1
fix a all nve

#Dumps output files and thermodynamic properties

dump 1 all atom 500 dump.waterNPT4
thermo 500

#Timestep is 1 fs (0.001 ps)

timestep 0.001

#Number of time steps (reduced)

#run 20000
compute mobile all temp
fix b all temp/berendsen 373 373 100

fix_modify b temp mobile
#run 20000
compute baro all pressure mobile

fix jsb water press/berendsen iso 1 1 1000

fix jsb water press/berendsen xyz 1 1 2000
fix_modify jsb press baro
thermo_style custom step pe ke etotal temp press vol xlo xhi ylo yhi zlo zhi
run 2000000

Advance thanks to all
Jesseela

If your box size changes with NPT (or press/berendsen) it is because
you have a large positive or negative pressure that is pushing
the box volume to change. In this case you have probably parameteriized
your water incorrectly, or are using a bad damping constant.

Steve

Units metal?
Matt

Quoting Steve Plimpton <[email protected]>:

Dear Developers and users
Thank you all for your suggestions.
Now I have converted the bond coefficients and angle coefficients to metal unit , but no change in volume of the simulation box occurred. But while changing the damping coefficient to high value (1000 to 6000), the volume could be reduced from 7 times to 3 times of the initial volume. My doubt is , what is the criteria in fixing the value of damping coefficient? How much it can be? Kindly help me to clear these doubts. Thanks in advance.
Jesseela

It's a time constant. The bigger you make it, the longer you need
to run to equilibrate the system. But you can't make it arbitrarily
small, or the system volume oscillates too wildly. For temperature,
a damping constant of about 100 timesteps is good, for pressure
about 1000. Those are just rules of thumb.

Steve

Dear Developers and users
Thank you all for your suggestions.
Now I have converted the bond coefficients and angle
coefficients to metal unit , but no change in volume of the simulation box
occurred. But while changing the damping coefficient to high value (1000 to
6000), the volume could be reduced from 7 times to 3 times of the initial
volume. My doubt is , what is the criteria in fixing the value of damping
coefficient? How much it can be? Kindly help me to clear these doubts.

my guess is that this is an indication that your initial coordinates
are very bad. have you tried running for a while at nvt or nve+langevin
first and then switch to npt. the problem with npt is: if your system first
"explodes", it can take a very long time until you reach equilibrium
density.

axel.

Dear developers and users
Thanks to all.
I have tried with nve+langevin for some 500 ps and then with berendsen barostat
for some 15000ps. Now the volume change is almost ok. but the system is not equilibrating to the
required pressure. Shall i get a solution for this?
jesseela

the pressure in MD simulations will fluctuate significantly, especially for small systems. Fluctuations on the order of hundreds of atmospheres are indeed possible. You had better to monitor the density or average pressure over some steps (~1000 steps) that these should be small fluctuation about desired value.

Farrokh

Why is the system not equilibrating to the target pressure? It could
be b/c your damping constant is bad. Or you are not measuring
temp and pressure in the thermo/barostat the same as you
are outputting. Or it is fluctuating and you need to average
to see that the ave is indeed the target.

Steve