PIGLET Thermostat with More Than 64 Beads

[hein...@uw.edu]

Hi,

I am in a situation where I would like to use the PIGLET thermostat to get the decreased number of beads compared to other thermostats, but the simulations we are doing are sufficiently low in temperature that the number of beads required is larger than 64 beads. This is only a problem because, as far as I know, the only place to get the covariance and friction matrices for the PIGLET thermostat is from gle4md.org, and the largest number of beads available in the library there is 64.

I am wondering if there is any other place that these matrices can be found for a larger number of beads, or if there is a simple way to expand the matrices provided in the library at gle4md to a larger number of beads?

Specifically, for water in the gas-phase as low as 5K, it seems that hbar*omega/kT=2000 would be necessary and up to 128 beads might be required to reach convergence of the zero-point energy. For temperatures around 50K, something closer to hbar*omega/kT=200 with about 100 beads would be required to reach convergence of the ZPE for some of the floppier possible systems.

I would love to know if there is a way to get the necessary matrices for these larger number of beads.

Thanks,

Joe

[Michele Ceriotti]

Dear Joe,

    Together with Felix Uhl from Dominik Marx group we did test PIGLET for ultra-low temperature conditions (see http://dx.doi.org/10.1063/1.4959602). He fitted matrices for extreme sets of conditions and perhaps he can add them to gle4md.org. However, we recently also implemented high-order PIMD, which might also help a lot reducing the number of beads, and can perhaps be more accurate in the end. Venkat Kapil would be the best person to give you advice, but in all cases you find it implemented and with examples in the current -dev master branch on github. 

What level of accuracy do you need, and what kind of observables you want to compute?

Michele

[hein...@uw.edu]

Thanks for the recommendations. We will only be calculating static properties (average energies, geometrical data, etc.), but ideally we would like to converge ZPE to ~0.1 kcal/mol or lower.

If those matrices could be added to gle4md.org that would be ideal. I have contacted Felix Uhl by email asking what he would think about sharing these. If he is fine with this will he need to be in touch with you about adding them to the library?

Thank you very much for the help. I think it will be very useful to many people if these matrices can be made available. I don't want to step on any toes though, so no pressure.

[Michele Ceriotti]

Look, I could generate equivalent matrices in an afternoon... but who has an afternoon free these days :-) Let's see what Felix says, it would be much nicer to have a set of parameters that has already been put to the test. However read carefully his article. At those ultra-low temperatures, PIGLET gives you very quickly qualitative accuracy (say 5-10%) but if you're after a 0.1kcal.mol (which I'd say depending on whether you count per atom or per molecule it 's between 5 and 1%) you might have a better luck with high-order PIMD. At least, you should give them a go, which you can do already now: just use <dynamics mode="sc"> and print out sc-op and sc-td estimators. Purely structural properties (pair correlation functions &C) can be readily computed by taking just the even beads from the output. 

[hein...@uw.edu]

Felix responded and it seems many of the matrices I would need for temperatures between 50K and 5K would not be suitable because of the maximum frequency they can accommodate. He also mentioned that the matrices he has which were used for around 1.5K would likely not be suitable for the 5K simulations. I am not sure why this is. Can't the parameters be scaled in some way for different temperatures? In any case, the max frequency for those might be very high.

I will look into running some high-order PIMD, but I still think it would be ideal to get some of these PIGLET matrices. Your point about the accuracy at very low temperatures is well taken. I do not know that this will be as much of a problem, however, as I would be content to hit temperatures of 50K and 25K and then look into the 5K data to see if it will be robust or not. In any case, the simulations done at 20K in the paper you linked do not seem to suffer from this same accuracy issue.

[Michele Ceriotti]

Hi Joe, 

   to clarify: matrices are built for a given value of beta hbar omega_max, and then depending on the temperature the actual cutoff varies. 

Now, this cutoff is not to be taken too strictly. If you're doing water and it says it is 3500cm^-1, then you should be fine as it's not like at 4000 it

will be total garbage. If it says it's 5000cm^1, you should be also fine. It gets problematic when the cutoff is *much* higher, since at that point

you might need a shorter time step to integrate accurately, as the noise will contain many high frequency components. 

So I'm actually convinced that matrices you have already online can get you quite far -- there are even some, marked "experimental" that get

down to hbar omega max=2000. 

My point was more that you probably won't need to go as high as 128 beads. PIGLET at those cryo T is best to get qualitative results. 

All of that said, keep an eye on GLE4MD, will be adding a few more matrices over the next few days; your milage may however vary, as I won't have 

time to do much testing.

Best

Michele