L'algorithme prédicteur-correcteur de Gear

L'erreur commise par l'algorithme de Verlet est en $O(\Delta t ^4)$ . Celle-ci peut être diminuée en augmentant par exemple l'ordre du développement limité utilisé [29,30,31].

Si nous prenons un développement de cinquième ordre, comme ce sera le cas dans nos simulations nous obtenons le système d'équations suivant:

$\vec r_i (t+ \Delta t)= \vec r_i (t) + \vec r^{'}_i (t) \frac{ ( \Delta t)^2} ... ...vec r^{iv}_i (t)\frac{( \Delta t)^5}{5!} +\vec r^{v}_i \frac{( \Delta t)^6}{6!}$

$\vec r^{'}_i (t + \Delta t) = \vec r{'}_i (t) + \vec r^{''}_i (t) \frac{ ( \De... ...{iv}_i (t) \frac{( \Delta t)^4}{4!} + \vec r^{v}_i (t)\frac{( \Delta t)^5}{5!}$

$\vec r^{''}_i (t + \Delta t) = \vec r{''}_i (t) + \vec r^{'''}_i (t) \frac{ (... ...v}_i (t) \frac{( \Delta t)^3}{3!} + \vec r^{v}_i (t) \frac{( \Delta t)^4}{4!}$

$\vec r^{'''}_i (t + \Delta t) = \vec r^{'''}_i (t) + \vec r^{iv}_i (t) \frac{ ( \Delta t)2} {2!} + \vec r^{v}_i (t) \frac{( \Delta t)^3}{3!}$

$\vec r^{iv}_i (t + \Delta t) = \vec r^{iv}_i (t) + \vec r^{v}_i (t) \frac{ ( \Delta t)^2}{2!}$

$\vec r^{v}_i (t + \Delta t) = \vec r^{v}_i (t)$

La précision du développement en série de Taylor dépend de l'ordre de celui-ci. Nous appellerons les valeurs obtenues à partir du développement "valeurs prédites". Nous pouvons corriger ces valeurs en remarquant que la force, c'est-à-dire le gradient du potentiel est fonction de la position. Nous calculons à nouveau la valeur du potentiel, au niveau de la position prédite. Ainsi, nous déduisons une nouvelle valeur de l'accélération et donc une position différente pour le pas de temps suivant et ainsi l'erreur commise:

$$\Delta \vec a(t+\Delta t) = \vec a_{corr}(t+\Delta t) - \vec a_{pred}(t+\Delta t)$$ (2.13)

L'erreur sur l'accélération est ensuite utilisée pour corriger la position prédite. Dans cette équation nous réservons l'indice "pred" pour la valeur prédite et l'indice "corr" pour la valeur corrigée.

$$\vec r_{corr}(t+\Delta t) = \vec r_{pred} (t+\Delta t) + a_0 \Delta \vec r(t+\Delta t)$$

$$\vec v_{corr}(t+\Delta t) = \vec v_{pred} (t+\Delta t) + a_1 \Delta \vec v(t+\Delta t)$$

$$\vec a_{corr}(t+\Delta t) = \vec a_{pred} (t+\Delta t) + a_2 \Delta \vec a(t+\Delta t)$$

$$etc ...$$ (2.14)

Ainsi, pour chaque composante du développement nous pouvons déduire une correction à ajouter à la prédiction. Cette correction est donnée par un ensemble de coefficients trouvés empiriquement pour avoir la meilleure stabilité de l'algorithme. Dans notre cas nous avons [29,30]:

$$a_0 = 3/16$$

$$a_1 = 251/360$$

$$a_2 = 1$$

$$a_3 = 11/18$$

$$a_4 = 1/6$$

$$a_5 = 1/60$$

L'algorithme de prédicteur-correcteur de Gear consiste donc à évaluer la position au pas de temps suivant, de la même façon qu'avec l'algorithme de Verlet, et de corriger l'erreur commise ( $F_{ij}(t) \neq F_{ij}(t+\Delta t)$ ) par des pondérations visant à optimiser la stabilité énergétique. Nous pouvons voir que cet algorithme est beaucoup plus précis que l'algorithme de Verlet (Fig 2.4).

Figure 2.4: L'erreur pour un algorithme de Gear est bien moins importante que celle d'un algorithme de Verlet. Le pas de temps ainsi que la précision sont en unités réduites.