Re: popc pressure profile

From: Marcos Sotomayor (sotomayo_at_ks.uiuc.edu)
Date: Tue Oct 10 2006 - 12:10:29 CDT

The membrane plugin of VMD provides membrane templates that are a bit too
ordered. I would suggest you do what Justin did during equilibration:

"After this initial construction stage, the headgroup atoms were fixed and
the rest of the atoms heated to 300 K for 500 ps to "melt" the aliphatic
carbon tails. Fixing the headgroups prevented oppositely charged groups in
neighboring lipids from bonding with each other, a condition that would
prevent proper solvation of the interface. "
  Gullingsrud/Schulten (Biophysical Journal 86:3496-3509 (2004))

Fixing the headgroups during part of the equilibration may help you get a
more disordered (real) membrane.

In addition, note that all the averages presented in Justin's paper are
performed over simulation with at least 1.5ns of equilibration and then
10~ns of dynamics. Also, how many lipids do you have in your membrane? 60?
Seems a bit small...

Marcos
PS: Why do you have PmeGridsizeZ 128 ? Is your system that large in the Z
direction? You may also want to use some of the NAMD features to eliminate
center of mass motion.

On Tue, 10 Oct 2006, Allan Haldane wrote:

> Hi,
>
> I am quite new to MD, so I thought I would try to reproduce the result
> Gullingsrud/Schulten (Biophysical Journal 86:3496-3509 (2004)) found for the
> lateral pressure profile of a POPC membrane. It doesn't look right, though. I
> have ideas of what might be wrong, but maybe someone can quickly put me in
> the right direction.
>
> Here is a plot of the pressure tensor components I get averaged over 5 ns
> (there is a bit of the water regions removed from the plot):
> http://www.ugrad.physics.mcgill.ca/~ahalda/upload/protein/NAMD/press26b1.png
>
> Here are the bonded + nonbonded components: (compare to simulation D4 in
> figure 6 of the paper)
> http://www.ugrad.physics.mcgill.ca/~ahalda/upload/protein/NAMD/profComponents.png
>
> I start with the popc membrane template (popc_box) that comes with vmd, and
> after adding some extra water (10 A above and below) I heat over 10 ps (by
> langevin dynamics) and then equilibrate it for 1 ns at NPT. N, P, T and V are
> indeed constant for most of this. Actually I calculate that popc_box has an
> area per lipid of 76 A^2 By the end of the equilibration which is high. (I
> calculate area/lipid as cell area/lipids in monolayer (30). Is it truly that
> simple?) I see one of Gullingsrud/Schulten's POPC runs was at around this
> area/lipid, though (simulation D4).
>
> I then do a pressure profile for 5 ns at NVT with nonbonded part done
> seperately (with NAMD 2.6.b1). (I also tried with 2.6.2 + Ewald and got quite
> different results...)
>
> I see two things happen during the pressure profile run: First, the global
> pressure increases by 53 bar over the 5 ns. It is initially -50 bar (unlike
> at the end of the NPT equilibration), but from a post on this list this is
> apparently normal for NVT simulations with liquid. Second, the bilayer moves
> quite a lot in the z (normal) direction, by a total shift of 4 A. (For the
> NAMD 2.6.2 run I did, it shifted by 10 A) This happens gradually and mostly
> continuously through the run.
>
> Here is my config file for the profile calculation. The equilibration uses
> the same parameters except that: langevinPiston on, langevinDamping 10,
> pressureProfile off and I output more frequently. My system is 50 by 46 by 69
> A. I do nonbonded seperately, as shown in the manual.
>
> Thanks!
> Allan
>
> -----------------------
>
> structure start.psf
> coordinates requil.pdb
> parameters par_all27_prot_lipid.inp
> parameters par-extraterms.inp
> paraTypeCharmm on
>
> extendedSystem requil.xsc
> velocities requil.vel
>
> outputEnergies 1000
> outputTiming 1000
> xstFreq 100000
> dcdFreq 10000
> wrapAll on
>
> restartname output/restart
> restartfreq 100000
>
> timestep 1
> nonBondedFreq 2
> fullElectFrequency 4
> stepsPerCycle 20
>
> switching on
> switchDist 8.5
> cutoff 10
> pairlistdist 11.5
>
> Pme on
> PmeGridsizeX 64
> PmeGridsizeY 64
> PmeGridsizeZ 128
>
> exclude scaled1-4
> 1-4scaling 1.0
>
> langevin on
> langevinDamping 1
> langevinTemp 310
> langevinHydrogen off
>
> langevinPiston off
> langevinPistonTarget 1.01325
> langevinPistonPeriod 200
> langevinPistonDecay 100
> langevinPistonTemp 310
>
> #these two only should affect equilibration:
> useGroupPressure yes #hydrogen excluded
> useFlexibleCell yes
>
> binaryoutput off
> outputname output/equil_out
>
> pressureProfile on
> pressureProfileSlabs 120
> pressureProfileFreq 10000
>
> run 5000000
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:44:03 CST