Hi Everyone,<div><br></div><div>I did another simulation where i found the the DPPC area per lipid is 0.61 nm^2 . Is this acceptable ?</div><div><br></div><div>I have seen issues like this on the mailing list before can one of the experts give me some hints.</div>
<div><br></div><div>Amit<br><br><div class="gmail_quote">On Wed, Nov 9, 2011 at 12:24 PM, Amit Choubey <span dir="ltr"><<a href="mailto:kgp.amit@gmail.com">kgp.amit@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex;">
<div><br></div>5Hello all,<div><br></div><div>I am trying to use CHARMM 36 for DPPC membrane simulation. I did the following so far:</div><div><br></div><div>1. Download pdb file containing 128 DPPC molecules from <a href="http://www.charmm-gui.org/?doc=archive&lib=lipid_pure" target="_blank">http://www.charmm-gui.org/?doc=archive&lib=lipid_pure</a></div>
<div><br></div><div>2. I separated one lipid molecule from the obtained pdb file and used </div><div>pdb2gmx -f 1dppc.gro -nochargegrp </div><div><br></div><div>The top file obtained was converted to itp file by commenting few things out.</div>
<div><br></div><div>3. I looked at a file named dppc_n128.inp in the downloaded dppc bilayer from CHARMM GUI for the box length. I used editconf to convert the pdb to gro and manually inserted the box size in the gro file. </div>
<div><br></div><div>4. Created a system.top file which looks like</div><div><div>#include "charmm36.ff/forcefield.itp"</div><div>#include "lipid.itp"</div><div>#include "charmm36.ff/tips3p.itp"</div>
<div><br></div><div>[ system ]</div><div>chol</div><div><br></div><div>[ molecules ]</div><div>DPPC 128</div><div>SOL 3659</div></div><div><br></div><div>5. I used following mdp settings</div><div><br></div><div><div>integrator<span style="white-space:pre-wrap"> </span>= md<span style="white-space:pre-wrap"> </span>; leap-frog integrator</div>
<div>nsteps<span style="white-space:pre-wrap"> </span>= 20000000<span style="white-space:pre-wrap"> </span>; 40 ns</div><div>dt<span style="white-space:pre-wrap"> </span>= 0.002<span style="white-space:pre-wrap"> </span>; 2 fs</div>
<div><br></div><div>; Output control</div><div>nstxout<span style="white-space:pre-wrap"> </span>= 100000<span style="white-space:pre-wrap"> </span>; save coordinates every 200 ps</div>
<div>nstvout<span style="white-space:pre-wrap"> </span>= 100000<span style="white-space:pre-wrap"> </span>; save velocities every 200 ps</div><div>nstenergy<span style="white-space:pre-wrap"> </span>= 10000<span style="white-space:pre-wrap"> </span>; save energies every 2 ps</div>
<div>nstlog<span style="white-space:pre-wrap"> </span>= 10000<span style="white-space:pre-wrap"> </span>; update log file every 2 ps</div><div><br></div><div>; Neighborsearching</div>
<div>ns_type<span style="white-space:pre-wrap"> </span>= grid<span style="white-space:pre-wrap"> </span>; search neighboring grid cels</div><div>nstlist<span style="white-space:pre-wrap"> </span>= 10<span style="white-space:pre-wrap"> </span>; 10 fs</div>
<div>rlist<span style="white-space:pre-wrap"> </span>= 1.0<span style="white-space:pre-wrap"> </span>; short-range neighborlist cutoff (in nm)</div><div>rlistlong <span style="white-space:pre-wrap"> </span>= 1.4</div>
<div>vdwtype<span style="white-space:pre-wrap"> </span>= switch</div><div>rvdw<span style="white-space:pre-wrap"> </span>= 1.2<span style="white-space:pre-wrap"> </span>; short-range van der Waals cutoff (in nm)</div>
<div>rvdw_switch<span style="white-space:pre-wrap"> </span>= 0.8</div><div><br></div><div>; Electrostatics</div><div>rcoulomb = 1.0</div><div>rcoulomb_switch<span style="white-space:pre-wrap"> </span>= 0.0</div>
<div>coulombtype<span style="white-space:pre-wrap"> </span>= PME<span style="white-space:pre-wrap"> </span>; Particle Mesh Ewald for long-range electrostatics</div><div>pme_order<span style="white-space:pre-wrap"> </span>= 6<span style="white-space:pre-wrap"> </span>; cubic interpolation</div>
<div>fourierspacing<span style="white-space:pre-wrap"> </span>= 0.15<span style="white-space:pre-wrap"> </span>; grid spacing for FFT</div><div><br></div><div>; Temperature coupling is on</div>
<div>tcoupl<span style="white-space:pre-wrap"> </span>= Nose-Hoover<span style="white-space:pre-wrap"> </span>; More accurate thermostat</div><div>tc-grps<span style="white-space:pre-wrap"> </span>= DPPC<span style="white-space:pre-wrap"> </span>SOL ; two coupling groups - more accurate</div>
<div>tau_t<span style="white-space:pre-wrap"> </span>= 0.2<span style="white-space:pre-wrap"> </span>0.2 ; time constant, in ps</div><div>ref_t<span style="white-space:pre-wrap"> </span>= 323.15 <span style="white-space:pre-wrap"> </span>323.15 ; reference temperature, one for each group, in K</div>
<div><br></div><div>; Pressure coupling is on</div><div>pcoupl<span style="white-space:pre-wrap"> </span>= Parrinello-Rahman<span style="white-space:pre-wrap"> </span>; Pressure coupling on in NPT</div>
<div>pcoupltype<span style="white-space:pre-wrap"> </span>= semiisotropic<span style="white-space:pre-wrap"> </span>; uniform scaling of x-y box vectors, independent z</div><div>tau_p<span style="white-space:pre-wrap"> </span>= 5.0<span style="white-space:pre-wrap"> </span>; time constant, in ps</div>
<div>ref_p<span style="white-space:pre-wrap"> </span>= 1.0<span style="white-space:pre-wrap"> </span>1.0<span style="white-space:pre-wrap"> </span>; reference pressure, x-y, z (in bar)</div>
<div>compressibility = 4.5e-5 4.5e-5<span style="white-space:pre-wrap"> </span>; isothermal compressibility, bar^-1</div><div><br></div><div>; Periodic boundary conditions</div><div>pbc<span style="white-space:pre-wrap"> </span>= xyz<span style="white-space:pre-wrap"> </span>; 3-D PBC</div>
<div><br></div><div>; Dispersion correction</div><div>DispCorr<span style="white-space:pre-wrap"> </span>= No<span style="white-space:pre-wrap"> </span>; account for cut-off vdW scheme</div>
<div><br></div><div><br></div><div>; Velocity generation</div><div>gen_vel = yes</div><div>gen_temp = 200.0</div><div>gen_seed = 173529</div><div><br></div><div>; COM motion removal</div>
<div>nstcomm = 1</div><div>comm-mode = Linear</div><div>comm-grps = DPPC SOL</div><div><br></div><div>; Energy monitoring</div><div>energygrps = DPPC SOL</div><div><br></div><div>; Bond parameters</div>
<div>constraint_algorithm = lincs ; holonomic constraints</div><div>constraints = hbonds ; all bonds (even heavy atom-H bonds) constrained</div><div>lincs_iter = 1 ; accuracy of LINCS</div><div>
lincs_order = 4 ; also related to accuracy</div></div><div><br></div><div>The resulting simulation gave me an area per lipid =0.591 nm^2. This is not quite right. Can somebody suggest what else could i do it get close to the 0.63 nm^2 value. </div>
<div><br></div><div>Also when i looked at trajectory file it seemed that the box length was not set up quite right because the downloaded pdb file looked ok with the membrane in the middle but at the next frame the membrane was at a totally different location. The leaflets were separated and the middle of the membrane was at the edge of the box.</div>
<div><br></div><div>Also when i used pdb2gmx -nochargegrp on the downloaded file it created a gro file whose system size was </div><div> 8.20170 8.64120 6.54740 in contrast to 6.40000 6.40000 6.38452 . </div>
<span class="HOEnZb"><font color="#888888">
<div><br></div><div>Amit</div>
</font></span></blockquote></div><br></div>