[gmx-users] Dubious results with NPT

[Michael Shirts]

No full solutions to this problem -- I'll write a few notes.

* With large systems run for relatively small amounts of time, the
ensemble could be statistically indistinguishble from the true
distribution even if the averages don't line up correctly. See:

http://pubs.acs.org/doi/abs/10.1021/ct300688p

For ways to examine the full distribution of volumes.

* Parrinello-Rahman isn't quite correct in GROMACS because virial and
kinetic energy are offset by 1/2 a step, though getting something more
like 1.2 atm when running at 1 atm is more common.

* Note that with systems that are heterogeneous in direction (not
uniform in x y and z), any errors in the pressure may be magnified.

* If you have a property that changes when the pressure changes from 1
to 5 atm, odds are that molecular simulation is not the right way to
study it -- the force fields are not that accurate.


On Tue, Feb 25, 2014 at 6:48 AM, Dr. Vitaly Chaban <vvchaban at gmail.com> wrote:
> in that case the average is reasonable, but you should state in the
> manuscript that all the reported results correspond to 5.15 bar.
>
>
> Dr. Vitaly V. Chaban
>
>
> On Tue, Feb 25, 2014 at 10:42 AM, sujithkakkat . <sujithks58 at gmail.com> wrote:
>> Hi,
>>
>>   In fact the average pressure 5.15 bar was obtained  for the case where the
>> reference pressure was 5 bar.  I had mentioned it in one of the previous
>> mails. But still the average value is not that good.
>>
>>>I had run an earlier simulation with the same parameters for a pure water
>>> system for 5ns and reference pressure 5bar, and things worked fine there
>>> with the average pressure at 5.15bar.
>>
>>
>> Sujith.
>>
>>
>> On Tue, Feb 25, 2014 at 2:29 PM, Dr. Vitaly Chaban <vvchaban at gmail.com>
>> wrote:
>>>
>>> your average pressure is the pressure you should report in the
>>> publication. if you got 5 bars instead of 1 bar, you should write
>>>
>>> "I simulated the system at 5 bars"
>>>
>>>
>>> Dr. Vitaly V. Chaban
>>>
>>>
>>> On Tue, Feb 25, 2014 at 6:01 AM, sujithkakkat . <sujithks58 at gmail.com>
>>> wrote:
>>> > On Mon, Feb 24, 2014 at 5:37 AM, sujithkakkat . <sujithks58 at gmail.com>
>>> > wrote:
>>> >> Hello,
>>> >>
>>> >>   Thank you both for  the comments. I am using gromos96 forcefield . I
>>> >> read
>>> >> a little bit  and as you said the nonbonded cutoff has to be higher.
>>> >> The tau_p=5ps was chosen , since the manual mentions that the value has
>>> >> to
>>> >> be raised by 4-5 times on going from berendsen to parrinello-rahman
>>> >> barostat, though I did not completely follow the reasons behind it. I
>>> >> will
>>> >> try with lower values.
>>> >>
>>> >>    I had run an earlier simulation with the same parameters for a pure
>>> >> water
>>> >> system for 5ns and reference pressure 5bar, and things worked fine
>>> >> there
>>> >> with the average pressure at 5.15bar.
>>> >
>>> >
>>> > No, this average pressure is not fine.
>>> >
>>> >
>>> >
>>> >  Why isn't it fine? I admit it is not very good. I guess due to huge
>>> > fluctuations the average pressure is tend to vary a bit to either side
>>> > of
>>> > the reference value, unless I continue the equilibration for a very long
>>> > time. Also the average was calculated with g_energy, which I guess
>>> > considered the entire simulation, where the first few time steps, when
>>> > the
>>> > system is away from equilibrium,  would have poor pressure values with
>>> > respect to reference , which would reflect in the average pressure
>>> > calculated.
>>> > Also using g_analyze to calculate average pressure between different
>>> > times
>>> > during the simulation did not give the same value.  So I am not sure if
>>> > I
>>> > should aim at an average pressure which is equal to the reference
>>> > pressure.
>>> >
>>> > Thanks,
>>> > Sujith.
>>> >
>>> >
>>> >
>>> > On Mon, Feb 24, 2014 at 4:26 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>> >>
>>> >>
>>> >>
>>> >> On 2/23/14, 11:37 PM, sujithkakkat . wrote:
>>> >>>
>>> >>> Hello,
>>> >>>
>>> >>>    Thank you both for  the comments. I am using gromos96 forcefield .
>>> >>> I
>>> >>> read
>>> >>> a little bit  and as you said the nonbonded cutoff has to be higher.
>>> >>> The tau_p=5ps was chosen , since the manual mentions that the value
>>> >>> has
>>> >>> to
>>> >>> be raised by 4-5 times on going from berendsen to parrinello-rahman
>>> >>> barostat, though I did not completely follow the reasons behind it. I
>>> >>> will
>>> >>> try with lower values.
>>> >>>
>>> >>>     I had run an earlier simulation with the same parameters for a
>>> >>> pure
>>> >>> water system for 5ns and reference pressure 5bar, and things worked
>>> >>> fine
>>> >>> there with the average pressure at 5.15bar. I guess the sigma  for the
>>> >>> case
>>> >>> of water was low and therefore the small cut-off of 0.8nm did not
>>> >>> matter.
>>> >>> However the case of cyclohexane alone remains to be tried.
>>> >>>
>>> >>> I guess Dr Vitaly was saying about using a switch/shift function.
>>> >>> I will try the simulation with the new settings and see.
>>> >>>
>>> >>
>>> >> The Gromos force fields were parametrized with a plain cutoff for van
>>> >> der
>>> >> Waals.  They do not require switching or shifting.  It's true that
>>> >> plain
>>> >> cutoffs are a bit harsh, but switching and shifting have their own
>>> >> issues,
>>> >> and one should always stick to the parametrization methods unless there
>>> >> is
>>> >> demonstrable evidence that the new method is better.  I have never
>>> >> known
>>> >> Gromos96 force fields to need anything other than plain cutoff for van
>>> >> der
>>> >> Waals to function as expected.
>>> >>
>>> >>
>>> >> -Justin
>>> >>
>>> >> --
>>> >> ==================================================
>>> >>
>>> >> Justin A. Lemkul, Ph.D.
>>> >> Ruth L. Kirschstein NRSA Postdoctoral Fellow
>>> >>
>>> >> Department of Pharmaceutical Sciences
>>> >> School of Pharmacy
>>> >> Health Sciences Facility II, Room 601
>>> >> University of Maryland, Baltimore
>>> >> 20 Penn St.
>>> >> Baltimore, MD 21201
>>> >>
>>> >> jalemkul at outerbanks.umaryland.edu | (410) 706-7441
>>> >> http://mackerell.umaryland.edu/~jalemkul
>>> >>
>>> >> ==================================================
>>> >> --
>>> >> Gromacs Users mailing list
>>> >>
>>> >> * Please search the archive at
>>> >> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
>>> >> posting!
>>> >>
>>> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>> >>
>>> >> * For (un)subscribe requests visit
>>> >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
>>> >> send
>>> >> a mail to gmx-users-request at gromacs.org.
>>> >
>>> >
>>> --
>>> Gromacs Users mailing list
>>>
>>> * Please search the archive at
>>> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!
>>>
>>> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>>>
>>> * For (un)subscribe requests visit
>>> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send
>>> a mail to gmx-users-request at gromacs.org.
>>
>>
> --
> Gromacs Users mailing list
>
> * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-request at gromacs.org.