[gmx-users] NVE Temperature Drift

[Johnny Lu]

Hi. I'm trying to run NVE on a
~170 amino acid protein
with a harmonic dihedral angle restraint of 4184 kJ/mol rad^2,
in a rhombic dodecahedron water box that is 15 angstroms away from the
protein
using single precision gromacs,
and 0.5 fs time step.

Since I only want to get many short stretches of NVE simulations, I
alternate npt and nve:

0. Steepest descent minimization with emtol = 1000.0.
1. [NVT equilibration, 2fs timestep 100,000 frames]
2. [NPT equilibration, 2fs timestep, 100,000 frames]
3. [NPT equilibration, 2fs timestep, 5,000 frames, Berendsen thermostat at
300K, position restraint?]
4. [NVE md,0.5 fs timestep, 35,000 frames, no thermostat]
5. go to step 3.

(all of these steps have the dihedral angle restraint)

In the NVE step, I get about 0.1% total energy decrease from start to
end(figure). (the total energy keep dropping).
But, the temperature still fluctuate in a range of ~ 5 Kelvin (figure
<http://oi60.tinypic.com/sni5at.jpg>). Is that normal? Does that happen in
real world? Where does this temperature fluctuation come from? How to
reduce that? Does the restraint have anything to do with this?
[image: Inline image 2]

[image: Inline image 1]

Below is the mdp file that I was using:

; Amber99SB MD
; Run parameters
integrator    = md        ; leap-frog integrator
nsteps        = 32768        ; 2 * 100000 = 200 ps, 1 ns
dt        = 0.0005        ; 1 fs
; Output control
nstfout        = 2
nstxtcout    = 2        ; xtc compressed trajectory output every 2 fs
nstlog        = 1000        ; update log file every 2 ps
xtc-grps=protein;

continuation    = yes        ; Restarting after NPT
constraint_algorithm = Lincs    ; holonomic constraints
constraints    = h-bonds    ; all bonds (even heavy atom-H bonds)
constrained
lincs_iter    = 2        ; accuracy of LINCS
lincs_order    = 4        ; also related to accuracy
; Neighborsearching
ns_type        = grid        ; search neighboring grid cells
nstlist        = 10        ; 10 fs
rlist        = 1.8        ; short-range neighborlist cutoff (in nm)
rcoulomb    = 1.5        ; short-range electrostatic cutoff (in nm)
rvdw        = 1.5        ; short-range van der Waals cutoff (in nm)
; Electrostatics

coulombtype    = PME        ; Particle Mesh Ewald for long-range
electrostatics
pme_order    = 6        ; cubic interpolation

fourierspacing    = 0.12        ; grid spacing for FFT
optimize_fft    = yes

tcoupl                   = no
pcoupl                   = no

compressibility = 4.5e-5    ; isothermal compressibility of water, bar^-1
; Periodic boundary conditions
pbc        = xyz        ; 3-D PBC
; Dispersion correction
DispCorr    = EnerPres    ; account for cut-off vdW scheme
; Velocity generation
gen_vel        = no        ; Velocity generation is off
verlet-buffer-drift=-1
cutoff-scheme = Verlet;

Thanks again.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20140825/ef7ec591/attachment.html>