AMBER Archive (2002)

Subject: Re: a little bit confused with leap-frog algorithm...

From: Margaret Cheung (cheung_at_physics.ucsd.edu)
Date: Sat Mar 30 2002 - 21:35:36 CST

Dear David,

I just looked up the velocity-verlet in amber7. It looks like
that only sander uses the new algorithm, whereas sander_classic
still uses the old algorithm. Please correct me if I'm wrong. Thank you
very much.

Sincerely,
Margaret S. Cheung
Biophysics Program
Physics Department 0319
University of California, San Diego
9500 Gilman Drive,
La Jolla, CA 92093-0319
http://www-physics.ucsd.edu/~cheung

On Thu, 28 Mar 2002, David Case wrote:

> On Thu, Mar 28, 2002, Margaret Cheung wrote:
> >
> >
> > If I understand the leap-frog correctly, the expression for the new
> > velocities, based on the old velocities:
> >
> > v(t+dt/2) = v(t-dt/2) + dt* f(t)/m
> >
> > The update of the new coordiates is:
> >
> > r(t+dt)=r(t)+dt*v(t+dt/2)
> >
> > Which implied that kinetic energy (obtained at time t+dt/2) and potential
> > energy (obtained at time t+dt) can not be obtained
> > at the same time.
> >
>
> For "strict" Verlet, this is true. There are various approaches to getting
> the KE at the same time as you have the PE (from the force() routine).
> Amber 4,5,6 used the following ansatz:
>
> KE(t) approx= [ KE(t-dt/2) + KE(t+dt/2) ]/2
>
> i.e. averaging the kinetic energy at the two points where we know the
> velocities. These values are called EOLD3 and ENEW in the code.
>
> Amber 7 uses a slightly better algorithm: do a "half-step" update to
> propagate the velocities from t-dt/2 to t; then use the V(t) values to
> compute the KE.
>
> Note that the method of estimating the KE has no effect on the trajectory;
> this is why I say that how to get the KE goes "beyond" the Verlet
> integrator.
>
> [Note that all of the above becomes more complicated when SHAKE is turned on.]
>
>
> > However, in the code (runmd.f)
> >
> > velocity is obtained by
> > V(I3) = (V(I3) + F(I3)*wfac)
>
> This is just your first equation, above.
>
> >
> > and the coordinates are later obtained by
> > X(I3) = X(I3)+V(I3)*DTX
>
> This is your second equation, above.
>
> ...hope this helps....dac
>
> --
>
> ==================================================================
> David A. Case | e-mail: case_at_scripps.edu
> Dept. of Molecular Biology, TPC15 | fax: +1-858-784-8896
> The Scripps Research Institute | phone: +1-858-784-9768
> 10550 N. Torrey Pines Rd. | home page:
> La Jolla CA 92037 USA | http://www.scripps.edu/case
> ==================================================================
>
>