Peter C. Lai wrote:
Hi
I haven't played with gromacs in years, so I'm pretty new at the ability
to play with large systems like lipid bilayers.
What is the "correct" or lets say canonical method to equilibrate a
membrane patch of say 9x9 nm POPC with sufficient vertical waters to take
care of pbc and long range interactions (> 1.5nm above/below the bilayer)
using Tom's Charmm36 FF port.
I am seeing a lot of different methods out there to: construct and
equilibrate such a system in prepration for g_embed or whatever.
I see from
http://www.mail-archive.com/gmx-users@gromacs.org/msg33812.html
the ability to use genbox to replicate Tieleman's 128 patch
but I also can build my bilayer inside VMD, which supports charmm36->pdb
atom names (and things seem to work fine)
I have yet to see anybody's run parameters and constraints for running
the equilibration, particularly using charmm36 (using vdwtype=switch).
Should I constrain all heavy atoms in the lipid during the NVT run, or
just the polar group (or even just P) and allow the tails to "melt" as seen
elsewhere? Should I only set Z axis constraints? How long should I run in
NVT before switching to NPT and for how long (10ns?) And should I change the
constraints moving from NVT to NPT?
As a point of clarity, "constraints" and "restraints" serve very different
functions in Gromacs:
http://www.gromacs.org/Documentation/Terminology/Constraints_and_Restraints
I presume you mean "restraints" in this context. I see no point in restraining
the lipids completely during equilibration. It sort of defeats the purpose.
Depending on how reasonable the starting configuration is, restraints on any
group may not be necessary at all. Poor starting configurations, particularly
with inadequate solvation within the interfacial region, often require vertical
position restraints on P atoms to prevent the membrane from separating.
Restraining the entire headgroup probably restricts lipid motion too much.
There is no standard timeframe for equilibration. Rotational relaxation occurs
within a very short amount of time, translational within 10-20 ns for most
lipids. You should pay close attention to diffusion constants, membrane
thickness, area per lipid, etc to judge the quality of equilibration.
How should the thermostat be coupled? System or separate the water from the
lipid?
This consideration is uniform among just about all systems and is a function of
the thermostat rather than the system itself. Couple water and lipids
separately, just as you would couple protein and solvent separately for most
thermostats.
COM motion should be removed separately for water and lipids as well, to prevent
lateral sliding that is hidden through global COM motion removal.
I experienced the same thing as Ng Hui Wen back in October with waters
being extremely dehomogenized (and diffusing out of the box) after 100ps
NVT with all heavy POPC constrained at 1000 kJ/nm for all axes. I am sure
that moving to NPT will restore the waters but still it was a bit scary to
see for the first run.
That is probably an artifact of restraining all lipid atoms and not allowing for
proper response between water and lipids.
If a lot of these approaches haven't been published, maybe we can just
stick them in the wiki to let everyone "pick their poison" when it comes
to equilibrating their bilayers.
I can think of at least a dozen plausible methods for equilibrating a lipid
bilayer system, making it probably more confusing to simply list them and tell
users to pick, especially if they don't know what they're doing. As with any
equilibration of any system, you have to judge (1) whether the desired
thermodynamic ensemble is achieved and (2) whether the properties of the system
of interest are established. With these general rules in mind, the user is free
to design their own procedure.
I wrote a membrane protein tutorial a long time ago that outlines many of these
principles and discusses several more.
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/membrane_protein/index.html
Even though your purposes may not involve a protein or the same force field, you
may find some of the general information useful.
-Justin
Thanks for any ideas/hints/suggestions
--
========================================
Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin
========================================
--
gmx-users mailing list gmx-users@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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
