On Mon, Dec 10, 2012 at 5:31 AM, Jan-Philip Gehrcke <jgehrcke.googlemail.com
> wrote:
> Hello,
>
> I am wondering how forces with jar=1 are actually calculated in case of
> a distance restraint between two atoms. The way I understand it, a
> 'correctional' force is dynamically calculated with each MD time step
> and applied in order to force the distance d between two atoms on a
> certain trajectory d(t), which is a linear relation in this case.
>
> In sander's nmr.F90, I've found `fcurr = -2.0*rk2*(rint-r2)` (line 3778
> for Amber 12).
>
The forces with jar==1 are calculated the same way as the forces with
nmropt==1 (e.g., as with umbrella sampling)--as in, they use the same
subroutine calls. The difference with jar=1 is that you also have to
accumulate the work done by the moving potential. Since work is not a
vector (it's the inner product of the force and position vectors), it has
no direction.
> First of all, I am wondering about the direction of the force. I'm not a
> Fortran expert, but assuming all of the variables in the code line above
> are scalars, only the force magnitude is calculated here. The direction
> should be given by the gradient of the according potential. Where/how is
> the direction of the force calculated?
>
The variables you quoted above are part of the work calculation, not the
force calculation. Hence, no direction ;). The forces and energies are
calculated in different subroutines on a per-restraint basis. If you are
looking for a distance, look for the subroutine "disnrg" in nmr.F90. This
is where the force (and its direction) are calculated. The direction comes
at the end, where dfr(m) is calculated. Note that xij is a vector (since
it is just xi-xj), which gives the direction and magnitude in each (x, y,
and z) direction. Also note that, as per Newton's 3rd law, that the force
exerted by the restraint on the first atom is the same, but opposite of
that on the second atom. In this subroutine, f() is the force array.
The code is a bit convoluted due to the numerous options and fairly complex
mathematical form of the restraint potential, but the ideas are still
there. You can look at angnrg, tornrg, etc. for how the forces are
calculated and accumulated for the different types of restraints.
HTH,
Jason
--
Jason M. Swails
Quantum Theory Project,
University of Florida
Ph.D. Candidate
352-392-4032
_______________________________________________
AMBER mailing list
AMBER.ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber
Received on Mon Dec 10 2012 - 09:00:03 PST
