Leapfrog-Verfahren

Das Leapfrog-Verfahren ist eine einfache Methode zur numerischen Integration einer gewöhnlichen Differentialgleichung vom Typ

$\ddot x = \dot v = a(x)$.

bzw. allgemeiner von konservativen Systemen $\ddot{\mathbf{x}} =- M^{-1}\nabla V(\mathbf{x})$ der klassischen Dynamik, die z.B. die Bewegung eines oder mehrerer Objekte in einem Potentialfeld V beschreiben,

$\ddot {\vec x}_k =-\frac{1}{m_k}\nabla_{{\vec x}_k} V\left({\vec x}_1,\dots,{\vec x}_N\right),\quad k=1,\dots,N.$

Die Leapfrog-Integration ist eine Methode zweiter Ordnung und liefert deshalb im Allgemeinen genauere Ergebnisse als zum Beispiel das eulersche Polygonzugverfahren, das nur von erster Ordnung ist. Außerdem ist es invariant unter Zeitumkehr und erhält in physikalischen Problemstellungen Größen wie den Impuls und Drehimpuls, die auch Erhaltungsgrößen des Originalsystems sind, exakt. Des Weiteren wird eine gestörte Energiefunktion in Ordnung 3 erhalten, während das Verfahren die globale Konvergenzordnung 2 hat.

Darstellung als Leapfrog-Verfahren

Die Leapfrog-Integration berechnet abwechselnd die Positionen x und die Geschwindigkeiten v zu unterschiedlichen Zeitpunkten, ähnlich wie beim Bockspringen (engl. leapfrog). Die Schritt-Gleichungen für das Verfahren lauten:

$\begin{align} x_{i+1} &= x_i + v_{i+1/2}\, \Delta t\\[.3em] v_{i+3/2} &= v_{i+1/2} + a(x_{i+1}) \,\Delta t \end{align}$,

mit den Startwerten x₀ und v_{1 / 2}.

Darstellung als Einschrittverfahren

Durch lineare Interpolation von Zwischenwerten kann das Leapfrog-Verfahren als Kombination der zwei Varianten des symplektischen Eulerverfahrens betrachtet werden:

$\begin{align} \text{(SE1)}\quad&\left\{\begin{array}{rcl} v_{i+1/2}&=&v_i+a(x_i)\, \Delta t/2\\[.5em] x_{i+1/2}&=&x_i+v_{i+1/2}\, \Delta t/2 \end{array}\right.\\[.7em] \text{(SE2)}\quad&\left\{\begin{array}{rcl} x_{i+1} &=& x_{i+1/2} + v_{i+1/2}\, \Delta t/2\\[.5em] v_{i+1} &=& v_{i+1/2} + a(x_{i+1}) \,\Delta t/2 \end{array}\right.\\ \end{align}$

Jeder einzelne Schritt und damit auch die Zusammensetzung ist eine symplektische Transformation und erhält daher Volumina im Phasenraum. Daraus folgt auch die exakte Erhaltung von Impuls und Winkelgeschwindigkeit, soweit das exakte System diese erhält.

Darstellung als Mehrschrittverfahren

Eliminiert man aus der Leapfrog-Version die Geschwindigkeitsberechnungen, so ergibt sich

$\begin{align} (x_{i+1}-x_i)-(x_i-x_{i-1}) &=(v_{i+1/2}-v_{i-1/2})\,\Delta t\\[.3em] x_{i+1}-2x_i+x_{i-1} &= a(x_i)\,\Delta t^2\\[.7em] \implies\quad x_{i+1}=2x_i-x_{i-1}&+a(x_i)\,\Delta t^2, \text{ mit Startwerten }x_0, x_1=x_0+v_0\,\Delta t+\tfrac12 a(x_0)\,\Delta t^2 \end{align}$

das Verlet-Verfahren das als symmetrische Diskretisierung von $\ddot x=a(x)$ auch direkt hergeleitet werden kann. Diese Diskretisierung hat einen lokalen Fehler von O(Δt⁴) und damit (wegen der zweifachen Integration) einen globalen Fehler der Größe $O(t_N^2\Delta t^2)$ für die Differenz x(t_N) − x_N zwischen exakter und Näherungslösung zur Endzeit t_N = t₀ + NΔt. Eine weitere Variante ist unter dem Namen Velocity-Verlet bekannt^[1].

Geschichte

Eine erste Beschreibung dieses Verfahrens wurde von Richard Feynman in Isaac Newtons Principia von 1687 in einem Argument zur Herleitung der Keplerschen Gesetze aus den Bewegungsgleichungen gefunden. Neben anderen wurden Varianten dieses Verfahrens 1860 von J. F. Encke und 1907 von C. Störmer verwendet.^[2]

Beispiel

Betrachtet man die Schwingungsgleichung $\ddot x+x=0$ mit der exakten Lösung $x(t)=C\,\cos (t+\phi)$, so erhält man für $\dot x=v,\;\dot v=-x$ in der Einschrittformulierung den Übergang

$\binom{x_{k+1}}{v_{k+1}}= \begin{pmatrix}1-\tfrac12 h^2&h\\-h(1-\tfrac14 h^2)&1-\tfrac12 h^2\end{pmatrix} \binom{x_k}{v_k}.$

Es stellt sich heraus, dass das modifizierte Energiefunktional $E_h(v,x)=\tfrac12\left(v^2+(1-\tfrac14h^2)x^2\right)$ exakt erhalten wird, genauer gilt E_h(v_n,x_n) = E(v₀,x₀). Damit hat der Fehler in der Energie eine globale Schranke der Ordnung O(h²). Die Näherungslösung verläuft für alle Zeiten auf den durch das konstante modifizierte Energieniveau definierten Ellipsen.

Quellen

1. ↑ Michael Griebel, Stefan Knapek, Gerhard Zumbusch, Attila Caglar: Numerische Simulation in der Moleküldynamik. Springer, 2004.
2. ↑ Ernst Hairer, Christian Lubich, Gerhard Wanner: Geometric numerical integration illustrated by the Störmer/Verlet method. In: Acta Numerica. 12, 2003, S. 399-450. doi:10.1017/S0962492902000144.