[gmx-users] Calculation of dihedral force constants

[Mark Abraham]

On 22/06/2012 11:45 PM, Covington, Cody Lance wrote:
> Dear gmx-users
>
> I am trying to calculate dihedral force constants from the diagonal elements in the hessian matrix in internal coordinates.
> The matrix is given in hartree/rad^2. Upon conversion to kJ/(mol*rad^2) the values seem an order of magnitude higher than I have seen in the force fields (~3.77 kJ/mol for CCCH dihedral but I get ~42 kJ/(mol*rad^2). The conversion I have used is 2625.5 kJ/(mol*hartree) but is there any consideration needed for the angular units?
> Another thought I had was,
> The GROMACS dihedral Potential Energy function #1 (eq 4.62) listed in the manual 4.5.4 is   Vd=k(1+Cos(nx-x0)) where x is dihedral angle.
> so the curvature of this function (second derivative) at the minimum is (n^2)*k
> so do I need to convert my force constants to this function by k_harm /n^2 = k ?

You should really be looking at the literature for the force field you 
are trying to extend or emulate, and see what data and method they used 
to fit these values. The GROMACS documentation is of necessity generic. 
Some force fields have a scaled 1,4 electrostatic interaction that will 
vary with the dihedral angle. There can be multiple dihedral angle 
functions around a single bond and these may add together.

Mark