Hello,
I am aware that:
compute temp/profile out bin
cannot correctly calculate temperature profiles when “fix shake” is used (problem when one atom in the fixed bond is in one bin, and the other atom is in an adjacent one).
I am very willing to write, and share a code to correctly calculate temperature profiles when “fix shake” is applied. Looking around on the LAMMPS mailing list, I have found that “compute temp/profile out bin” does not remove the degrees of freedom affected by “fix shake”. The only problem is that I do not really know what this implies/changes in the calculation of each atom’s KE:
KE = 3/2 k * mass * (vx^2 +vy^2 + vz^2)
If anyone could provide any helpful suggestions, once I finish the code, I will gladly share it.
Thanks a lot
Hi Thomas,
The problem is not in the definition of KE, it's in the definition of the number of degrees of freedom. Specifically, when a (potentially complicated) molecule straddles one of the bin boundaries for your profile, how many of its degrees of freedom go into each bin?
Consider I have a molecule of hydrogen gas (H2), the bond length of which I
constrain with fix shake. Why not only consider the barycenter of the
constrained bond length ? Even if the molecule of H2 straddles two bins,
could we not take the barycenter of the H-H bond, and put all of the DoF
associated with that bond in one of the two bins?
This reasoning could be generalized to a more complicated molecule, I am
assuming ?
2013/8/7 Niall Jackson <[email protected]...>
Absolutely. Or, for something simple like H2, you could put half of the degrees of freedom on each atom. You just need to make sure you're consistent in how you partition the DoFs. Obviously, if your thermostat is working according to different rules, then you'll probably find discrepancies between your input temperatures and your temperature profile.