On 2/23/14, 8:30 AM, sujith wrote:
Dear GROMACS users, I am new to GROMACS, and recently started using the version 4.6.5. I have seen a lot of NPT related issues raised earlier in this forum, but in my case the error looks much more severe. I am following Justin Lemkul's tutorial (http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/biphasic/01_genconf.html) on Biphasic system, containing cyclohexane and water.Everything went well till the NPT simulation to bring system to reference pressure of 1 bar. Here are the details of the system , .mdp file , what I did and where I find the problem. SYSTEM : 512 cyclohexane + 2720 water molecules. CURRENT STATUS: (1) Energy minimization : energy converged. (2) NVT (tcoupl=berendsen , tau_t=0.1ps / ref_t=298K) completed and target pressure of 298 K attained. (3) NPT (tcoupl=nose-hoover, tau_t=0.5ps / pcoupl=berendsen, tau_p=1ps / ref_p = 1bar , total time=15 ns) completed with the following result: Energy Average Err.Est. RMSD Tot-Drift ------------------------------------------------------------------------------- Pressure 1.08924 0.67 114.66 -1.74033 (bar) Since berendsen barostat doesn't generate the true ensemble, an equilibration with pcoupl=parrinello-rahman is performed in the next step. (4) NPT (pcoupl=parrinello-rahman, tau_p=5ps, total time=5ns) ----> PROBLEM FACED: The average pressure too high as shown below which I feel is not going to improve. Energy Average Err.Est. RMSD Tot-Drift ------------------------------------------------------------------------------- Pressure 23.6651 0.53 144.464 -1.00634 (bar) I am aware of the fact that pressure is subject to large fluctuations in small sized systems and that this may affect the average value of the pressure. But , here the average pressure looks too large to be ignored. The pressure-vs-time graph doesn't show any upward trend, and the pressure looks like fluctuating about a central value. Here is the .mdp file. Only changes made from the previous .mdp for berendsen pressure coupling , are with pcoupl and tau_p. ; 7.3.3 Run Control integrator = md tinit = 0 dt = 0.002 nsteps = 2500000 comm_mode = Linear nstcomm = 1 comm_grps = CHX SOL ; 7.3.8 Output Control nstxout = 2500 nstvout = 2500 nstfout = 2500 nstlog = 2500 nstenergy = 100 nstxtcout = 1000 xtc_precision = 1000 xtc_grps = System energygrps = System ; 7.3.9 Neighbor Searching nstlist = 1 ns_type = grid pbc = xyz rlist = 0.8 ; 7.3.10 Electrostatics coulombtype = PME rcoulomb = 0.8 ; 7.3.11 VdW vdwtype = cut-off rvdw = 0.8
What force field are you using? If it's Gromos96 like my tutorial, the value of rvdw is much too short and can lead to artifacts. Using 0.8 for rlist/rcoulomb is fine, though 0.9 is more common. rvdw should be 1.4.
DispCorr = EnerPres ; 7.3.13 Ewald fourierspacing = 0.12 pme_order = 4 ewald_rtol = 1e-5 ; 7.3.14 Temperature Coupling tcoupl = nose-hoover tc_grps = CHX SOL tau_t = 0.5 0.5 ref_t = 298 298 ; 7.3.15 Pressure Coupling pcoupl = parrinello-rahman pcoupltype = isotropic tau_p = 5.0 compressibility = 4.5e-5 ref_p = 1.0
The value of tau_p seems a bit long to me; does changing it to 1.0 improve the results?
; 7.3.17 Velocity Generation gen_vel = no ; 7.3.18 Bonds constraints = all-bonds constraint_algorithm = LINCS continuation = yes lincs_order = 4 lincs_iter = 1 lincs_warnangle = 30 "npt.mdp" 63L, 4080C I guess there is something seriously wrong in the choice of methods/parameters in the .mdp file, which I cant figure out. Kindly go through and let me know your comments. I would be happy to give any further details. Any help would be appreciated.
An even better test is to simplify the system. Run pure water or pure cyclohexane with the existing settings, then try the modifications above. That should help root out whether you're having a problem with the .mdp settings, force field, or maybe even the combination of both.
-Justin -- ================================================== Justin A. Lemkul, Ph.D. Ruth L. Kirschstein NRSA Postdoctoral Fellow Department of Pharmaceutical Sciences School of Pharmacy Health Sciences Facility II, Room 601 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.