AMBER Archive (2003)

Subject: AMBER: PMEMD and EAMBER

From: Kristina Furse (kristina.e.furse_at_vanderbilt.edu)
Date: Thu Aug 07 2003 - 15:32:45 CDT

Hello-

We've been trying out PMEMD, using it for some potential of mean force
calculations (using NMR-style torsional constraints to bias the sampling, then
a WHAM method to get the PMF curve).

I ran some initial test runs to compare pmemd with sander in AMBER7 (running
on an SGI origin). While comparing the test runs, I noticed something about
the computation of EAMBER in minimization and dynamics output (below). Looking
at the first step of minimization output for the very same input coordinates
run in sander (amber7) and pmemd (amber7_compat=1), the energies are identical
(well, to be exact, EEL is different by 0.0001) except for EAMBER, which is
off by 2.5778kcal, which is exactly the restraint energy minus the 0.0001. I
summed the energies by hand and it looks like PMEMD may be subtracting the
restraint energy twice, or else not including it in the total potential energy
calculation then still subtracting it to get EAMBER?

I took a look at some dynamics output (same input coordinates as before) from
sander7 and pmemd, and I think the same thing may be going on there. For
pmemd, it seems like the number for EPtot already does not include restraint
energy, then EAMBER has restraint energy subtracted again.

Has anyone else noticed this? Did I miss something somewhere, or am I just
adding wrong...

Thanks!

Kristina

Minimization output:

Sander (amber 7):

   NSTEP ENERGY RMS GMAX NAME NUMBER
      1 -3.5193E+05 1.6203E+01 1.3187E+02 CD 17121

 BOND = 3402.1061 ANGLE = 9311.0284 DIHED = 11245.0567
 VDWAALS = 48314.8917 EEL = -472531.9934 HBOND = 0.0000
 1-4 VDW = 4345.3096 1-4 EEL = 43982.3959 RESTRAINT = 2.5779
 EAMBER = -351931.2050
==============================================================================
=
                      NMR restraints for step 1
 Energy (this step): Bond = 2.578 Angle = 0.000 Torsion = 0.000
 Energy (tot. run) : Bond = 2.578 Angle = 0.000 Torsion = 0.000

 DEVIATIONS: Target=(r2+r3)/2 Target = closer of r2/r3
            This step Entire run This step Entire run
           ave. rms ave. rms ave. rms ave. rms
 Bond 0.253 0.254 0.253 0.254 0.253 0.254 0.253 0.254
 Angle 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
 Torsion 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
==============================================================================
=

PMEMD (amber7_compat=1):

   NSTEP ENERGY RMS GMAX NAME NUMBER
      1 -3.5193E+05 1.6203E+01 1.3187E+02 CD 17121

 BOND = 3402.1061 ANGLE = 9311.0284 DIHED = 11245.0567
 VDWAALS = 48314.8917 EEL = -472531.9933 HBOND = 0.0000
 1-4 VDW = 4345.3096 1-4 EEL = 43982.3959 RESTRAINT = 2.5779
 EAMBER = -351933.7828
==============================================================================
=
                      NMR restraints for step 1
 Energy (this step): Bond = 2.578 Angle = 0.000 Torsion = 0.000
 Energy (tot. run) : Bond = 2.578 Angle = 0.000 Torsion = 0.000

 DEVIATIONS: Target=(r2+r3)/2 Target = closer of r2/r3
            This step Entire run This step Entire run
           ave. rms ave. rms ave. rms ave. rms
 Bond 0.253 0.254 0.253 0.254 0.253 0.254 0.253 0.254
 Angle 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
 Torsion 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
==============================================================================
=

Molecular Dynamics Output:

Sander (amber7)

 NSTEP = 100 TIME(PS) = 3150.100 TEMP(K) = 296.72 PRESS = 32.9
 Etot = -287279.7774 EKtot = 64194.7299 EPtot = -351474.5072
 BOND = 3498.2008 ANGLE = 9330.0143 DIHED = 11174.5893
 1-4 NB = 4342.5808 1-4 EEL = 44008.7977 VDWAALS = 48274.1435
 EELEC = -472109.4246 EHBOND = 0.0000 RESTRAINT = 6.5909
 EAMBER (non-restraint) = -351481.0982
 EKCMT = 25330.3166 VIRIAL = 24605.6087 VOLUME = 1019086.6920
                                                Density = 1.0541
 Ewald error estimate: 0.4495E-04

PMEMD (amber7_compat=1)

 NSTEP = 100 TIME(PS) = 3150.100 TEMP(K) = 296.72 PRESS = 32.8
 Etot = -287292.8061 EKtot = 64195.1417 EPtot = -351487.9478
 BOND = 3498.2009 ANGLE = 9330.0138 DIHED = 11174.5889
 1-4 NB = 4342.5814 1-4 EEL = 44008.7981 VDWAALS = 48276.0446
 EELEC = -472118.1755 EHBOND = 0.0000 RESTRAINT = 6.5910
 EAMBER (non-restraint) = -351494.5388
 EKCMT = 25330.5359 VIRIAL = 24608.2331 VOLUME = 1019085.5924
                                                Density = 1.0541
 Ewald error estimate: 0.4460E-04

**In case you wonder about the difference in restraint energies, the dynamics
runs have NMR-style constraints on the torsion we're driving, as well as a
"virtual bond" to a heme group accomplished using an NMR-style distance
constraint. The minimizations just have the virtual bond.

Finally, for good measure, dynamics input files:

sander (amber7)

Biased PMF: Cox dimer w/radical, explicit water, 9.0 cut
 &cntrl

  nmropt = 1,
  ntx = 5, irest = 1, ntrx = 1, ntxo = 1,
  ntpr = 100, ntwx = 100, ntwv = 0, ntwe = 0,

  ntf = 2, ntb = 2,
  cut = 9.0,

  ibelly = 0, ntr = 0,

  imin = 0,
  nstlim = 10000,
  nscm = 100,
  t = 0.0, dt = 0.001,

  temp0 = 298.0, tempi = 298.0,
  ig = 71277, heat = 0.0,
  ntt = 1,
  tautp = 1.0,
  vlimit = 20.0,

  ntp = 1, pres0 = 1.0, comp = 44.6,
  taup = 1.0,

  ntc = 2, tol = 0.00001, watnam = 'SPC ',

 &end
 &ewald

 nbflag = 1, skinnb = 2.0, use_pme = 1,

 &end
 &wt type='END' &end
DISANG=tor0.rst

pmemd input file just has the additional flag, amber7_compat=1,

****************************************************
Kristina E. Furse
Department of Chemistry
Center for Structural Biology
Vanderbilt University
email: kfurse_at_structbio.vanderbilt.edu

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu