[gmx-users] System does not heat up to 300K.

[Marc Hömberger]

Hi,

I have a system setup (minimized etc.) and wanted to heat this system
gradually with simulated annealing. I set up the mdp file as you can see
below (simulated annealing part is marked bold). But when the heating is
finished, the log file tells me that the temperature of the system is only
1K (see below).

Can someone tell me what I am doing wrong?


Best,
Marc

;************ PREPROCESSING *********************************************
;include                =
;define                 =
;************ RUN CONTROL ***********************************************
dt                      = 0.002         ; timestep, 2 fs
integrator              = md            ; leap-frog integrator
nsteps                  = 500000     ; 2 * 500000 = 1000 ps
tinit                   = 0             ; starting time, ps
init_step               = 0             ; starting step
comm-mode               = linear        ; remove COM translation
nstcomm                 = 10            ; freq steps comm-mode
comm-grps               = system        ; groups COM removal ; protein
nstcalcenergy           = -1            ; freq calc ener; default = nstlist
;************ DYNAMICS **************************************************
;### Initial velocity assignment
gen_vel                 = yes           ; assign velocities from Maxwell
distribution
gen_temp                = 300           ; temperature for Maxwell
distribution
gen_seed                = -1            ; random generator seed; -1 from
job ID
;************ OUTPUT CONTROL ********************************************
;### xtc
nstxtcout               = 2500          ; .xtc freq write coor
xtc-precision           = 1000          ; .xtc precision
xtc_grps                = system        ; .xtc saved groups
;### trr
nstxout                 = 500000        ; .trr freq write coor
nstvout                 = 500000        ; .trr freq write vel
nstfout                 = 0             ; .trr freq write forces
;### log
nstlog                  = 2500          ; .log freq write ener
;### edr
nstenergy               = 2500          ; .edr freq write ener; =
i*nstcalcenergy
energygrps              = protein non-protein   ; .edr saved groups
;************ NON-BONDED INTERACTIONS ***********************************
ns_type                 = grid          ; search neighbor grid cells
nstlist                 = 10            ; freq update NB list

vdwtype                 = switch        ; switching function for VdW
coulombtype             = PME           ; PME for long-range electrostatics
rvdw_switch             = 0.9           ; switching dist VdW
rvdw                    = 1.0           ; cutoff VdW
rcoulomb                = 1.2           ; cutoff elec
rlist                   = 1.2           ; cutoff NB list; must = rcoulomb
when PME

epsilon-r               = 1             ; dielectric const
DispCorr                = EnerPres      ; Long-range VdW dispersion
correction

;### PME & Ewald
ewald_rtol              = 1.0e-5        ; direct space tolerance;
1.0e-6/2.2e-5
fourierspacing          = 0.12          ; max grid spacing for FFT
pme_order               = 6             ; interpolation order; depends on
cutoff
ewald_geometry          = 3d            ; Ewald in 3d
epsilon_surface         = 0             ; dipole correction for 3d Ewald
optimize_fft            = yes           ; start-up FFT grid optimization
;************ BOUNDARY CONTROL ******************************************
pbc                     = xyz           ; 3-D PBC
;************ CONSTRAINTS CONTROL ***************************************
constraints             = all-bonds     ; all bonds (even heavy atom-H
bonds) constrained
continuation            = no            ; constraints applied to initial
conf
constraint_algorithm = lincs    ; holonomic constraints
lincs_iter              = 1             ; accuracy of LINCS
lincs_order             = 4             ; also related to accuracy
lincs_warnangle         = 30            ; max allowed bond angle before
warning
;************ TEMPERATURE CONTROL ***************************************
tcoupl                  = nose-hoover   ; Nose-Hoover T-coupling; v-rescale
tc-grps                 = system        ; two coupling groups - more
accurate
tau_t                   = 0.6           ; coupling time constant at equil,
in ps
ref_t                   = 300           ; reference temperature, one for
each group, in K
;************ PRESSURE CONTROL ******************************************
pcoupl                  = Parrinello-Rahman     ; Pressure coupling on in
NPT
pcoupltype              = isotropic     ; uniform scaling of box vectors
tau_p                   = 1.0           ; coupling time constant at equil,
in ps
ref_p                   = 1.0           ; reference pressure, in bar
compressibility         = 4.5e-5        ; isothermal compressibility of
water at 1 atm, bar^-1
;************ SIMULATED ANNEALING ******************************************
*annealing               = single        ; type of annealing (single
sequence of annealing points)*
*annealing_npoints       = 2             ; # of temperature points*
*annealing_time          = 0 300000      ; time window in which temperature
shall be increased*
*annealing_temp          = 0 300         ; temperature at timepoints*


LOG-FILE:
Writing checkpoint, step 500000 at Thu Jun  6 21:23:55 2013


*Current ref_t for group System:      1.0*
   Energies (kJ/mol)
            U-B    Proper Dih.  Improper Dih.      CMAP Dih.          LJ-14
    8.64012e+03    5.98740e+03    4.75546e+02   -1.72734e+03    4.25420e+03
     Coulomb-14        LJ (SR)  Disper. corr.   Coulomb (SR)   Coul. recip.
    4.78689e+04    1.75858e+05   -5.83865e+03   -1.21361e+06   -1.12292e+05
 Position Rest.      Potential    Kinetic En.   Total Energy    Temperature
    7.52594e+00   -1.09037e+06    1.35115e+05   -9.55257e+05    2.17729e+02
 Pres. DC (bar) Pressure (bar)   Constr. rmsd
   -1.35977e+02    1.38325e+02    2.75899e-05