Dear Justin, Thanks a lot for pointing out the issues. I now understand why there were such high oscillations.
Could you please also tell me if there are any ideal values for pme_order and fourier spacing with respect to the cut offs' value of 1.4. Does the following Note imply I can raise the fourier grid spacing to 0.25? NOTE 2 [file sim-new.mdp]: The optimal PME mesh load for parallel simulations is below 0.5 and for highly parallel simulations between 0.25 and 0.33, for higher performance, increase the cut-off and the PME grid spacing Thank you again, Nisha On Thu, Jun 30, 2016 at 6:55 PM, Justin Lemkul <jalem...@vt.edu> wrote: > > > On 6/30/16 9:16 AM, NISHA Prakash wrote: > >> Dear Justin, >> >> Thank you for your reply. >> It is a protein carbohydrate system. Including the solvent, the number of >> atoms is 43499. >> I have minimized the system for 200 ps followed by NPT and NVT simulations >> for 200 ps respectively >> >> > Given that your temperature output started from 0 K, then you did not > continue from the equilibration properly by supplying the checkpoint file > to grompp -t. This is important to get right, otherwise you're basically > starting over from some random point (an equilibrated structure without any > velocities likely isn't a physically realistic state). > > Below is the .mdp file. >> >> >> ; VARIOUS PREPROCESSING OPTIONS >> title = REMD Simulation >> define = -DPOSRES >> >> >> ; RUN CONTROL PARAMETERS >> integrator = md-vv ; velocity verlet algorithm - >> tinit = 0 ; >> dt = 0.002 ; timestep in ps >> nsteps = 5000000 ; >> simulation-part = 1 ; Part index is updated automatically on >> checkpointing >> comm-mode = Linear ; mode for center of mass motion removal >> nstcomm = 100 ; number of steps for center of mass motion >> removal >> comm-grps = Protein_Carb Water_and_Ions ; group(s) for >> center of mass motion removal >> >> > In a solvated system, you should not be separating these groups. This > could explain the sudden jump in temperature - you could have things > clashing badly over the course of the simulation. > > > >> ; ENERGY MINIMIZATION OPTIONS >> emtol = 10 ; Force tolerance >> emstep = 0.01 ; initial step-size >> niter = 20 ; Max number of iterations in relax-shells >> fcstep = 0 ; Step size (ps^2) for minimization of >> flexible constraints >> nstcgsteep = 1000 ; Frequency of steepest descents steps >> when >> doing CG >> nbfgscorr = 10 >> >> >> ; OUTPUT CONTROL OPTIONS >> nstxout = 50000 ; Writing full precision coordinates >> every >> ns >> nstvout = 50000 ; Writing velocities every nanosecond >> nstfout = 0 ; Not writing forces >> nstlog = 5000 ; Writing to the log file every step 10ps >> nstcalcenergy = 100 >> nstenergy = 5000 ; Writing out energy information every >> step 10ps >> nstxtcout = 2500 ; Writing coordinates every 5 ps >> xtc-precision = 1000 >> xtc-grps = Protein_Carb Water_and_Ions ; subset of >> atoms for the .xtc file. >> energygrps = Protein_Carb Water_and_Ions ; Selection of >> energy groups >> >> >> ; NEIGHBORSEARCHING PARAMETERS >> nstlist = 10 ; nblist update frequency- >> ns-type = Grid ; ns algorithm (simple or grid) >> pbc = xyz ; Periodic boundary conditions: xyz, >> no, >> xy >> periodic-molecules = no >> rlist = 1.4 ; nblist cut-off >> rlistlong = -1 ; long-range cut-off for switched >> potentials >> >> >> ; OPTIONS FOR ELECTROSTATICS >> coulombtype = PME ; Method for doing electrostatics >> rcoulomb = 1.4 ; >> epsilon-r = 1 ; Relative dielectric constant for the >> medium >> pme_order = 10; >> >> >> ; OPTIONS FOR VDW >> vdw-type = Cut-off ; Method for doing Van der Waals >> rvdw-switch = 0 ; cut-off lengths >> rvdw = 1.4 ; >> DispCorr = EnerPres; Apply long range dispersion >> corrections for Energy and Pressure >> table-extension = 1 ; Extension of the potential lookup tables >> beyond the cut-off >> fourierspacing = 0.08 ; Spacing for the PME/PPPM FFT grid >> >> > This small Fourier spacing, coupled with the very high PME order above, is > going to unnecessarily slow your system down. Is there some reason you > have set these this way? > > >> ; GENERALIZED BORN ELECTROSTATICS >> gb-algorithm = Still ; Algorithm for calculating Born radii >> nstgbradii = 1 ; Frequency of calculating the Born radii >> inside rlist >> rgbradii = 1 ; Cutoff for Born radii calculation >> gb-epsilon-solvent = 80 ; Dielectric coefficient of the implicit >> solvent >> gb-saltconc = 0 ; Salt concentration in M for Generalized >> Born models >> >> >> ; Scaling factors used in the OBC GB model. Default values are OBC(II) >> gb-obc-alpha = 1 >> gb-obc-beta = 0.8 >> gb-obc-gamma = 4.85 >> gb-dielectric-offset = 0.009 >> sa-algorithm = Ace-approximation >> sa-surface-tension = -1 ; Surface tension (kJ/mol/nm^2) for the >> SA >> (nonpolar surface) part of GBSA - default -1 >> >> > Implicit solvent should not be used if you have explicit solvent, though > it looks like these options are probably off since the default for the > implicit-solvent keyword is "no," but be aware that these are extraneous. > > >> >> ; Temperature coupling >> tcoupl = nose-hoover >> nsttcouple = 10 ; >> nh-chain-length = 10 >> tc-grps = Protein_Carb Water_and_Ions ; Groups to >> couple separately >> tau-t = 10 10; Time constant (ps)- >> > > tau-t = 10 is far too permissive and will allow the temperature to > oscillate very widely. A value of 1 or so should generally be used with > N-H. > > -Justin > > > ref-t = 270.0 270.0; reference temperature (K) >> >> >> ; pressure coupling >> pcoupl = no ;- >> >> >> ; GENERATE VELOCITIES FOR STARTUP RUN >> gen-vel = no >> gen-temp = 270.0 >> gen-seed = 173529 >> >> >> ; OPTIONS FOR BONDS >> continuation = yes ; constrain the start configuration >> >> constraints = all-bonds >> constraint-algorithm = lincs ; Type of constraint algorithm- >> lincs-order = 4 >> lincs-iter = 1 >> lincs-warnangle = 30 >> >> >> Thank you for your help. >> >> Nisha >> >> >> >> On Thu, Jun 30, 2016 at 6:21 PM, Justin Lemkul <jalem...@vt.edu> wrote: >> >> >>> >>> On 6/30/16 8:46 AM, NISHA Prakash wrote: >>> >>> Dear all, >>>> >>>> I have conducted a 10ns REMD simulation for a protein ligand complex >>>> with >>>> the temperature range - 270 to 350 K, however the temperature >>>> distribution >>>> plot of the replicas show that the sampling has occurred at higher >>>> temperatures as well that is beyond 350K - >>>> Below is an excerpt from the temperature xvg file >>>> >>>> >>>> @ title "Gromacs Energies" >>>> @ xaxis label "Time (ps)" >>>> @ yaxis label "(K)" >>>> @TYPE xy >>>> @ view 0.15, 0.15, 0.75, 0.85 >>>> @ legend on >>>> @ legend box on >>>> @ legend loctype view >>>> @ legend 0.78, 0.8 >>>> @ legend length 2 >>>> @ s0 legend "Temperature" >>>> 0.000000 0.000000 >>>> 10.000000 350.997864 >>>> 20.000000 353.618927 >>>> 30.000000 350.068481 >>>> 40.000000 353.921753 >>>> 50.000000 359.485565 >>>> 60.000000 353.463654 >>>> 70.000000 352.015778 >>>> 80.000000 350.657898 >>>> 90.000000 351.927155 >>>> 100.000000 354.539429 >>>> 110.000000 354.287720 >>>> 120.000000 349.436096 >>>> 130.000000 352.960541 >>>> 140.000000 351.631317 >>>> 150.000000 354.217407 >>>> 160.000000 350.185852 >>>> 170.000000 350.294434 >>>> 180.000000 350.980194 >>>> 190.000000 350.914429 >>>> .... >>>> .... >>>> 470.000000 349.224060 >>>> 480.000000 350.819458 >>>> 490.000000 348.541748 >>>> 500.000000 350.393127 >>>> 510.000000 398.775208 >>>> 520.000000 444.802856 >>>> 530.000000 470.899323 >>>> 540.000000 466.652740 >>>> 550.000000 465.600677 >>>> 560.000000 469.225555 >>>> 570.000000 470.548370 >>>> 580.000000 470.011566 >>>> 590.000000 470.643951 >>>> 600.000000 472.433197 >>>> 610.000000 470.451172 >>>> 620.000000 469.991699 >>>> 630.000000 469.073090 >>>> 640.000000 467.259521 >>>> 650.000000 464.561798 >>>> 660.000000 468.416901 >>>> 670.000000 468.754913 >>>> 680.000000 469.259613 >>>> 690.000000 467.641144 >>>> 700.000000 468.542328 >>>> >>>> >>>> Temperature coupling was done using Nose hoover algorithm. >>>> >>>> Does this imply the sampling is wrong or insufficent? >>>> Any help / suggestion is appreciated. >>>> >>>> >>>> How large is your system, and what is it? What were your (full) .mdp >>> settings? The fact that your temperature started at 0 K and ramped up >>> suggests that you did not equilibrate prior to the run, did not generate >>> appropriate velocities, or did not continue properly. The sudden jump in >>> temperature later suggests instability, and could be due to incorrect >>> settings. N-H allows for large oscillations, but I wouldn't expect a >>> stable system to that degree. >>> >>> -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. >>> >>> > -- > ================================================== > > 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. > -- 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.