Dear professor,
i want to obtain some transform properties of CO2 use the classical EPM2 model,
a model with rigid bond lengths but with a flexible bond angle potential(harmoonic)
i always cannot obtain the right density and i have read some advice about fix shake(which is not fit to linear molecule but i only shake the bond of CO2 not angle).could you give me some suggestion,thank you !
here are my script.
variable t equal 423
#setup
units real
dimension 3
atom_style full
neighbor 1.7 bin
neigh_modify delay 0 every 1 check yes
read_data data.1000co2
mass 1 12.0107 # C
mass 2 15.9994 # O
group Cgroup type 1
group co2 type 1 2
pair_style lj/cut/coul/long 12.0 8.5
pair_modify tail yes
kspace_style pppm 1.0e-5
pair_coeff 1 1 0.05587 2.757
pair_coeff 1 2 0.0945 2.8917
pair_coeff 2 2 0.1599 3.033
bond_style harmonic
bond_coeff 1 0 1.149
angle_style harmonic
angle_coeff 1 14.77 180
velocity all create $t 22222222 dist gaussian
fix 1 all shake 0.000001 100 0 b 1
fix 3 all nve
fix 4 all temp/berendsen $t t (100*dt)
thermo_style custom step temp ke pe etotal press vol
timestep 1
thermo 10000
run 1000000
unfix 3
unfix 4
reset_timestep 0
fix 5 all npt temp $t t (100dt) iso 246.73 246.73 $(1000dt)
thermo_style custom step temp ke pe etotal press vol
timestep 1
thermo 10000
run 1000000
unfix 5
reset_timestep 0
fix 6 all npt temp $t t (100dt) iso 246.73 246.73 $(1000dt)
thermo_style custom step temp ke pe etotal press density vol
dump 3 Cgroup custom 50 CARBON.lammpstrj id type vx vy vz
dump_modify 3 sort id
timestep 1
thermo 20000
run 3000000