# [lammps-users] Temperature calculation

Hello,

I am new to LAMMPS and am trying to set up a simulation.
I have a mix of water and point-atoms, with some rigid carbon walls.

I set the temperature to a constant 298K using fix nvt:
fix 2 movingatoms nvt temp 298.0 298.0 20.0

When the simulation runs however, I get output indicating that the temperature is only 200K, like this:

---------------- Step 1200 ----- CPU = 128.6769 (sec) ----------------
TotEng = -24676.1103 KinEng = 2724.7333 Temp = 221.1690
PotEng = -27400.8436 E_bond = 0.0000 E_angle = 0.0000
E_dihed = 0.0000 E_impro = 0.0000 E_vdwl = 3303.2261
E_coul = 58147.8695 E_long = -88851.9391 Press = -485.7289

And the temperature does not vary significantly from around 220K. It should be 298K though.
I tried setting the dampening time both to 20.0 and 200.0, but it makes no difference.

The water is done with shake:
fix 1 water shake 0.0001 10 500 b 1 a 1

I’m guessing it could be something to do with the degrees of freedom, but I can’t work out how to fix it.
My frozen atoms are specified like this:
group carbon type 1
fix 3 carbon setforce 0.0 0.0 0.0

I tried telling LAMMPS just to calculate the temperature of the moving atoms by using:
group movingatoms type 2 3 4 5
compute matemp movingatoms temp

but the output is the same.

I would appreciate any insight into what I’m missing here.

Best regards,

Hen

Hello,

I am new to LAMMPS and am trying to set up a simulation.
I have a mix of water and point-atoms, with some rigid carbon walls.

I set the temperature to a constant 298K using fix nvt:
fix 2 movingatoms nvt temp 298.0 298.0 20.0

When the simulation runs however, I get output indicating that the
temperature is only 200K, like this:

[...]

I'm guessing it could be something to do with the degrees of freedom, but I
can't work out how to fix it.
My frozen atoms are specified like this:
group carbon type 1
fix 3 carbon setforce 0.0 0.0 0.0

I tried telling LAMMPS just to calculate the temperature of the moving atoms
by using:
group movingatoms type 2 3 4 5
compute matemp movingatoms temp

but the output is the same.

I would appreciate any insight into what I'm missing here.

you are missing: thermo_modify temp matemp

axel

Axel,
Great, it works now.
Thank you!

Hello,

Can anyone tell me if it’s possible to vary the velocity of some atoms (a wall) in a strict way, depending on the force?

I calculate the force on a wall using:
variable wf1 equal fcm(wall1,z)

Further down, I specify the velocity of the wall to be one of two values:
if {wf1} <= -40.0 then "fix 10 wall1 move linear 0.0 0.0 0.0 units box" if {wf1} > -40.0 then “fix 11 wall1 move linear 0.0 0.0 0.001 units box”

I have monitored the value of wf1, and -40 is an appropriate value.
It seems that this is only read once though, at the beginning of the simulation, since the wall velocity doesn’t change after initial setting.
I can’t figure out how to make this be read at every timestep. If anyone could offer any advice, I would very much appreciate it.

Thank you.

It sounds like you want the fix to evaluate
and change the wall velocity during a run.
The only way to do this is to use the
variable option within fix move.

And assign a variable that evaluates the
fcm() and sets the velocity accordingly.

Something like:
variable fcm equal fcm(wall1,z)
variable vz equal (v_fcm > -40.0) * 0.001
fix 10 wall1 move variable NULL NULL NULL NULL NULL v_vz

Steve

Dear Steve,

I was able to successfully implement the moving wall, however it seems it’s better to move based on an average of, for example, over the last 20 steps, rather than instantaneous values.

I spent some time looking into this, but it’s not clear how to do this.
fix ave/time applies to a group of atoms, whereas I want to simply average a number.
I also considered compute, and I saw in the mail archive [1] that you added something for variables, but I can’t get it to work. It also requires a group of atoms to act on, while I just want to average a number.

Would anyone be able to point me in the right direction?

Many thanks.

You can average anything you want, like a variable, with fix ave/time.
The group will just be ignored.

Steve

Hello,
Thanks for your help. I got a step closer however I’m having trouble using the averaged variable.

If you do this:
fix 9 all ave/time 1 1000 1000 v_wf
variable wfav equal f_9

Then it produces the error: Fix in variable not computed at compatible time

So until the 1000th step, I need it to ignore the variable command.

I understand this, but I can’t see how to ensure that the variable is only evaluated in concert with the fix.
I can’t see a timestep variable to use, so when I tried this:

variable wfdav equal 0.0
fix 9 all ave/time 1 1000 1000 v_wfd
if timestep equal 1000 “variable wfdav equal f_9”

then the “if” command produced an error. So I guess “timestep” is not a real variable.

Furthermore, there appears to be no option within the “variable” command to only run it on certain steps.
I’m really stumped on this one. I could start trying to build something into the source code, but I’ve seen references to this problem within the manual so it must be a common issue - unfortunately, although it’s referenced many times within the manual, it doesn’t say anywhere how actually to deal with such a case.

Would anyone be able to give any suggestions?

Furthermore, once I reach step 1001, can I refer to the fix, and if so, will the fix yield the same value at step 1000 (ie, the average of 1-1000) or is it a moving average (ie, step 2-1001)?

Thank you again very much for your help.

Well, the fix as you've written it only
produces a value every 1000 steps,
so you can't access it on any arbitrary
step. The fix ave/time command has many
options that alter how it averages, delays
the start of it, alter when it does output
etc. I suggest you read its doc page
carefully. The error messages will
tell when you do something wrong.

Steve

Hello. I have a short question about fix deform.

I’m trying to get my simulation box to increase / decrease with size, depending on the value of a variable. I first tried this:

variable boxvel equal 0.00001
fix boxbigger all deform 1 z vel v_boxvel units box remap none

but nothing happens. When I do any of the following however, it does work:

fix boxbigger all deform 1 z vel \${boxvel} units box remap none # based on initial value of boxvel only
OR
fix boxbigger all deform 1 z vel 0.00001 units box remap none

Does this mean it’s not possible to assign a variable velocity to this fix? If anyone could offer any insight regarding this, I’d very much appreciate it!
Best regards.

Does this mean it's not possible to assign a variable velocity to this fix?

If you mean specify a constant like 0.00001 to a variable and use
that in the fix deform command, then the way you did it with \${boxvel} is
the way to do it.

If you mean assign a formula to a variable (see the variable doc page)
and use that as the time-dependent box deformation parameter, then no,
the fix deform command cannot do that. Several commands like
fix setforce do allow this by using v_boxvel in the command, but
not fix deform. If it did, it would be in the doc page description
of the command. What happened when you did this is LAMMPS just
read v_boxvel as a number and got 0.0.

Steve

Dear Steve,

Thanks very much for clarifying that. I appreciate it.

The way I got around this, for anyone who might have such a problem in the future, was to run lots of mini simulations in a loop. At the beginning of the new mini simulation I refer to the variable using \${boxvel}, the value of which is determined during the previous simulation.

Thanks again, and happy christmas / new year.