themerature distribution in spce water

Dear Lammps users

I am trying to simulate a 3d water spce flow between two graphene sheets looking to the temperature distribution in the water domain. I set the temperature of the walls via nvt thermostat and used nve ensemble for water molecules. to extract the temperature distribution in the domain I used the fix ave/spacial command and computed the temperature from T=2/3*Ke/kb, where Ke is the per atom kinetic energy and kb is boltzman constant. the obtained temperature distribution is kind of weird. the temperature of the fluid is in the order of 3000 k!!
It would be kind of you to help me through that. the temp distribution between two walls is:

144.777
3002.74
3078.13
3166.06
3244.73
3291.92
3308.02
3330.33
3314
3276.31
3223.01
3122.99
3101.85
114.177

which the first and the last are the wall temperatures.
I’m using “lammps-22Feb13” version. and my input script is as below:

----------------- Init Section -----------------

dimension 3
boundary p p f
units real
atom_style full
pair_style hybrid airebo 2.5 0 1 lj/charmm/coul/long 9.0 10.0 10.0 lj/cut 8.5
bond_style hybrid harmonic
angle_style hybrid harmonic
kspace_style pppm 0.0001
kspace_modify slab nozforce #overlap no
pair_modify mix arithmetic
neighbor 3.0 bin
neigh_modify check yes

----------------- Atom Definition Section -----------------

read_data “system.data”
group topwallatoms type 1
group lowwallatoms type 2
group wallatoms union topwallatoms lowwallatoms
group spce type 3 4

----------------- Variables Section -----------------

variable kb equal 1.38066e-23
variable An equal 6.0221415e23
variable convertor equal 4184/${An}

----------------- Settings Section -----------------

pair_coeff * * airebo CH.airebo C C NULL NULL
pair_coeff 1 2 none
pair_coeff 3 3 lj/charmm/coul/long 0.1553 3.166
pair_coeff 4 4 lj/charmm/coul/long 0.0 2.058
pair_coeff 12 3 lj/cut 0.196844494 3.283 8.5 #LJ using Mixing rules
pair_coeff 1
2 4 lj/cut 0 3.4 1.0 #no interaction between Ar and #Hydrojen
bond_coeff 1 harmonic 1000.0 1.0
angle_coeff 1 harmonic 1000.0 109.47

----------------- Run Section -----------------

minimize 1.0e-5 1.0e-7 1000 10000
fix fShakeSPCE spce shake 0.0001 10 100 b 1 a 1
velocity spce create 115 1234 dist gaussian mom yes rot yes
fix 11 wallatoms setforce NULL NULL 0.0

timestep 0.1 #femtosecond (e-15)

thermo 1000
thermo_style custom step etotal vol
dump 1 all custom 1000 traj_npt.lammpstrj id mol type x y z ix iy iz

fix 5 spce nve
fix 3 topwallatoms nvt temp 140.0 140.0 10
fix 4 lowwallatoms nvt temp 100.0 100.0 10

run 1000000

write_restart restart2.dat
compute KErn all ke/atom
variable tempe atom (2/3)(c_KErn)(v_convertor)/v_kb
variable vxx atom vx*100000
fix valave all ave/spatial 1 10 10 x 0 1 y 0 1 z 0 0.05 v_vxx v_tempe units reduced norm sample ave running file Vel20.dat overwrite # ---- (end of examples) ----

run 1000000

Dear Lammps users

I am trying to simulate a 3d water spce flow between two graphene sheets looking to the temperature distribution in the water domain. I set the temperature of the walls via nvt thermostat and used nve ensemble for water molecules. to extract the temperature distribution in the domain I used the fix ave/spacial command and computed the temperature from T=2/3*Ke/kb, where Ke is the per atom kinetic energy and kb is boltzman constant. the obtained temperature distribution is kind of weird. the temperature of the fluid is in the order of 3000 k!!
It would be kind of you to help me through that.

Fat chance.

Unless you can prove, that there is a bug in lammps, you will have to debug your input by yourself like everybody else.

As usual, complex inputs like yours have to be built and tested in stages and also each piece and subsystems should be checked.

the temp distribution between two walls is:

I suggest you dump the velocities of all the flow
atoms into a dump file. Then you can post-process
to see if you can compute the temperature you
expect. There are several issues to consider
a) flow velocity vs thermal velocities

b) SPC/E may be for constrained molecule, e.g. SHAKE, so

degrees of freedom are important
c) rotational vs translational motion

Once you are sure you know how to compute T correctly, then
you can see if you can get LAMMPS to do the same thing,

Steve