Re: inconsistence in NAMD and AMBER energy

From: Floris Buelens (floris_buelens_at_yahoo.com)
Date: Tue Jul 21 2009 - 10:57:41 CDT

Hi,

> Now I need to think about how to get away with it. And I suggest all
people
> that are or want to do NAMD using AMBER parameters need to keep
in
> mind that VDW energy is NOT going to be the same if you are
using
> periodic boundary and NAMD 2.7b1, or older version.

the above is correct but I think it's fair to say that for most applications the impact will be minimal (mainly you'll see pressure and equilibrium density being slightly affected by the choice of cutoff; dynamics are actually not affected at constant vollume). But of course consistency is a good thing and the solution is to use the latest CVS code (or future releases after 2.7b1).
Best,

Floris

________________________________
From: Yi Shang <mirandaisbest_at_gmail.com>
To: Thomas C. Bishop <bishop_at_tulane.edu>
Cc: floris_buelens_at_yahoo.com; gouthambs_at_gmail.com; namd-l_at_ks.uiuc.edu
Sent: Tuesday, 21 July, 2009 16:56:04
Subject: Re: namd-l: inconsistence in NAMD and AMBER energy

Hi all,
Thank you for replying.
Thomas, are you using periodic boundary condition or not? Because if you look at my email. The energy gap is there only if I use periodic boundary condition. I was not using any SHAKE, temperature or pressure control.

Goutham, I am using same amber parameter file for both simulations when I did comparison, so I believe force parameters are the same.
 
Floris, you are right, by turning of vdw correction in amber simulation (vdwmeth = 0), I get following energies, notice that vdw energy is exactly the same as NAMD!! now all energies are consistent...
 NSTEP = 0 TIME(PS) = 0.000 TEMP(K) = 0.00 PRESS = 0.0
 Etot = -57081.5138 EKtot = 0.0000 EPtot = -57081.5138
 BOND = 341.5077 ANGLE = 716.0017 DIHED = 2372.7669
 1-4 NB = 1076.3715 1-4 EEL = 9100.8801 VDWAALS = 6132.6047
 EELEC = -76821.6464 EHBOND = 0.0000 RESTRAINT = 0.0000
 Ewald error estimate: 0.2976E-03

Now I need to think about how to get away with it. And I suggest all people that are or want to do NAMD using AMBER parameters need to keep in mind that VDW energy is NOT going to be the same if you are using periodic boundary and NAMD 2.7b1, or older version.

 
On Tue, Jul 21, 2009 at 10:33 AM, Thomas C. Bishop <bishop_at_tulane.edu> wrote:

I have been able to obtain the exact same energies from NAMD and amber for
>several time steps using suggested parameters in the NAMD manual.
>
>For this it is easiest to start with a temperature of 0 but you can also read
>in a velocity file s.t. they have same coords AND momentum
>
>Obviously you have to turn off anything that relies on random number generator
>>so NVE it primary target here.
>
>bonds, angle, dihedrals, improper agree w/out any problem
>
>vdw and and elec are tricky.
>first start with a simple cut-off to make sure your parameters are consistent.
>then turn on long range interactions. PAY CAREFUL ATTENTION TO NAMD SWITCHING
>>OPTIONS! see manual on swithcing functions too. This introduces a shift in
>the elec and vdw terms that will make it look like NAMD and AMBER are
>radically different but in terms of forces they are nearly same. make
>>switching agree and energies agree also.
>
>finally amber includes a long range vdw correction that will affect the
>trajectory over time but in first few steps ~10 this is minimal. W/ cutoff
>they should give almost same result for much longer
>
>Please post what you finally come up w/
>
>Tom
>
>
>
>On Tuesday 21 July 2009, Floris Buelens wrote:
>> Hi,
>>
>> I'm not sure about exactly how AMBER handles this issue, but I suspect the
>> difference might arise from the neglect of van der Waals interactions
>>> beyond the cutoff. I think AMBER applies an analytical correction by
>> default for these interactions (search the manual for 'vdwmeth'), while
>> NAMD as of 2.7b1 doesn't yet. I would suggest re-running your analysis
>>> using an increasing series of cutoff values, I think you should see the
>> NAMD value converging towards that of AMBER. Alternatively, the option for
>> analytical corrections for long-range van der Waals has now been added to
>>> the latest development code, which you can access through CVS as described
>> on the NAMD website, if you use that code with the option "LJCorrection on"
>> I hope you shouldn't see such a big discrepancy. Best wishes,
>>>
>> Floris
>>
>>
>> ________________________________
>> From: Yi Shang <mirandaisbest_at_gmail.com>
>> To: namd-l_at_ks.uiuc.edu
>>> Sent: Tuesday, 21 July, 2009 1:31:47
>> Subject: namd-l: inconsistence in NAMD and AMBER energy
>>
>>
>> Hi all,
>> I did some comparison of energy calculated by NAMD 2.7b1 and AMBER 10. I
>>> found something intriguing. During simulation, VDW energy from AMBER
>> simulation is constantly much lower from NAMD. So I tried to compare point
>> energy first. I turned off SHAKE, temperature, and pressure control and
>>> only looked at basic energy calculation of one single structure: bonds,
>> angles, dihedrals, electrostatic, and vdw energies. For non-periodic
>> boundary, I had no problem, energies are exactly the same. For periodic
>>> boundary condition here comes the problem:
>> -electro energy off by ~14kcal/mol, which is ok, parameters might not be
>> identical in two methods. -vdw energy off by > 600 kcal/mol. That's
>> terrible.
>>>
>> I tested these inputs on two different periodic systems, and I get the same
>> result (same energy gap). So I don't think it's system specific problem.
>> Here is my inputs and outputs from one system's comparison (I combined
>>> AMBER VDWAALS and 1-4 NB energies to compare to NAMD vdw energy, likewise
>> for electrostatic interaction):
>> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>> NAMD config.in:
>>> amber on
>> parmfile box.top
>> ambercoor box.crd
>> binaryoutput no
>> outputname md6
>> outputpressure 1
>> outputenergies 1
>> exclude scaled1-4
>>> 1-4scaling 0.83333
>> timestep 2
>> temperature 0
>> cutoff 8
>> switching off
>> PME yes
>> PMEGridSpacing 1.0
>> cellbasisvector1 79.6614 0 0
>>> cellbasisvector2 0 60.2065 0
>> cellbasisvector3 0 0 63.7403
>> cellorigin 0 0 0
>> run 0
>> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>> amber md.in:
>>> &cntrl
>> imin = 0, nstlim = 1, dt = 0.002,
>> irest = 0, ntx = 1,
>> tempi = 0.0, temp0 = 0.0,
>> ntwx = 0, ntwe = 0, ntwr = 0, ntpr = 1,
>> cut = 8.0,
>> /
>>
>> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>>> namd energies:
>> ETITLE: TS BOND ANGLE DIHED IMPRP
>> ELECT VDW BOUNDARY MISC KINETIC
>> TOTAL TEMP POTENTIAL TOTAL3
>>> TEMPAVG PRESSURE GPRESSURE VOL UME PRESSAVG
>> GPRESSAVG
>> ENERGY: 0 341.5078 716.0018 2372.7673 0.0000
>> -67739.0356 7208.9795 0.0000 0.0000 0.0000
>>> -57099.7791 0.0000 -57099.7791 -56529.2524 0.0000
>> 11740.3005 -276.1608 305707.0 250 11740.3005
>> -276.1608
>>
>> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>>> amber energies:
>> NSTEP = 0 TIME(PS) = 0.000 TEMP(K) = 0.00 PRESS =
>> 0.0 Etot = -57686.8091 EKtot = 0.0000 EPtot =
>> -57686.8091 BOND = 341.5077 ANGLE = 716.0017 DIHED =
>>> 2372.7669 1-4 NB = 1076.3715 1-4 EEL = 9100.8801 VDWAALS
>> = 5527.3094 EELEC = -76821.6464 EHBOND = 0.0000
>> RESTRAINT = 0.0000 Ewald error estimate: 0.2976E-03
>>>
>> %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
>>
>> It's strange because in NAMD seems most important parameter I would play
>> around with VDW energy is switching, and it's turned off. Same cutoff is
>>> set in two inputs, but still, I get very different results. any suggestions
>> would be appreciated! Thanks!
>
>
>
>--
>**********************
>Thomas C. Bishop *
>>Office: 504-862-3370 *
>Fax: 504-862-8392 *
>**********************
>

-- 
Miranda
      

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:53:03 CST