Here is the tutorial I followed, just to verify everything works on my computer.the purpose of this tutorial is to study the self-assembling of a dipalmitoyl-phosphatidylcholine (DPPC) bilayer from a random configuration of lipids and water in the simulation box.

xeon:MD/bilayer-lipidome-tutorial/spontaneous-assembly/initial_assembly

dppc_single.gro:  coordinates for a single DPPC molecule.

water.gro martini model of water box.

we start with a box of randomly filled lipids and waters.

gmx_mpi insert-molecules -ci dppc_single.gro -box 7.5 7.5 7.5 -nmol 128 -try 500 -o 128_noW.gro

Perform a short energy minimization of the system containing only the lipids:

gmx_mpi grompp -f minimization.mdp -c 128_noW.gro -p dppc.top -o dppc-min-init.tpr

gmx_mpi mdrun -deffnm dppc-min-init -v -c 128_minimized.gro

now we solvate the box with CG waters.

gmx_mpi solvate -cp 128_minimized.gro -cs water.gro -o waterbox.gro -maxsol 768 -radius 0.21

now minimize the solvated box.

gmx_mpi grompp -f minimization.mdp -c waterbox.gro -p dppc.top -o dppc-min-solvent.tpr

gmx_mpi mdrun -deffnm dppc-min-solvent -v -c minimized.gro

now let's run a 30ns md simulation

gmx_mpi grompp -f martini_md.mdp -c minimized.gro -p dppc.top -o dppc-md.tpr

gmx_mpi mdrun -deffnm dppc-md -v

You may want to check the progress of the simulation to see whether the bilayer has already formed before the end of the simulation by gmx_mpi view,

gmx_mpi view -f dppc-md.xtc -s dppc-md.tpr

or you can view the trajectory by vmd

./do_vmd.sh(or vmd dppc-md-start.gro dppc-md-viz.xtc, then source cg_bonds.tcl, cg_bonds -top dppc.top, if the script doesn't work )

 

Bilayer thickness:

The bilayer thickness can be obtain from the distance between the headgroup peaks(NH3) in the density profile, before we get density profile, we need to rotate the box and make the bilayer in X-Y plane.

make_ndx_mpi -f dppc-md.gro -o index.ndx

>a N*

>q

editconf_mpi -f dppc-md.gro -rotate 50 0 0 -o dppc_aligned.gro

g_density_mpi -f  dppc_aligned.gro -n index.ndx -s dppc-md.tpr -o head_density.xvg(choose group N*)

From the profile, we estimate the thickness of bilayer approximately equal to 13.8-9.6=4.2nm .

Lateral diffusion:

Before calculating the lateral diffusion, remove jumps over the box boundaries (trjconv -pbc nojump). Then, calculate the lateral diffusion using g_msd. Take care to remove the overall center of mass motion (-rmcomm*), and to fit the line only to the linear regime of the mean-square-displacement curve (-beginfit and -endfit options of g_msd). To get the lateral diffusion, choose the -lateral z option.

trjconv_mpi -pbc nojump -f dppc-md.xtc -o dppc_nojump.xtc

g_msd_mpi -f dppc_nojump.xtc -s dppc-md.tpr -lateral z 

Order parameters

(second-rank) order parameter is defined as:

P₂ = 1/2 (3 cos²<θ> − 1),

where θ is the angle between the bond and the bilayer normal. P₂ = 1 means perfect alignment with the bilayer normal, P₂ = −0.5 anti-alignment, and P₂ = 0 random orientation.

do-order.py dppc-md.xtc dppc-md.tpr 0 300000 100 0 0 1 1280 DPPC

will for example read 0 to 300 ns from the trajectory dppc-md.xtc. The simulation has 1280 DPPC lipids and the output is the P2, calculated relative to the Z-axis, and average over results from every 100th equally-spaced snapshot in the trajectory.

Abs average order parameters for carbon chains <Sn> = 0.054

Here is another tutorial of DPPC by Justin Lemkul.