[gmx-users] Nstlist and constrain simulations

[Szilárd Páll]

--
Szilárd


On Tue, Feb 2, 2016 at 4:23 PM, Mark Abraham <mark.j.abraham at gmail.com> wrote:
> Hi,
>
> On Tue, Feb 2, 2016 at 1:11 PM Michail Palaiokostas Avramidis <
> m.palaiokostas at qmul.ac.uk> wrote:
>
>> Hi Mark and thank you for your answer.
>> Please see below :)
>>
>> On 01/02/16 18:28, Mark Abraham wrote:
>> > Hi,
>> >
>> > On Mon, Feb 1, 2016 at 6:42 PM Michail Palaiokostas Avramidis <
>> > m.palaiokostas at qmul.ac.uk> wrote:
>> >
>> >> Dear GMX users,
>> >>
>> >>
>> >> I would like to ask about your opinion on the size of the neighbour list
>> >> (nstlist).
>> >>
>> >>
>> >> I am running constraint simulations
>> >
>> > What do you mean by a "constraint simulation?"
>> Sorry, I should have been clearer. I am using the z-constraint method in
>> which I constrain the permeant along the z-axis in various positions and
>> I record the constraint force. Similar to umbrella sampling but with
>> constrained distances between the permeant and membrane.
>>
>
> Using which GROMACS feature? e.g. freeze groups and NPT are usually an
> unhappy combination.
>
>
>> >> of a permeant along a lipid bilayer. In the initial setup the system
>> runs
>> >> on NPT and nstlist 10 (which GROMACS changes to 40 automatically). With
>> >> this setup my simulation is very fast but always crashes at various
>> points
>> >> with the same way always. Initially it gives some LINCS warnings and
>> then
>> >> it gives an error that a particle "communicated to PME rank 2 are more
>> than
>> >> 2/3 times the cut-off out of the domain decomposition cell".
>> >>
>> > Usually this indicates your setup is unstable e.g. see
>> > http://www.gromacs.org/Documentation/Terminology/Blowing_Up. I'll
>> proceed
>> > by assuming you're very confident that your membrane system has had a
>> > suitable few dozen+ nanoseconds of equilibration ;-)
>> Yes absolutely confident. The system does not crash with unconstrained
>> simulations. The membrane-solvent system is well equilibrated for 1μs
>> and then, after I introduce the permeant, I perform an energy
>> minimization and a short 500ps NPT-Berendsen equilibration to relax the
>> pressure that the permeant might have introduced.
>>
>
> I'd be skeptical that 500ps was enough.
>
>
>> > When I see the last pdb steps before crash, sometimes there seem to be
>> >> overlaps between the permeants and the lipids and then system explosions
>> >> and PBChaos. Other times it is not so dramatic.
>> >>
>> >>
>> >> So the question is, should I be conservative with the nstlist? When I
>> run
>> >> NPT simulations and change the nstlist to 1, the system does not fail.
>> >> Alternatively, if I change to NVT and nstlist to 10, again the system
>> runs
>> >> successfully. But NVT and 40 crashes.
>> >>
>> > This might all be a wild goose chase. If you are pulling permeants into
>> the
>> > membrane, then that creates pressure on the membrane and perhaps
>> > destabilizes the pressure coupling. If that's with Parrinello-Rahman then
>> > it can easily oscillate out of control, ending up violating the
>> assumptions
>> > under which the code is written.
>> Yes, this is why I run an small equilibration in the beginning. And to
>> be honest, technically, I do not pull the permeants. For each position,
>> I add it to the system in VMD, and then run the
>> minimization-equilibration-production sequence.
>> >
>> > Since GROMACS was doing the nstlist update automatically, I thought it
>> was
>> >> a "safe" option. Based on the above,however, I believe that since the
>> >> simulation runs with constraints, the system is more sensitive and thus
>> the
>> >> nstlist should be smaller. Can anyone validate that my hypothesis is
>> >> correct?
>> >>
>> > On the limited evidence available, I think your observations are more
>> > likely to be symptoms of the problem than the problem itself. How does
>> e.g.
>> > one permeant molecule behave?
>> Once again sorry for not being clear. The plural in permeants comes from
>> the fact that I test several molecules (e.g. water, ammonia etc) but it
>> is only one per simulation/system.
>>
>> So, you think that the updating of neighbour list, shouldn't affect the
>> problem?
>>
>
> Yep. You can force mdrun to use an nstlist of your choosing with mdrun
> -nstlist n. I would expect you to observe similar instability in such
> cases. The only case in which I think nstlist could matter is when
> something else is so badly wrong that particles diffuse further than the
> automated buffer anticipates, in which case any reduction of nstlist that
> keeps the simulation running is merely deferring the problem to analysis
> time...

Good point. However, I think that could be the case only if the error
in the estimate would cause several orders of magnitude increase in
drift. The larger nstlist, the more the Verlet buffer calculation
tends to over-estimate the buffer (i.e. err on the safe side). Hence,
the missed interactions would need to have a very large drift
contribution to be able to cancel the effect of the buffer
estimation's.

@Michail: have you looked at energy conservation? I suggest to try
multiple nstlist values and plot the energy drift. If what you see is
really related to a _direct_ effect of the search frequency, that
should manifest in a buffer under-estimate and consequently high
non-bonded drift contribution which is relatively easy to diagnose.

--
Szilárd

>
> Mark
>
>
>> Kind Regards,
>> Michail
>>
>> > Mark
>> >
>> > Thanks in advance.
>> >>
>> >> Kind Regards,
>> >>
>> >> Michail
>> >>
>> >>
>> >> -------------------------------------------------------------------
>> >> Michail (Michalis) Palaiokostas
>> >> PhD Student
>> >> School of Engineering and Materials Science
>> >> Queen Mary University of London
>> >> -------------------------------------------------------------------
>> >> --
>> >> 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.