Dear all,
I am using Lammps-1Feb2014 and I am having a problem with dpd/tstat: Including dpd/tstat with hybrid/overlay option changes the E_coul energy (typically to very large positive number for condensed phases). However all the other energy terms remain the same. Here I copied my input files for just 3 water molecules with and without the thermostat.
With dpd/tstat
----------------- Init Section -----------------
units real
atom_style full
pair_style hybrid/overlay lj/charmm/coul/long 9.0 10.0 10.0 dpd/tstat 300 300 9.0 34387
bond_style hybrid harmonic
angle_style hybrid harmonic
kspace_style pppm 0.0001
----------------- Atom Definition Section -----------------
LAMMPS Description
9 atoms
6 bonds
3 angles
0 dihedrals
0 impropers
2 atom types
1 bond types
1 angle types
0.0 41.000 xlo xhi
0.0 41.000 ylo yhi
0.0 101.000 zlo zhi
Masses
1 15.9994
2 1.008
Atoms
1 1 1 -0.8476 5.690 12.751 11.651
2 1 2 0.4238 4.760 12.681 11.281
3 1 2 0.4238 5.800 13.641 12.091
4 2 1 -0.8476 15.551 15.111 7.030
5 2 2 0.4238 14.981 14.951 7.840
6 2 2 0.4238 14.961 15.211 6.230
7 3 1 -0.8476 17.431 6.180 8.560
8 3 2 0.4238 17.761 7.120 8.560
9 3 2 0.4238 17.941 5.640 9.220
Bonds
1 1 1 2
2 1 1 3
3 1 4 5
4 1 4 6
5 1 7 8
6 1 7 9
Angles
1 1 2 1 3
2 1 5 4 6
3 1 8 7 9
----------------- Settings Section -----------------
bond_coeff 1 harmonic 1000.0 1.0
angle_coeff 1 harmonic 1000.0 109.47
pair_coeff * * dpd/tstat 2.4
pair_coeff 1 1 lj/charmm/coul/long 0.1553 3.166
pair_coeff 2 2 lj/charmm/coul/long 0.0 2.058
group spce type 1 2
fix fShakeSPCE spce shake 0.0001 10 100 b 1 a 1
----------------- Run Section -----------------
communicate single vel yes
thermo_style multi
thermo 15
unfix fShakeSPCE
minimize 1.0e-5 1.0e-7 1000 10000
The output for step 0 is like this:
---------------- Step 0 ----- CPU = 0.0000 (sec) ----------------
TotEng = -140.5061 KinEng = 0.0000 Temp = 0.0000
PotEng = -140.5061 E_bond = 0.0789 E_angle = 0.1960
E_dihed = 0.0000 E_impro = 0.0000 E_vdwl = -0.0008
E_coul = -33.7868 E_long = -106.9934 Press = 2.3176
For the case without dpd/tstat, I only change the init and setting sections:
----------------- Init Section -----------------
units real
atom_style full
pair_style hybrid lj/charmm/coul/long 9.0 10.0 10.0
bond_style hybrid harmonic
angle_style hybrid harmonic
kspace_style pppm 0.0001
----------------- Settings Section -----------------
bond_coeff 1 harmonic 1000.0 1.0
angle_coeff 1 harmonic 1000.0 109.47
pair_coeff 1 1 lj/charmm/coul/long 0.1553 3.166
pair_coeff 2 2 lj/charmm/coul/long 0.0 2.058
group spce type 1 2
fix fShakeSPCE spce shake 0.0001 10 100 b 1 a 1
And I get these energy terms:
---------------- Step 0 ----- CPU = 0.0000 (sec) ----------------
TotEng = 0.2434 KinEng = 0.0000 Temp = 0.0000
PotEng = 0.2434 E_bond = 0.0789 E_angle = 0.1960
E_dihed = 0.0000 E_impro = 0.0000 E_vdwl = -0.0008
E_coul = 106.9627 E_long = -106.9934 Press = 1.5867
I don’t know what I am doing wrong here. I appreciate any help on this problem.
Best regards,
Mohammad
Postdoctoral fellow
Department of Chemical and Petroleum Engineering
University of Wyoming
Phone: 307-761-0316