problem with temperature control

Hi all,

I met a problem of temperature control when constructing the liquid/vapor interface. I used the “fix spring/self” command to fix the atoms of upper and lower walls, and used “fix langevin” command to control the two walls at different temperatures (1.1 and 0.75) after equilibration. Then I found that the temperature of liquid between two walls is much higher than both walls, which is not reasonable. Can someone tell me what the problem is? Thank you in advance.

The following is my input file:

2-d LJ flow simulation

dimension 3
boundary p p s
atom_style atomic
neighbor 1.0 bin
neigh_modify delay 5

create geometry

lattice fcc 0.05
region box block 0 50 -13 13 0 140 units box
create_box 3 box
mass 1 1.0
mass 2 4.0
mass 3 4.0

LJ potentials

pair_style lj/cut 3.0
pair_coeff * * 1.0 1.0 3.0

define groups

region 1 block 0.5 49.5 -12.5 12.5 1.3 138 units box
create_atoms 1 region 1 units box

lattice fcc 2.6

region 2 block 0.5 49.5 -12.5 12.5 0 1.3 units box
create_atoms 2 region 2 units box
group lower region 2

region 3 block 0.5 49.5 -12.5 12.5 138 140 units box
create_atoms 3 region 3 units box
group upper region 3

group boundary union upper lower
group flow subtract all boundary

initial velocities

velocity flow create 1.1 1081275 units box

fix 1 all nve
fix 2 flow langevin 1.1 1.1 0.02 699483
fix 8 lower langevin 1.1 1.1 0.02 699483
fix 9 upper langevin 1.1 1.1 0.02 699483

fix 105 lower spring/self 3200
fix 106 upper spring/self 3200

compute Tupper upper temp/com
compute Tlower lower temp/com
compute Tflow flow temp/com

thermo_style custom step temp epair emol etotal press vol c_Tupper c_Tlower c_Tflow

Run

timestep 0.001
thermo 100

dump 1 all atom 5000 dump.flow
dump_modify 1 scale no

run 20000

unfix 2
unfix 8
unfix 9

fix 8 lower langevin 0.75 0.75 0.02 699483
fix 9 upper langevin 1.1 1.1 0.02 699483

run 100000

Best regards
Wenjing Zhou

Don't know - don't see anything obviously wrong with
your script. I would simplify the model, starting
with the basics, and add features one at a time.

Steve

Thanks, Steve. Another question: I know that using the combination of “fix temp/rescale” and “compute temp/region” commands I can thermostat atoms in a geometric region, but I want the temperature to be calculated after subtracting the spatially-averaged velocity in the region during the thermostat. What kind of command should I use? or I have to write my own “compute” and “fix”? Thank you.

Best regards
Wenjing Zhou

In LAMMPS lingo, both the COM velocity and the
region definition are "biases" on the temperature.
There is a compute temp/com and temp/region that
does this. You're asking for a temperature that
subtracts both biases. Can't think of a way to
do this at the moment. Maybe there would be a clever
way to apply multiple biases if the computes were
restructured a bit. But for now I think you'd need
to write your own compute which applied both biases.

Steve