Actually, if you read on you will find that both variants are represented on the wikipedia page.
Basic Verlet
The basic second order central difference quotient discretization for second order ODE x''(t)=a(x(t)) is
xn+1 - 2*xn + xn-1 = an*dt^2
Note that there are not velocities in the iteration and also not in the acceleration function a(x). This is because Verlet integration is only then superior to other integration methods when the dynamical system is conservative, that is, -m*a(x) is the gradient of some potential function, and potential functions are static objects, they depend only on the position, not on time and not on velocity. Many frictionless mechanical systems fall into this category.
Velocity Verlet
Now set, using first order central difference quotients, the velocity at time tn to
vn*(2*dt) = xn+1 - xn-1
and add and subtract this to and from the first equation to obtain
-2*xn + 2*xn-1 = -2*vn*dt + an*dt^2
2*xn+1 - 2*xn = 2*vn*dt + an*dt^2
or
vn = (xn - xn-1)/dt + 0.5*an*dt
xn+1 = xn + vn*dt + 0.5*an*dt^2
This is one variant to write the velocity-Verlet algorithm.
(updated to have all state variables correspond to the same time before and after the iteration step)
Using the equations of the previous step from n-1 to n, one can replace xn-1 by vn-1 and an-1 in the velocity computation. Then
vn = vn-1 + 0.5*(an-1 + an)*dt
To avoid having two instances of any of the vectors x,v,a, one can arrange the update process so that everything is in place. Assume that at entry of the iteration step, the stored data corresponds to (tn-1,xn-1,vn-1,an-1). Then the next state is computed as
vn-0.5 = vn-1 + 0.5*an-1*dt
xn = xn-1 + vn-0.5*dt
Do collision detection with xn and vn-0.5
an = a(xn)
vn = vn-0.5 + 0.5*an*dt
Do statistics with xn and vn
or as code
v += a*0.5*dt;
x += v*dt;
do_collisions(x,v);
a = eval_a(x);
v += a*0.5*dt;
do_statistics(x,v);
Changing the order of these operations will destroy the Verlet scheme and significantly alter the results, rotation of the oparations is possible, but one has to be careful about the interpretation of the state after the iteration step.
The only initialization needed is the computation of a0 = a( x0 ).
Leapfrog Verlet
From the formulas of velocity Verlet one can see that for the updates of the position one does not need the velocities vn, but only the half point velocities vn+0.5. Then
an = a(xn)
vn+0.5 = vn-0.5 + an*dt
xn+1 = xn + vn+0.5*dt
or in code
a = eval_a(x);
v += a*dt;
x += v*dt;
Again, the order of these operations is fundamentally important, changes will lead to strange results for conservative systems.
(Update) However, one may rotate the execution sequence to
x += v*dt;
a = eval_a(x);
v += a*dt;
This corresponds to the iteration of the triples (tn,xn,vn+0.5) as in
xn = xn-1 + vn-0.5*dt
an = a(xn)
vn+0.5 = vn-0.5 + an*dt
The initialization only needs to compute
v0+0.5 = v0 + 0.5*a( x0 )*dt
(end update)
Any statistics that one computes using xn and vn-0.5 or vn+0.5 will be off by an error proportional to dt, since the time indices do not match. Collisions can be detected with the position vector alone, but in the deflection the velocities need to be sensibly updated too.
parse error in constructor in data/newtype declaration: (Double, Double, Double)I'm going to take this to another thread. – Gadmann