[gmx-users] Reference pressure and pressure fluctuations in an NPT simulation
Mark.Abraham at anu.edu.au
Sat May 14 03:38:53 CEST 2011
On 14/05/2011 11:02 AM, Andrew DeYoung wrote:
> I have been running simulations of 254 SPC/E water molecules using the OPLS
Bear in mind that this is a tiny system, and that temperature and
pressure are macroscopic quantities. Reliable measurements of them need
long times and lots of statistics.
> force field. As Nilesh mentioned earlier today, I am using the NPT
> ensemble. I am setting the reference pressure at 1 bar, but when I do a
> sequence of minimization, equilibration, and dynamics steps, and then
> finally use g_energy to determine the average pressure, I find an average
> pressure of ~2 bar, rather than 1 bar.
> Earlier today, Justin pointed out that this probably means that the system
> has not been equilibrated long enough, so I did another run with longer
> equilibration. However, I get a similar, puzzling result of ~2 bar.
> Finally, I tried changing the temperature and pressure coupling methods;
> this time, I get a seemingly more reasonable result for average pressure,
I think that is partly fortuitous, and partly the nature of Berendsen
> but, as I will describe below and in my PDF file, still somewhat puzzling to
> me (I am new to the field of molecular dynamics).
> If you have time, could you please look at my PDF file at
> http://www.andrew.cmu.edu/user/adeyoung/may13.pdf that summarizes what I
> have tried? In case it is helpful, also, here is a text description of what
> I have tried:
> (i) T coupling = v-rescale; P coupling = parrinello-rahman. Reference
> pressure = 1 bar. 1 ns equilibration; 2 ns dynamics. Gromacs tells me that
> the average pressure is ~2.77 bar, but when I use g_energy to extract
> pressure as a function of time to an xvg list and then use software such as
> Mathematica or Matlab to compute the mean of the list, I find that the
> average pressure is ~1.24 bar.
See second para of g_energy -h.
> (ii) Use the same parameters as before, except equilibrate for 5 ns instead
> of 1 ns. Gromacs says that the average pressure is ~2.82 bar, but when I
> extract the pressure data and compute the mean of the list, I find that the
> average pressure is ~6.81 bar.
> (iii) Use the same parameters as in (ii), except use berendsen temperature
> coupling and berendsen pressure coupling. Gromacs says that the average
> pressure is ~1.00 bar, but when I extract the pressure data and compute the
> mean of the list, I find that the average pressure is ~3.19 bar.
> In my PDF file at http://www.andrew.cmu.edu/user/adeyoung/may13.pdf, I have
> plotted pressure vs time for the dynamics runs of (i), (ii), and (iii), and
> these plots show VERY large magnitude oscillations.
Looks normal. I see that Justin's reply has dealt with the issues below.
> If you have time, I have three questions about these results:
> (1) Does it seem reasonable that I obtain an average pressure closer to my
> reference pressure (1 bar) when I use berendsen temperature coupling and
> berendsen pressure coupling -- instead of v-rescale temperature coupling and
> parrinello-rahman pressure coupling?
> (2) When I use g_energy to extract pressure vs time, and then compute the
> mean of the list of pressures, why is my answer so different from Gromacs'
> "black box" calculation (that is, using g_energy to have Gromacs simply
> print the average pressure over the dynamics run) of the average pressure?
> For all of my runs, I am using a step size of 1 fs (=0.001 ps). However,
> when I extract the pressures (using the default settings), I get a pressure
> value for every 0.1 ps, not 0.001 ps. Could it be that the difference in
> average pressures is due to the fact that "under the hood" Gromacs is using
> the pressure data at every step (0.001 ps), instead of every 0.1 ps? But,
> even if this is true, I still wouldn't necessarily expect such a big
> disagreement between the two calculations.
> (3) When I plot pressure vs time (as I have done in my PDF file), why is
> there such a large magnitude of pressure fluctuation? For example, in my
> runs, the maximum pressure is on the order of 3000 bar, whereas the minimum
> pressure is on the order of -3000 bar. These pressures are unreasonably
> large in magnitude (despite the fact that the average pressure nevertheless
> turns out to be of the correct order of magnitude in the long run). Is this
> true? Also, is "negative pressure" unphysical? Or, does "negative
> pressure" correspond to "compression" and "positive pressure" corresponds to
> "expansion," or something like this?
> Thank you very much for your time. I truly appreciate it.
> Andrew DeYoung
> Carnegie Mellon University
More information about the gromacs.org_gmx-users