# [gmx-users] Reg: NPT Equilibration of water

Tsjerk Wassenaar tsjerkw at gmail.com
Sun Oct 31 10:17:54 CET 2010

```Hey,

There's quite a number of people that should be rereading their statistics... :p

First of all, there are instantaneous measurements, averages and
fluctuations. If the statistics (mean/fluctuation) are 7 +/- 500, then
that doesn't mean that the average of 7 is an estimate that may be off
by 500, and the true average can lie between -493 and 507. It means
that the average is 7, but that any instantaneous measurement may be
expected to lie between those values, and that two measurements, say
-300 and +300, can not be expected to be from different distributions.
So a pressure of 300 bar measured in one frame does not suggest that
the pressure of the system is higher than that of a system of the same
constitution with a pressure measured in one frame of -300 bar.

Now the real issue is that the average value (7) is an estimate of the
true mean, which is expected to be 1, as that is given as target
value. The average was determined over a short time, 100ps. Every
interval of 100ps will yield a different value for the observation
'average over 100ps', and these averages will have a certain
distribution around an average 'average over 100ps' observations. If
you do ten of these then the average over averages will still only be
an average over 1 ns, which will be a better estimator of the true
mean value, but is still an estimate itself. And this can be continued
ad infinitum. The average observed over a certain time interval will
converge to the true mean with increasing interval length.

The standard deviation of 'average pressures observed over 100ps' is
likely more than a few bar if the fluctuation in instantaneous values
is so huge. In that regard, an observed average pressure of 7 bar over
100ps is probably consistent with a system at a mean pressure of 1
bar.

Hope it helps,

Tsjerk

On Fri, Oct 29, 2010 at 4:58 PM, Mark Abraham <Mark.Abraham at anu.edu.au> wrote:
> On 29/10/2010 7:50 PM, vinothkumar mohanakrishnan wrote:
>
> Hi Mark
>
> How come 7 +/- 500 bar is approximately 1 bar. Can you explain it more
> clearly? It will be of more useful to me to understand the concept.
>
> It's a pretty basic one in scientific measurement. See
> http://en.wikipedia.org/wiki/Error_bar and links thereon, or an introductory
> physics text.
>
> When we make a measurement, that's only an estimate of the true value. When
> someone cites a range of uncertainty, they are acknowledging that there is
> an appreciable chance that the true value, of which they have taken an
> error-prone measurement, lies somewhere in a given range. So if 7 +/-500
> acknowledges that the true value could lie anywhere from about -493 to 507,
> then that's not much less useful than an observation of 1 +/- 500, which
> acknowledges that the true value could lie anywhere from about -499 to 501.
> The details of that chance will vary with exactly what variation quantity is
> being reported (standard deviation, standard error, maximum deviation,
> quartile, etc.).
>
> Here, g_energy is effectively reporting that the instantaneous pressure
> values varied over a range of many hundreds of bar, and that you haven't
> made enough observations to report a reliable average.
>
> Mark
>
>
> One more thing, density (from g_energy command) of water is found to low
> (expected 1000) after NPT is equilibration. why? given below is my average
> density value. Can i proceed with this density value?
>
> Energy                      Average       RMSD     Fluct.      Drift
> Tot-Drift
> -------------------------------------------------------------------------------
> Density (SI)                 979.37    14.3238    14.2431 -0.0525784
> -5.25789
>
> Regards
> Vinoth
>
> On Fri, Oct 29, 2010 at 2:06 PM, Mark Abraham <Mark.Abraham at anu.edu.au>
> wrote:
>>
>> On 29/10/2010 6:39 PM, vinothkumar mohanakrishnan wrote:
>>
>> I am using the semiisotropic pressure scaling because i want the box size
>> to remain the same as that of the original box size in X and Y axis and want
>> to change it only on the Z axis. I am doing this because at a latter stage i
>> want to create an interface with organic solvent where i need the cross
>> sections (X and Y axis length) of both the box should be the same.
>>
>> what i should do now  get the pressure of 1 bar?
>>
>> 7 +/- 500 *is* approximately 1 bar, and hardly any better an approximation
>> than 1 +/- 500. If you want lower fluctuations, use a larger system and run
>> for a much longer time. Or, since you'll have to re-equilibrate once you
>> combine the solvent boxes, don't bother.
>>
>> Mark
>>
>> Has no one has got the pressure of water to be close to 1 bar in GROMACS
>> till now?
>>
>> Regards
>> Vinoth
>>
>> On Fri, Oct 29, 2010 at 12:54 PM, David van der Spoel
>> <spoel at xray.bmc.uu.se> wrote:
>>>
>>> On 2010-10-29 09.08, vinothkumar mohanakrishnan wrote:
>>>>
>>>> Hi Mark
>>>>
>>>> I read the link before posting the question. even though the fluctation
>>>> is between 500-600 (as said in the link) bar the average pressure is
>>>> around 7.4 bar. my concern is the average pressure?. Because at a latter
>>>> stage i am going to combine this water box with another organic solvent
>>>> and i want to avoid any complications there.
>>>
>>> Think again of what you just wrote. Your value is 7 +/- 500. In fact
>>> according to normal statistical rules you should round the value of 7 to 0.
>>> If you want more accurate number you could increase the box size by a factor
>>> of 100.
>>>
>>> By the way, why are you using semiisotropic pressure scaling in a water
>>> box?
>>>
>>>>
>>>> Regards
>>>> Vinoth
>>>>
>>>> On Fri, Oct 29, 2010 at 12:22 PM, Mark Abraham <Mark.Abraham at anu.edu.au
>>>> <mailto:Mark.Abraham at anu.edu.au>> wrote:
>>>>
>>>>    On 29/10/2010 5:48 PM, vinothkumar mohanakrishnan wrote:
>>>>
>>>>        Hi Gromacians
>>>>
>>>>        I want to do equilibration of water (spc model) first in the NVT
>>>>        ensemble and then in the NPT ensemble to maintain a temperature
>>>>        of 300K and a pressure of 1 bar respectively. The NVT
>>>>        equilibration works fine and the average temperature turns out
>>>>        to be 299.229 K.
>>>>
>>>>        Energy                      Average       RMSD        Fluct.
>>>>              Drift            Tot-Drift
>>>>        ------------------------------
>>>>        -------------------------------------------------
>>>>        Temperature                 299.229    10.3092    10.2463
>>>>          0.0393873    3.93877
>>>>        Heat Capacity Cv:       12.494 J/mol K (factor = 0.00118698)
>>>>
>>>>          when i do NPT equilibration i am not getting the desired
>>>>        pressure as 1 bar or atleast close to 1 bar (between 1-1.4 bar).
>>>>        In the mdp file i used semiisotropic pressure coupling type
>>>>        because i want to fix the length of the box same on two axis as
>>>>        that of the original box size and want to change it only on one
>>>>        axis. can any one tell me why iam not getting the desired
>>>>        pressure of 1 bar?.
>>>>
>>>>
>>>>    This looks normal for a smallish water system over 100ps. See
>>>>    http://www.gromacs.org/Documentation/Terminology/Pressure
>>>>
>>>>    Mark
>>>>
>>>>
>>>>
>>>>        Energy                      Average       RMSD     Fluct.
>>>>          Drift  Tot-Drift
>>>>
>>>>  -------------------------------------------------------------------------------
>>>>        Pressure (bar)               7.4339    574.052    573.574
>>>>        0.812129    81.2138
>>>>
>>>>        Given below is my mdp file (NPT equilibration). any help is
>>>>        highly appreciated.
>>>>
>>>>        title           = DCE NVT equilibration
>>>>        cpp           = usr/bin/cpp
>>>>        integrator   = md
>>>>        nsteps       = 100000
>>>>        dt              = 0.001
>>>>        nstxout       = 100
>>>>        nstvout       = 100
>>>>        nstenergy   = 100
>>>>        nstlog        = 100
>>>>        ns_type     = grid
>>>>        nstlist        = 1
>>>>        rlist            = 1.0
>>>>        coulombtype    = PME
>>>>        rcoulomb         = 1.0
>>>>        vdwtype          = Cut-off
>>>>        rvdw               = 1.0
>>>>        pme_order       = 4
>>>>        fourierspacing  = 0.16
>>>>        pbc                 = xyz
>>>>        tcoupl              = V-rescale
>>>>        tc-grps             = system
>>>>        tau_t                = 0.1
>>>>        ref_t                 = 300
>>>>        pcoupl              = berendsen
>>>>        pcoupltype        = semiisotropic
>>>>        tau_p                = 0.5
>>>>        ref_p                 = 1.0 1.0
>>>>        compressibility   = 0.0 4.5e-5
>>>>        DispCorr            = Enerpres
>>>>        gen_vel             = yes
>>>>        gen_temp          = 300
>>>>        gen_seed          = 173529
>>>>
>>>>        Regards
>>>>        Vinoth
>>>>
>>>>
>>>>    --
>>>>    gmx-users mailing list gmx-users at gromacs.org
>>>>    <mailto:gmx-users at gromacs.org>
>>>>    http://lists.gromacs.org/mailman/listinfo/gmx-users
>>>>    Please search the archive at
>>>>    http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>>>>    Please don't post (un)subscribe requests to the list. Use the www
>>>>    interface or send it to gmx-users-request at gromacs.org
>>>>    <mailto:gmx-users-request at gromacs.org>.
>>>>
>>>>
>>>
>>>
>>> --
>>> David van der Spoel, Ph.D., Professor of Biology
>>> Dept. of Cell & Molec. Biol., Uppsala University.
>>> Box 596, 75124 Uppsala, Sweden. Phone:  +46184714205.
>>> spoel at xray.bmc.uu.se    http://folding.bmc.uu.se
>>> --
>>> gmx-users mailing list    gmx-users at gromacs.org
>>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>>> Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>>> Please don't post (un)subscribe requests to the list. Use the www
>>> interface or send it to gmx-users-request at gromacs.org.
>>
>>
>>
>> --
>> gmx-users mailing list    gmx-users at gromacs.org
>> http://lists.gromacs.org/mailman/listinfo/gmx-users
>> Please search the archive at
>> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
>> Please don't post (un)subscribe requests to the list. Use the
>> www interface or send it to gmx-users-request at gromacs.org.
>
>
>
> --
> gmx-users mailing list    gmx-users at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/Search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
>

--
Tsjerk A. Wassenaar, Ph.D.

post-doctoral researcher
Molecular Dynamics Group
* Groningen Institute for Biomolecular Research and Biotechnology
* Zernike Institute for Advanced Materials
University of Groningen
The Netherlands

```