[lammps-users] Inconsistent heat flux and thermal conductivity with shake and fix rigid/nve/small

Dear All,

I’m trying to calculate the thermal conductivity using the Green Kubo method for an ensemble of rigid bodies (fix rigid/nve/small), for which I get very strange results. To test this, I considered an ensemble of dimers and got very different results (1-2 orders of magnitude) when using shake vs. fix rigid/nve/small.

To further investigate this, I tried a system with one dimer with the same initial conditions. I get the same PE (identically 0) and KE (which is now conserved) and similar stress/atom elements. In this case, the stress contributions come only from the constraint fixes. However, the heat flux and thermal conductivity are very different. I’ve attached the input files, and appended the salient features to the end of this email. I’m using LAMMPS version 23Jun18, and I get the same results for the single dimer using 3Mar20.

Why is there a difference between shake and fix rigid/nve/small? Is it possible to calculate the correct heat flux and thermal conductivity with fix rigid/nve/small?

Any help would be much appreciated.



in.dimer_rigid (1.68 KB)

mol.dimer (69 Bytes)

in.dimer_shake (1.65 KB)

mol.dimer_shake (91 Bytes)

Are the instantaneous values from compute heat/flux (variables JQxyz) different on every timestep from the beginning?
E.g. use thermo 1