Details on fix npt (nph) drag


from what I could figure out before writing the thread you already mentioned, the drag parameter is not just a proper damping term on the PDE representing the controller chains, but rather a scaling factor on its output, as you assume. I haven't found any equations or at least comments on it anywhere either.

2) In case there is no such thing in writing, I guess the answer has to be looked for in the code. Looking through the source files, I am guessing that this "dragging" is done in fix_nh.cpp. A pdrag_factor is created that seems to remain 1 without drag (i.e.,drag=0.0). Otherwhise, 1 is reduced by a factor depending on the time step, Pdamp?, drag, and the number of chains (e.g., line 689).

I think this was just supposed to be an additional PT1 with a time constant of pdamp/chain (in time steps).

Next, pdrag_factor is multiplied to etap_dot and omega_dot, which should be the barostat velocities (omega_dot) and the thermostat velocities belonging to the barostat (etap_dot).

Depending on whether the MTK formalism is used, it may be multiplied in again, and at that point, I still think something is double-counted, leading to the missing RT/vol pressure I observed.
From what I could find, the dynamics is unaffected (apart from the distributions being "narrower" by a deterministic factor, as expected), and in theory one could just pre-add this factor to the setpoint to get correct enthalpies. I have no idea what that does to other fixes that use pressures...

> When dealing with constant pressure simulations, it might be helpful to use drag to settle your system into the right average pressure and then
> sampling with a proper barostat without drag once you reach a steady state. That in itself might
> already significantly reduce the fluctuations. If they are still too bad, you can also experiment
> with the number of barostats through the pchain option. As far as I know, barostat chains will
> sample the right ensemble.

At least for metallic systems with a bulk modulus in the order of >10GPa, that doesn't help at all. The catastrophic fluctuations are just intrinsically there, even when starting from already stabilized states. The parameters only make them have shorter or longer time constants. The only real correct solution is using absurdly large systems - in my experience: if under NPT the system is down to 0.1 mK variance, the barostat might at least be in the right order of magnitude.

It always seemed odd to me, because the 1994 MTK article shows much larger volume fluctuations (corresponding to lower pressure fluctuations) than I have ever seen in Lammps, for a softer system. In consequence, the variance-based methods for Cp, beta etc. don't really work, as is regularly discussed on this list.