LAMMPS WWW Site - LAMMPS Documentation - LAMMPS Mailing List Archives
[lammps-users] viscosity calculation thermostat problem
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[lammps-users] viscosity calculation thermostat problem


From: Dezhao Huang <dhuang2@...46...>
Date: Thu, 19 Apr 2018 17:17:52 -0400

Dear Lammps users,
Sorry to ask an old question. I have searched the email list but still can not get a good explanation.
I am currently trying to use the shearing the wall method to calculate the shear viscosity of water which is confined between two SAM-functionalized gold substrate. If I only use fix nve for the water, I can get a linear velocity profile from the simulation, but the temperature of my water increases a lot, which I think is not right.
The problem is when I try to use (fix nve + fix langevin) for the water, I could not get the linear velocity profile. The velocity points for my water is just randomly dispersed around 0.
I would be really grateful if you can give me some hints or explanation.

Below is my lammps script (fix nve + fix langevin for my water)

read_restart       npt.restart
restart               1000 left.restart right.restart                               
pair_style          hybrid lj/class2/coul/long 10.0 10.0   lj/cut/coul/long 10.0 10.0    morse 8      
pair_modify        mix arithmetic  
kspace_style       pppm 1.0e-5                                    
bond_style          hybrid class2 harmonic                                               
angle_style         hybrid class2 harmonic                                             
dihedral_style     class2                                             
improper_style    class2  

include           'pair_coeff'
include           'bond_CH3'

group              SAM1 id <> 1 1280
group              SAM2 id <> 6081 7360
group              Gold_1 id <> 1281 6080
group              Gold_2 id <> 7361 12160
group              Water id <> 12161 18160
group              AuS type 3 5
group              Au union  Gold_1 Gold_2
group              bulk1 union Gold_1 SAM1
group              bulk2 union Gold_2 SAM2
group              wall  union bulk1 bulk2
group              moving union Water bulk1

#---------------  eighborlist------------------------------------
neighbor           2.5 bin
neigh_modify       delay 0 every 1 check yes
#--------------variable definition and initilization------------
variable            vtimpstep equal 0.25
timestep           ${vtimpstep}       
variable           Pdamp  equal  ${vtimpstep}*1000
variable           Tdamp  equal  ${vtimpstep}*100

variable        shearrate equal 0.0015

variable           atom2Pa  equal  101325
variable           fs2s     equal  1.0e-15
variable           ly        equal    34
variable           npt_run_time_1  equal   3000 
variable           npt_run_time_2  equal   3000
variable           runtime         equal   2000000
variable           PaS2mPaS        equal   1000

variable           Text      equal 300
variable           Pext      equal 1.0
variable           length    equal 38
fix                  constrain Water shake 1.0e-4 100 0 b 6 a 10

fix                   2 Water nve

compute            Water_thermo  Water temp/profile 0 1 1 x 20

fix                    1 Water  langevin 300 300 ${Tdamp} 123456 tally yes
fix_modify         1 temp Water_thermo


compute            1  Water stress/atom Water_thermo
compute            11 Water reduce sum c_1[1]
compute            12 Water reduce sum c_1[2]
compute            13 Water reduce sum c_1[3]
compute            14 Water reduce sum c_1[4]
compute            15 Water reduce sum c_1[5]
compute            16 Water reduce sum c_1[6]

variable           waterpress  equal  -(c_11+c_12+c_13)/(3*lx*ly*38)

fix                   move  bulk1 move linear ${shearrate} 0.0 0.0 units box

compute          layers Water chunk/atom bin/1d z lower 0.02 units reduced
fix       chunkavg Water ave/chunk 10 500 5000 layers vx c_1[4] c_1[5] c_1[6] file profile.wall

variable           viscxz  equal    ((c_15*${atom2Pa}/(lx*ly*${length}))/(${shearrate}/(${length}*${fs2s})))*${PaS2mPaS}

fix           visave  all ave/time 10 100 1000  v_viscxz  ave running file avg_vis.profile

thermo                1000
thermo_style    custom step c_Water_thermo temp press v_waterpress pxy c_14 c_15 c_16 f_1  v_viscxz  f_visave
thermo_modify      temp Water_thermo

dump               2 all custom 4000 shearCH3.lammpstrj id type x y z vx vy vz fx fy fz

run                  4000000

write_data         shear.data
write_restart      shear.restart



Thanks!

******************************************************************
Dezhao Huang
2nd year Ph.D. student
University of Notre Dame
Department of Aerospace and Mechanical engineering
150D CoMSEL Fitzpatrick Hall
Notre Dame, IN 46556-5710