Hi,
I'm trying to simulate Tip3p-Ew rigid water model using DPD thermostat in NVE ensemble. I'm using shake algorithm to constraint bonds and angles. I'm getting a very high pressure (10^14) even at 1000 steps and temperature is not controlling. The temperature increases with time. The system is supposed to be at 300 K and 1 atm. I did tests changing the DPD cutoff and the damping constant. but any of them not worked. I used the following LAMMPS input script. It would be really appreciated if anyone tell me if I did anything wrong in my LAMMPS input and any suggestions to resolve this problem.
LAMMPS input script
log lammps.shear
units real
newton on
neigh_modify delay 5 every 1 one 2500
comm_modify vel yes
boundary p p p
atom_style full
bond_style harmonic
angle_style harmonic
dihedral_style opls
special_bonds lj 0.0 0.0 0.5
pair_style hybrid/overlay dpd/tstat 1.0 1.0 12.0 34387 lj/cut/coul/long/kk 12.0
kspace_style pppm 1.0e-4
kspace_modify mesh 30 30 30
read_data nvt_3ns_tip3p_correct_density_deform_28_water_512_300K_0.4ns.data
kspace_style pppm 1.0e-4
group water type 1:2
pair_coeff 1 1 lj/cut/coul/long/kk 0.102 3.188 12.0
pair_coeff 1 2 lj/cut/coul/long/kk 0.051 1.594 12.0
pair_coeff 2 2 lj/cut/coul/long/kk 0.0 0.5 12.0
pair_coeff 1 1 dpd/tstat 0.01 12.0 #dpd cutoff
pair_coeff 1 2 dpd/tstat 0.01 12.0
pair_coeff 2 2 dpd/tstat 0.01 12.0
#bond coeff
bond_coeff 1 119 0.9572
#angle coeff
angle_coeff 1 101 104.52
fix 1 all nve
fix 2 water shake/kk 0.001 20 0 b 1 a 1
thermo 1000
thermo_style custom step temp press vol density cpu
timestep 1.0
restart 1000 restart.A restart.1
dump 2 all custom 10000 water_atom.1 id type xs ys zs ix iy iz mol
run 10000000
exit
Thank you in advance!