On 2/24/16 7:25 AM, Gmx QA wrote:
Dear list I have a simulation of a system of 64 POPC lipids per leaflet using the Charmm36 lipid parameter set, and I am calculating the electron density across the bilayer to compare to literature values. Unfortunately I cannot quite get the same results as eg. in Klauda et al ( http://www.sciencedirect.com/science/article/pii/S0005273614002193) I have made the starting structure using the charmm-gui, and gone through all the equilibration steps as suggested in the gromacs tar-ball from that website. My production mdp file is the following: integrator = md dt = 0.002 nsteps = 50000000 nstlog = 10000 nstxout = 500000 nstvout = 500000 nstfout = 0 nstcalcenergy = 10 nstenergy = 1000 nstxtcout = 500000 ; cutoff-scheme = Verlet nstlist = 20 rlist = 1.2 coulombtype = pme rcoulomb = 1.2 vdwtype = Cut-off vdw-modifier = Force-switch rvdw_switch = 1.0 rvdw = 1.2 ; tcoupl = Nose-Hoover tc_grps = POPC SOL tau_t = 1.0 1.0 ref_t = 303.15 303.15 ;nsttcouple = 1 ; pcoupl = Parrinello-Rahman pcoupltype = semiisotropic tau_p = 5.0 5.0 compressibility = 4.5e-5 4.5e-5 ref_p = 1.0 1.0 ; constraints = h-bonds constraint_algorithm = LINCS continuation = yes ; nstcomm = 10 comm_mode = linear comm_grps = POPC SOL ; refcoord_scaling = com Simulations are run with gromacs 5.0.4, and currently 100 ns long. Before running vmx density, I center my trajectory like this: $ gmx trjconv -f run.trr -s em.tpr -o run_center.trr -center selecting POPC as the group to center on. Then, my command for computing the electron density is the following (also using gmx density 5.0.4): $ gmx density -f run_center.trr -s em.tpr -o density_traj_center.xvg -dens electron -ei electrons.dat -center Again, selecting POPC to center the density profile on. The resulting profile can be found here (i hope that works) http://www.filedropper.com/compareelectrondensity As you can see, the depth of the central trough is not as deep as the literature value, nor is the height of the density the same even though the peak-peak distance appears to be similar, and values in the other regions also match _fairly_ closely. The average APL value over the entire trajectory is 64.0 +- 1.33, which is in line with the reported 64.7 from Klauda, but the electron density worries me a bit. Perhaps anyone (Justin?) can shed any light on this?
Why? The results match almost perfectly. The magnitude of error is within a few % and if your APL matches then I'm sure everything is fine. We rigorously tested the lipid inputs, cutoff setup, etc. and everything agrees quite well between lipid simulations in CHARMM, GROMACS, AMBER, NAMD, and OpenMM. There are minor variations, as would be expected in comparing any two software programs, but nothing is really aberrant.
http://dx.doi.org/10.1021/acs.jctc.5b00935 -Justin -- ================================================== Justin A. Lemkul, Ph.D. Ruth L. Kirschstein NRSA Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 629 University of Maryland, Baltimore 20 Penn St. Baltimore, MD 21201 jalem...@outerbanks.umaryland.edu | (410) 706-7441 http://mackerell.umaryland.edu/~jalemkul ================================================== -- Gromacs Users mailing list * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before posting! * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists * For (un)subscribe requests visit https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or send a mail to gmx-users-requ...@gromacs.org.
