Funny output from temp/profile with SHAKE

LAMMPS LAMMPS Mailing List Mirror
1

Dear all,

I am running a simulation to see streaming velocity profiles on Drude polarizable molecules with rigid SHAKE C-H bonds in a nano-channel. Things seem to be working (amazingly), including a temp/profile compute which I have set to thermostat my system at 300K via an NVT Nose-Hoover integrator for the Drude cores (with a separate NVT Nose-Hoover integrator for the Drude particles), but outputting the values of the temp/profile thermostat gives temperatures around 260K. (I have attached the relevant parts of my input script below.)

What is odd is that my system’s dynamics appear to be correctly 300K. A temp/drude compute during the simulation shows the cores’ temperature at 300K, as well as a temp/com compute rerun after the simulation on velocity dump files. I think this is related to SHAKE: when I rerun on the velocity dump files with fix shake my temperature computes are around 300K, but when I rerun on the velocity dump files without fix shake all the computes show ~ 260K temperatures, which is the same range as my in-simulation temp/profile values and also theoretically sensible (if the temp computes think there are more DOFs than there actually were during the simulation, it will calculate a lower temperature).

Is this correct / incorrect behavior, and what should I do about it (if anything)? My system shows similar dynamics under a temp/com thermostat, so I don’t think the temp/profile thermostat is doing anything wrong internally.

Cheers,

Shern

[Input script snippets:

group sol type 2:29

group MOVATOMS type 2:16

group DRUDES type 17:29

group SHAKECH type 4 5 6 7 15 16

fix SHAKE SHAKECH shake 0.0001 20 2000 t 15 16

compute TATOMS MOVATOMS temp/profile 1 0 0 z 5 out bin # various values of xflag, yflag, zflag, and nbin don’t change the behavior

fix DIRECT all drude/transform/direct

fix NVT1 MOVATOMS nvt temp 300. 300. 100

fix_modify NVT1 temp TATOMS

fix NVT2 DRUDES nvt temp 1. 1. 20

fix OUTATOMS sol ave/time 1000 1 1000 c_TATOMS[*] mode vector ave one file out_tstat.txt

fix INVERSE all drude/transform/inverse

]

2

fix shake removes degrees of freedom. and the number of degrees of freedom has to be known when you compute temperature. however, not all temperature computes can determine how many degrees of freedom have to be subtracted.
also when you use rerun without shake, you add degrees of freedom, that were not present and that lowers your temperature.

axel.