Lesson 1.5 · 14 min

Euler vs velocity Verlet

There are many ways to turn F = ma into steps. The obvious one, Euler, slowly pumps energy into the system. Velocity Verlet costs the same but keeps energy steady and can even run backwards. Leapfrog, used by GROMACS, moves the atoms along exactly the same path; only the way it reports the energy differs.

Why the simple Euler update makes energy grow
How velocity Verlet works: half kick, move, half kick
Two properties that make Verlet the MD standard: no energy drift and time reversibility
Leapfrog: the same positions as Verlet, with velocities at half steps, and what that does to the reported energy

The idea

In Lesson 1.3 we updated the velocity first and then moved the atom. The most obvious rule does it the other way round. It moves the atom using its old velocity, then updates the velocity with the old force. This is the explicit Euler method.

For an orbit, Euler always steps along the tangent of the circle, which points slightly outward. Every step lands a bit too far out, so the orbit spirals outward and the energy grows without limit. A simulation like that heats itself up and is useless.

Velocity Verlet splits the step into three moves:

Half kick: update the velocity with half of the current force.
Move: update the position with that half-updated velocity.
Half kick: compute the force at the new position, and use it for the other half of the velocity update.

The step is symmetric in time. Because of that:

The energy wiggles slightly but does not drift, even over millions of steps.
It is reversible: flip all velocities and the system retraces its path back to the start.
It needs only one force calculation per step, the same as Euler, because the new force is reused in the next step.

Leapfrog is the same idea written with velocities at half steps: kick the velocity a whole step, from v(t − Δt/2) to v(t + Δt/2), then move the atom with it. Starting from the same point, it puts every atom at exactly the same positions as velocity Verlet. Its catch: it never has the velocity at time t, so it must estimate the kinetic energy, and different estimates report slightly different energies (section 5 of the maths). Neither estimate drifts.

Key idea. A good integrator is not just "accurate for one step". It must respect the physics over millions of steps. Verlet does: same cost as Euler, far better behaviour.

The algorithm

Velocity Verlet in a real MD loop. Note that the force from the end of one step is reused at the start of the next:

F = compute_forces(pos)                 # once, before the loop
for step in range(n_steps):
    vel += 0.5 * dt * F / mass          # half kick
    pos += dt * vel                     # move
    F = compute_forces(pos)             # the only force call this step
    vel += 0.5 * dt * F / mass          # half kick

For comparison, explicit Euler and Euler–Cromer (Lesson 1.3):

# explicit Euler: energy grows
F = compute_forces(pos)
pos, vel = pos + dt * vel, vel + dt * F / mass

# Euler–Cromer: stable but first-order and not reversible
F = compute_forces(pos)
vel += dt * F / mass
pos += dt * vel

Leapfrog (GROMACS integrator = md) keeps the velocity half a step ahead of the position:

v_half = vel - 0.5 * dt * F / mass      # v(−dt/2), once, before the loop
for step in range(n_steps):
    F = compute_forces(pos)
    v_new = v_half + dt * F / mass      # whole kick: v(t − dt/2) → v(t + dt/2)
    KE = 0.5 * mass * (v_half**2 + v_new**2) / 2   # GROMACS: average of the two half-step KEs
    pos += dt * v_new                   # move with the half-step velocity
    v_half = v_new

Positions are the same as velocity Verlet's (the average of the two half-step velocities is exactly Verlet's v(t)). What is not fixed is the kinetic energy at time t: averaging the two half-step velocities reproduces Verlet's energy; averaging the two half-step kinetic energies, as GROMACS does, reports a slightly different one.