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.
