Re: [AMBER] Temperature regulation in HREMD

From: dani sd <daniandphy.yahoo.com>
Date: Wed, 15 Aug 2012 11:39:59 -0700 (PDT)

Hi Brian,

I think the problem is coming from HREMD code. After each  successful exchange, there is a relaxation time, which is related to the natural fluctuation of the system. In the current implementation, the system will be out of equilibrium for a while after each successful exchange. This problem is bold especially in the big systems because of  large relaxation times.
We solved the problem by modifying the code, however the modification is not uploaded to the amber repository. I will send the modified HREMD code to your email and you can recompile "sander.MPI".
Hope it solves the problem.

Best,
Danial


Danial Sabri Dashti
Quantum Theory Project &
Physics Department
University of Florida
Ph.D. Graduate Student


________________________________


From: Brian Radak <radak004.umn.edu>
Date: Wed, Aug 15, 2012 at 9:39 AM
Subject: [AMBER] Temperature regulation in HREMD
To: AMBER Mailing List <amber.ambermd.org>


I am running HREMD and have observed what seems to be unexpected behavior
in the instantaneous temperature. That is, although all of the replicas are
propagated using Langevin dynamics at 300K, the average temperatures are
generally around 270-280K. There is also a large variance in the total
energy (about 1000 kcal/mol), although perhaps this is ok and I am just not
understanding the output (do the mdout files track the replicas or track
the coordinates?)

I am using a patched version of AMBER 12. The system is a histidine residue
in a ~50 A box of TIP4P-Ew water. There are 12 replicas (although I also
observed this with 2 replicas) which consist of the atomic charges being
linearly scaled between the protonated and deprotonated forms (this was
accomplished using parmed). No changes to the atom types were performed.

All replicas had the same mdin:
 &cntrl
   nstlim = 250, numexchg = 1000, dt = 0.002,
   ntpr = 250, ntwx = 250,
   ntwr = 250000,
   ioutfm = 1, iwrap = 1, ig = -1,
   ntc = 2, ntf = 2,
   cut = 9,
   ntt = 3, temp0 = 300.0, tempi = 300.0, gamma_ln = 1.0
   /
 &ewald
  dsum_tol=1.0e-6
 /

The starting structures were generated by 200 ps of constant pressure
Langevin dynamics with no exchanges and then all of the box coordinates
changed to be the same. I thought it might have been a mistake not to use
irest = 1 when starting REMD, but am reluctant to blame that for the
observed behavior.

The free energies seem to be well converged based on the last 2 entries
from the remlog, but I'm sure this could just be misleading:

# Rep#, Neibr#, Temp0, PotE(x_1), PotE(x_2), left_fe, right_fe, Success,
Success rate (i,i+1)
# exchange      999
     1    12    300.00 -44753.10 -44571.23     -2.76      0.00    F
0.66
     2     3    300.00 -44301.78 -44699.62      0.00     -2.87    T
0.37
     3     2    300.00 -44696.59 -44298.66      2.82      0.00    T
0.40
     4     5    300.00 -44632.76 -44140.71      0.00     -0.41    F
0.40
     5     4    300.00 -44142.14 -44631.30      0.54      0.00    F
0.43
     6     7    300.00 -44282.42 -44425.70      0.00      1.75    F
0.52
     7     6    300.00 -44427.33 -44282.56     -1.74      0.00    F
0.45
     8     9    300.00 -44584.41 -44984.53      0.00      3.62    F
0.52
     9     8    300.00 -44988.08 -44586.89     -3.60      0.00    F
0.50
    10    11    300.00 -44692.44 -44763.92      0.00      5.55    T
0.53
    11    10    300.00 -44769.07 -44697.58     -5.52      0.00    T
0.71
    12     1    300.00 -44618.67 -44713.88      0.00    -33.86    F
0.00
# exchange     1000
     1     2    300.00 -44275.48 -44811.39      0.00     -1.91    T
0.66
     2     1    300.00 -44809.19 -44273.45      1.89      0.00    T
0.37
     3     4    300.00 -44737.41 -44249.38      0.00     -1.63    F
0.40
     4     3    300.00 -44248.15 -44735.31      1.57      0.00    F
0.40
     5     6    300.00 -43832.07 -43837.66      0.00      0.62    T
0.44
     6     5    300.00 -43838.34 -43833.08     -0.71      0.00    T
0.52
     7     8    300.00 -43958.10 -44146.05      0.00      2.64    T
0.46
     8     7    300.00 -44148.87 -43960.69     -2.73      0.00    T
0.52
     9    10    300.00 -44596.36 -44997.54      0.00      4.56    F
0.50
    10     9    300.00 -45002.94 -44600.83     -4.63      0.00    F
0.53
    11    12    300.00 -44792.55 -44199.67      0.00      3.14    F
0.71
    12    11    300.00 -44202.75 -44795.50     -3.16      0.00    F
0.00


Am I doing something ill advised here or making some other silly error?

Thanks,
Brian

--
================================ Current Address =======================
 Brian Radak                                             :     BioMaPS
Institute for Quantitative Biology
 PhD candidate - York Research Group       :     Rutgers, The State
University of New Jersey
 University of Minnesota - Twin Cities         :     Center for Integrative
Proteomics Room 308
 Graduate Program in Chemical Physics     :     174 Frelinghuysen Road,
 Department of Chemistry                          :     Piscataway, NJ
08854-8066
 radak004.umn.edu                                 :
radakb.biomaps.rutgers.edu
====================================================================
Sorry for the multiple e-mail addresses, just use the institute appropriate
address.
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
-- 
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Wed Aug 15 2012 - 12:00:07 PDT
Custom Search