**From:** Jim Phillips (*jim_at_ks.uiuc.edu*)

**Date:** Mon Jan 29 2007 - 16:55:40 CST

**Next message:**Richard Wood: "a philosophical question..."**Previous message:**Shirley Hui: "Questions re: protein-DNA complex"**In reply to:**Floris Buelens: "FAO developers: vdW energies, forces and derivatives"**Messages sorted by:**[ date ] [ thread ] [ subject ] [ author ] [ attachment ]

Hi,

The interpolation tables of U(r^2) are used to calculate dU/dr^2. This is

then used to calculate forces without calling sqrt() with the formula:

force = -dU/dr = -2 r dU/dr^2

The force vector is calculated simply by multiplying the vector between

the two interacting atoms by -2 dU/dr^2. Note that to save operations the

factor of 2 is already included in some of the table coefficients.

If the actual force was interpolated from the table it would be necessary

to call sqrt() to calculate the unit vector between atoms, which is very

expensive compared to multiplication and addition.

-Jim

On Mon, 29 Jan 2007, Floris Buelens wrote:

*> Hello NAMD users and developers,
*

*>
*

*> I have a question about the relationship between Lennard-Jones energy and force in NAMD. From an old NAMD programming guide (v1.5, 1998), I read:
*

*>
*

*> vdW energy = (A / r^12) - (B / r^6) - the basic Lennard-Jones potential function;
*

*>
*

*> It's then stated that the force is obtained by differentiating that function:
*

*>
*

*> vdW force = (12A / r^13) - (6B / r^7)
*

*>
*

*> However I can't seem to make this agree with the implementation I find in the 2.6 source... ComputeNonbondedUtil.C pre-calculates tables of r^2 versus various nonbonded energy parameters at startup. For the vdW forces, the calculation and differentiation all appear to be based on r^2, not r - I interpret:
*

*>
*

*> r_2 = 1 / r^2
*

*> r_6 = r^(-6) = (r_2) ^ 3;
*

*> r_12 = (r_6) ^ 2;
*

*> van der Waals energy = A * r_12 - B * r_6
*

*>
*

*> I'm clear on this so far, but then the routine stores vdwa_gradient and vdwb_gradient as:
*

*>
*

*> vdwa_gradient = -6 * r_2 * r_12 = -6 * r^(-14) [this corresponds to dE / dr^2 of E = (r^2)^-6 -> -6*(r^2)^-7]
*

*> vdwb_gradient = -3 * r_2 * r_6 = -3 * r^(-8) [this corresponds to dE / dr^2 of E = (r^2)^-3 -> -3(r^2)^-4]
*

*>
*

*>
*

*> (assuming for convenience that switching function is not active.)
*

*>
*

*> These gradient values are then later (in ComputeNonbondedBase2) multiplied by A and B respectively to give the vdW force.
*

*>
*

*> My understanding and the NAMD programming guide suggest we need to calculate dE / dr to get the force, while as I understand it, NAMD 2.6 calculates dE / dr^2. (12A / r^13) is not the same as (6A / r^14), and (6B / r^7) is not the same as (3B / r^8) - I don't understand the discrepancy. I need to get this straight in order to implement a modified vdW potential for free energy calculations. Can anybody help?
*

*>
*

*> Thanks in advance!
*

*>
*

*> Floris Buelens
*

*> Department of Crystallography, Birkbeck College, London
*

*>
*

*>
*

*>
*

*>
*

*>
*

*> p.s. - for comparison, the equivalent formula for electrostatics follows my expectations:
*

*> energy = Cqiqj / r = Cqiqj * r^-1
*

*> force = -1 * (Cqiqj) ^ -2 = -Cqiqj / r^2
*

*> ... from the programming guide, corresponds to:
*

*>
*

*> fast_energy = 1.0/r
*

*> fast_gradient = -1.0/r2
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

*>
*

*> ___________________________________________________________
*

*> New Yahoo! Mail is the ultimate force in competitive emailing. Find out more at the Yahoo! Mail Championships. Plus: play games and win prizes.
*

*> http://uk.rd.yahoo.com/evt=44106/*http://mail.yahoo.net/uk
*

*>
*

**Next message:**Richard Wood: "a philosophical question..."**Previous message:**Shirley Hui: "Questions re: protein-DNA complex"**In reply to:**Floris Buelens: "FAO developers: vdW energies, forces and derivatives"**Messages sorted by:**[ date ] [ thread ] [ subject ] [ author ] [ attachment ]

*
This archive was generated by hypermail 2.1.6
: Wed Feb 29 2012 - 15:44:20 CST
*