Dear LAMMPS users,
I am doing a testing run on Tip4p water with kspace_style pppm/tip4p at T=298K and P=1bar via fix npt for N=216 water molecules from a initial density of 0.9872 gcm-3 (the orientation of water molecules in the initial configuration are completely artificial, i.e. I have all of them lied up in the same direction, but I assume for liquid water the system should equilibrate itself to natural orientation)
I am using a time step of 2.0fs.
For two different runs, both the one with Pdamp = 50 and the one with Pdamp=1000 has been able to equilibrate the total potential energy to around the same as literature values and with each other (~-10kcal/mol) for the first 400,0000 steps, and the three radial distribution functions also match well with reference radial distribution functions from literature and with each other.
However, strange thing starts to happen after ~400,0000 steps (~10ns), the run with Pdamp = 50 starts to equilibrate to lower potential energy
(~-12kcal/mol) with lower density/large volume and seem to find its new equilibrium there, while the run with Pdamp = 1000 remains the same, with the pressures of both runs are still oscillating around the desired P=1.
Hence I wonder whether the change for the run with Pdamp = 50 is due to the fact that Pdamp is too small? Because from LAMMPS manual, it says that:
“If Pdamp is too small, the pressure and volume can fluctuate wildly; if it is too large, the pressure will take a very long time to equilibrate. A good choice for many models is a Pdamp of around 1000 time steps”
However, in my case the pressure and volumes fluctuation for both runs are around the same, actually the pressure fluctuation for the run with Pdamp = 50 is smaller.
Below is my input if further information on the input setting is needed, otherwise please ignore it:
replicate 1 1 1
group all type 1 2
group O type 1
group H type 2
neighbor 2.0 bin
neigh_modify delay 0 every 1 check yes one 4000
pair_style lj/cut/tip4p/long 1 2 1 1 0.15 8.5 8.5
kspace_style pppm/tip4p 1.0e-5
pair_modify shift yes
pair_modify tail no
pair_coeff 1 1 0.155 3.1536
pair_coeff 1 2 0.0000 0.0000
pair_coeff 2 2 0.0000 0.0000
bond_coeff 1 10000.00 0.9572
angle_coeff 1 10000.0 104.52
velocity all create 298 427 dist gaussian mom yes
restart 1000 old_config_1.dat old_config_2.dat
fix 1 all shake 0.00001 60 0 b 1 a 1
fix 5 all npt temp 298.0 298.0 100 iso 1 1 1000 drag 0.0
thermo_style custom step temp pe evdwl ecoul press vol etail
dump id all xyz 1000 dump.xyz
dump_modify id element O H