Hi all,
I’m a little stuck on deriving an angular force constant for the amber forcefield.
Consider in amber, that the potential energy from the angle ABC is: U_{ABC} = K(theta - theta_eq)^2 The K is the force constant with the units kcal/(mol rad^2) which reduced to kcal/mol.
Using Gaussian, I am able to get out the Force Constant Matrix (Hessian) from the check point file which is an nxn matrix (n = number of atoms), second derivative of the potential energy relative to the cartesian displacement. The Hessian matrix is in using the Hartree/bohr-radius^2. This can be converted to kcal/(mol ang^2).
Now using the attached maths (see the file) I can pull out the angular force constant simply by taking two 3x3 subsets from the hessian of the displaced atom and the atom experiencing the force. But the problem I face is that I do not know how I can get the obvious displacement units, into radians. If you look at equation (9) there is even a distance_AB and distance_CB term, which is in Angstrom.
How can I get these calculations to reflect the same units at those in the amber force field for angular potential energy?
Thanks
Anthony
Anyone got any thoughts?
Can anyone explain how the distance term (see the image) in the final equations loses it's Angstrom^2 unit? Given that it must, if it were to fit into the Amber force field for a bond angle i.e., K (kcal/(mol rad^2)) * (Angle(rad) - Angle_eq(rad))^2
?
I am really looking for the same answer and I really couldn't understand why all "forcefield parametrize/development articles" don't explain step by step how to develop a forcefield from quantum mechanic calculations. In all articles, some points are missing and you can not even follow half of the instructions.
Hi Anthony Nash and Barış Kurt
I am going through the same process of deriving force fields from frequency analysis in Gaussian.
I am a bit confused about the angular force constant.
K_thetha in eq9 is in Hartree/radians^2 or Hartree/degree^2?
One way the distance term (Angstrom^2 or borh^2) in the denominator of eq9 is lost is that the sum term in the denominator has a Angstrom^-2 or borh^-2 component.
kcal/(mol ang^2).
Eq9 is then similar to:
1/hartree.rad^-2 = 1/borh^2*(hartree.rad^-2.borh^-2) + 1/borh^2*(hartree.rad^-2.borh^-2).
I don't know if that makes sense.
Vuyani Moses any idea?
- Badges
- Science topic
- Similar topics
- Chemistry
- Theoretical Chemistry
- Computational Chemistry
- Molecular Dynamics