1) make random box of 128 DSPC molecules from single DSPC molecule:
 genbox -ci ~/mdtutorial/lipids/dspc_single.gro -nmol 128 -box 7 7 7 -try 100 -o 128_noW.gro 

2) minimise 
grompp -f em -c 128_noW.gro -p dspc.top -maxwarn 10
mdrun -v -c 128_minimised.gro

3) add 6 CG waters (or 24 all-atom waters) per lipid
 genbox -cp 128_minimised.gro -cs ~/mdtutorial/water.gro -o waterbox -maxsol 768

4) minimise again
grompp -f em.mdp -c waterbox.gro -p ../dspc.top -maxwarn 10
mdrun -v -c minimised

5) run MD for ca. 25 ns
cd md
grompp -f ../md.mdp -c ../minimised.gro -p ../../dspc.top -maxwarn 10
mdrun (ca 30 min on 2 CPUs)

6) check if you got a bilayer. If yes, check if the formed membrane is
normal to the z-axis (i.e., membrane in the xy-plane). If the latter
is not the case, rotate the system (with editconf). In case you did
not get a bilayer at all, continue with the pre-formed one that you find
here: dspc_bilayer.gro

7) continue simulation for another 10 ns at zero surface tension
(semiisotropic pressure coupling)

8) From this simulation, calculate properties such as area per lipid,
lateral diffusion of the lipids, and bilayer thickness. For the first,
you can simply divide the area of your simulation box (Box-X and -Y from
g_energy) by half the number of lipids in your system. For the latter
two, use the gromacs tools g_msd and g_density, respectively. Before
calculating the lateral diffusion, remove jumps over the box
boundaries (trjconv -pbc nojump). Then, when using g_msd, take care to
remove 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.

In general, for the analysis, you might want to discard the first few
ns of your simulation as equilibration time.

Compare your results to those from small-angle neutron
scattering (SANS) experiments (Balgavy et al., Biochim. et
Biophys. Acta 2001, 1512, 40-52): 
thickness = 4.98  0.15 nm; 
area per lipid = 0.65  0.05 nm^2

9) Next, introduce a double bond in each lipid tail and investigate
the effect on the bilayer properties. To do so, you can replace the
DSPC lipids by DOPC (compare the MARTINI topologies of these two
lipids). A simple way is to open the confout.gro from the last DSPC
simulation with your favorite text editor and to search-replace all
"DSPC" by "DOPC". Replace DSPC by DOPC in your .top and .mdp files
accordingly, and grompp does the rest for you (you can ignore the
"atom name does not match" warnings of grompp).

10) run another 10-ns MD simulation and compare the above properties for the
two lipid types (saturatated DSPC and unsaturated DOPC). How did the properties
change? Do the observed changes match your expectations? Why/why not?

11) Now, starting again from the final snapshot from the DSPC
simulation, change the head groups from PC to PE. Again, an easy way
to do so is to search-replace DSPC by DSPE in the .gro, .top, and .mdp
files, and run grompp. Also here, run an MD simulation and compare the
above properties of DSPC against DSPE. Discuss.
