Équation de Kepler

En astronomie, l'équation de Kepler est une formule liant, dans une orbite, l'excentricité e et l'anomalie excentrique E à l'anomalie moyenne M. L'importance de cette équation est qu'elle permet de passer des paramètres dynamiques du mouvement d'un astre (l'anomalie moyenne) aux paramètres géométriques (l'anomalie excentrique). Cette équation a été établie par Kepler dans le cas des orbites elliptiques, en analysant les relevés de position de la planète Mars effectués par Tycho Brahe. Elle fut ensuite généralisée aux autres formes d'orbites (paraboliques, hyperboliques, quasi-paraboliques, rectilinéaires) à l'aide des principes de la mécanique newtonienne.

Présentations de l'équation de Kepler

L'équation de Kepler en tant que telle est celle établie par Kepler pour les orbites elliptiques. Elle peut cependant être déclinée en plusieurs formes pour couvrir tous les cas d'orbites.

Cas de l'orbite elliptique

L'équation de Kepler en orbite elliptique est :

$E-e\cdot \sin (E)=M$

avec l'anomalie moyenne Modèle:Mvar définie par :

$M=n\cdot (t-t_0)$

avec Modèle:Mvar le moyen mouvement :

$n=\frac{2\cdot \pi }{T}$

Modèle:Mvar le temps et Modèle:Math étant l'instant du passage au périastre. Modèle:Mvar est la période orbitale.

Le moyen mouvement peut être aussi exprimé par :

$n=\sqrt{\frac{G\cdot (m_1+m_2)}{a^3}}$

où

• Modèle:Mvar est la constante universelle de la gravitation
• Modèle:Math et Modèle:Math les masses des deux corps
• Modèle:Mvar le demi-grand axe de l'ellipse $a=\frac{p}{1-e^2}$ ; Modèle:Mvar étant le paramètre de l'ellipse

L'équation de Kepler, associée au lien entre l'anomalie excentrique Modèle:Mvar et l'anomalie vraie Modèle:Mvar

$\tan (v/2)=\sqrt{\frac{1+e}{1-e}}\cdot \tan (E/2)$

permet de déterminer la position au cours du temps d'un astre sur son orbite.

Cas de l'orbite hyperbolique

En cas d'orbite hyperbolique (Modèle:Math), on peut démontrer analytiquement une relation équivalente à l'équation de Kepler :

$e\cdot \sinh (H)-H=M$

où Modèle:Math désigne le sinus hyperbolique.

M est défini de la même manière que dans le cas elliptique, avec l'expression du moyen mouvement suivant :

$n=\sqrt{\frac{G\cdot (m_1+m_2)}{(-a{)}^3}}$

L'argument H n'est plus un angle comme c'est le cas de E dans le mouvement elliptique. H est dans ce cas liée à l'anomalie vraie v par :

$\tan (v/2)=\sqrt{\frac{e+1}{e-1}}\cdot \tanh (H/2)$

Cas de l'orbite parabolique

L'équation de Kepler n'est pas définie dans le cadre du mouvement parabolique (Modèle:Math). Elle est remplacée par l'équation de Barker.

$M=s+\frac{s^3}{3}$

avec

$s=\tan (v/2)$

$M=n\cdot (t-t_0)$ et

$n=2\cdot \sqrt{\frac{G\cdot (m_1+m_2)}{p^3}}$

Cette équation cubique peut être résolue de manière analytique par la méthode de Cardan.

Expression universelle de l'équation de Kepler

Par un changement de variable, les équations de Kepler elliptiques, paraboliques et hyperboliques peuvent être regroupées en une seule équation "universelle"Modèle:,. Une des expressions possible est :

$qx+e{x}^3{c}_3(\alpha x^2)=\sqrt{G\cdot (m_1+m_2)}(t-t_0)$

avec le périastre Modèle:Math, et Modèle:Math. Modèle:Mvar est positif pour les orbites elliptiques, nul pour les orbites paraboliques et négatif pour celles hyperboliques. La nouvelle variable Modèle:Mvar est définie par :

$x=E\cdot \sqrt{a}$

et la fonction Modèle:Math est une des fonctions de Stumpff, qui s'écrit dans le cas général :

$c_k(t)=\sum _{n=0}^{+\mathrm{\infty }}\frac{(-1{)}^n}{(2n+k)!}{t}^n$

La détermination de Modèle:Mvar d'après l'équation universelle permet de déterminer la position du corps sur son orbite Modèle:Math par :

$X=a(\cos E-e)=q+a\cdot \cos \left(\frac{x}{\sqrt{a}}\right)=q-x^2\sum _{n=0}^{+\mathrm{\infty }}\frac{(-1{)}^n}{(2n+2)!}{(\alpha x^2)}^n=q-x^2{c}_2(\alpha x^2)$

$Y=a\sqrt{1-e^2}\sin E=\sqrt{a(1-e^2)}\sqrt{a}\sin \left(\frac{x}{\sqrt{a}}\right)=x\cdot \sqrt{q(1+e)}\sum _{n=0}^{+\mathrm{\infty }}\frac{(-1{)}^n}{(2n+1)!}{(\alpha x^2)}^n=x\cdot \sqrt{q(1+e)}{c}_1(\alpha x^2)$

Les fonctions Modèle:Math et Modèle:Math étant définies de la même manière que Modèle:Math plus haut.

Orbites rectilinéaires

Les orbites rectilinéaires sont des cas limites des autres orbites, en faisant tendre la distance au périastre Modèle:Mvar vers zéro tout en conservant le demi grand axe Modèle:Mvar constant : l'orbite tend alors vers un segment ou une demi-droite. Pour le cas des orbites elliptiques et hyperboliques, cela suppose de faire tendre l'excentricité Modèle:Mvar vers 1, car demi grand axe Modèle:Mvar, excentricité Modèle:Mvar et périastre Modèle:Mvar sont liés par Modèle:Math. Il existe donc trois types d'orbites rectilinéaires : elliptiques, paraboliques et hyperboliques. En pratique, une partie seulement de ces orbites est décrite par l'astre, aboutissant soit à une collision, soit à une évasion. Certaines comètes kamikazes détectées par les observatoires solaires spatiaux (SoHO, SDO…) ou bien des tronçons d'orbites de sondes interplanétaires sont proches des orbites rectilinéaires.

Pour l'orbite rectilinéaire elliptique, l'équation de Kepler devient :

$E-\sin (E)=M$

avec l'anomalie moyenne Modèle:Mvar définie par :

$M=n\cdot (t-t_0)$

L'anomalie vraie n'ayant plus de signification pour une orbite rectilinéaire, la position de l'astre est définie par sa distance le séparant de l'astre principal r:

$r=a(1-\cos (E))$

Pour l'orbite rectilinéaire hyperbolique, l'équation de Kepler devient :

$\sinh (H)-H=M$

et la position de l'astre :

$r=a(1-\cosh (H))$

Modèle:Mvar étant négatif pour les orbites hyperboliques

Enfin, pour l'orbite rectilinéaire parabolique :

$M=\frac{s^3}{6}$

avec

$M=n\cdot (t-t_0)$ et

$n=\sqrt{G\cdot (m_1+m_2)}$

et la position de l'astre :

$r=\frac{s^2}{2}$

Résolution de l'équation de Kepler

L'équation de Kepler

$E-e\cdot \sin (E)=M$

permet de calculer de manière directe la date (liée à Modèle:Mvar) correspondant à une position donnée (liée à Modèle:Mvar), par exemple déterminer la date des équinoxes. Par contre le problème inverse, déterminer la position d'une planète pour une date donnée, nécessite la détermination de Modèle:Mvar, connaissant Modèle:Mvar et Modèle:Mvar. Ce problème ne peut pas être résolu de manière simple.

Résoudre l'équation de Kepler, c'est trouver Modèle:Math :

• comme série de Fourier puisque c'est une fonction périodique impaire de Modèle:Mvar
• comme série entière de Modèle:Mvar, si Modèle:Math, rayon de convergence de la série.
• comme une valeur numérique avec un nombre de chiffres (Modèle:Mvar), pour un temps de calcul Modèle:Math optimisé.

Série de Fourier

C'est LagrangeModèle:, qui trouve l'expression, bien que le nom Modèle:Math soit associé au nom de Bessel.

• Modèle:Math est une fonction impaire périodique de Modèle:Mvar :

$E-M=e\cdot \sin E=2\cdot \sum _{n=1}^{+\mathrm{\infty }}\frac{J_n(ne)}{n}\sin (nM)$

où Modèle:Math est la fonction de Bessel de Modèle:1re d'ordre Modèle:Mvar.

Série entière de l'excentricité

C'est encore Lagrange qui trouve la solution, que Laplace complétera en donnant le rayon de convergence. Ces travaux inspireront Cauchy, qui fondera la théorie des séries analytiques pour résoudre ce problème épineux ; celui-ci verra son aboutissement avec les travaux de Puiseux.

L'application du théorème d'inversion de série de Lagrange fournit :

$E-M=\sum _{n=1}^{+\mathrm{\infty }}\frac{e^n}{n!}\cdot a_n(M)$

avec

$a_n(M)=\frac{{\mathrm{d}}^{n-1}}{\mathrm{d}{M}^{n-1}}({\sin }^{n}M).$

Le rayon de convergence minimum de la série, qui dépend de Modèle:Mvar, est atteint pour Modèle:Math, et vaut Modèle:Math tel qu'indiqué par Laplace (1823) et démontré par Cauchy et Puiseux :

$e_0=2\frac{\sqrt{x(1-x)}}{1-2x}$ et Modèle:Mvar tel que $\frac{1-x}{x}=\exp \left(\frac{2}{1-2x}\right)$.

Ceci rend cette formule inapplicable pour déterminer la position des comètes, dont l'excentricité est souvent voisine de 1.

Les premiers termes sont :

$E-M=e\cdot \sin M+e^2\cdot \sin M\cdot \cos M+e^3\cdot (\sin M\cdot {\cos }^{2}M-\frac{{\sin }^{3}M}{2})+\dots$

Note : il est possible d'obtenir ce développement en série en remplaçant dans la série de Fourier précédente les fonctions de Bessel par leur développement limité :

$J_n(x)=(x/2{)}^n\sum _{p=0}^{\mathrm{\infty }}\frac{(-1{)}^p}{p!(n+p)!}(x/2{)}^{2p}$

On obtient alors le développement limité beaucoup plus simplement que par la méthode de l'inversion de série :

$E-M=$

$(e-\frac{e^3}{8}+\frac{e^5}{192}+\cdots )\cdot \sin M+(\frac{e^2}{2}-\frac{e^4}{6}+\frac{e^6}{48}+\cdots )\cdot \sin 2M+(\frac{3}{8}{e}^3-\frac{27}{128}{e}^5+\cdots )\cdot \sin 3M$

$+(\frac{e^4}{3}-\frac{4}{15}{e}^6+\cdots )\cdot \sin 4M+(\frac{125}{384}{e}^5+\cdots )\cdot \sin 5M+(\frac{27}{80}{e}^6+\cdots )\cdot \sin 6M$

Il est à noter que bien que la série de Fourier converge pour Modèle:Math, et que les développements des fonctions de Bessel aient un rayon de convergence infini, le résultat après réorganisation des termes ne converge que pour Modèle:Math

Cas des comètes : Modèle:Math

Le premier à se confronter au problème est Horrocks, puis surtout Halley, pour les calculs sur sa comète d'excentricité Modèle:Math.

Plusieurs solutions ont été proposées en modifiant légèrement l'équation de Barker (Modèle:Math). La solution proposée par Bessel (1805) couvre le domaine Modèle:Math. Gauss s'illustra en donnant une belle solution pour Modèle:Math.

Une généralisation de l'équation de Barker est un développement en série convergeant d'autant plus rapidement que l'excentricité Modèle:Mvar est proche de 1, ce qui s'avère bien adapté aux cas des comètes (cette série s'applique aussi aux orbites légèrement hyperboliques) :

$n(t-t_0)=S+(1-2\gamma )\frac{S^3}{3}-\gamma (2-3\gamma )\frac{S^5}{5}+{\gamma }^2(3-4\gamma )\frac{S^7}{7}-\dots =S+\sum _{p=1}^{\mathrm{\infty }}(-\gamma {)}^{p-1}(p-(p+1)\gamma )\frac{S^{2p+1}}{2p+1}$

qui a pour rayon de convergence : $|\gamma |S^2<1$

avec Modèle:Math

$n=\frac{k}{2}\frac{\sqrt{1+e}}{q^{3/2}}$

$\gamma =\frac{1-e}{1+e}$

Modèle:Math étant l'anomalie vraie, Modèle:Math la constante gravitationnelle de Gauss, Modèle:Math et Modèle:Math étant respectivement l'excentricité et le périastre de l'orbite, Modèle:Math le temps et Modèle:Math étant l'instant du passage au périastre.

Lorsque Modèle:Math, la série se réduit à l'équation de Barker.

Calcul numérique

L'équation de Kepler peut être résolue à l'aide d'un algorithme de recherche d'un zéro d'une fonction. Les méthodes du type encadrement, méthode de dichotomie, méthode de la fausse position nécessitent un encadrement de départ dans lequel la racine est présente. Du fait de la périodicité et la parité de l'équation de Kepler, il est toujours possible de ramener l'intervalle de départ à Modèle:Math. Ceci fournit un encadrement de départ pour ces méthodes, mais il est aisé d'en trouver de plus fins.

Les méthodes du type point fixe nécessitent une estimation de départ de la racine, le germe de la méthode Modèle:Math, pour lancer les calculs : il en existe de multiples dans la littérature, dont la plus simple est Modèle:Math.

La méthode de type point fixe la plus simple, celle utilisée par Kepler, est :

$E_0=M$

$E_{n+1}=M+e\cdot \sin E_n\text{ pour }n=0,1\dots$

converge lentement lorsque Modèle:Mvar est proche de 1. Il est alors avantageux de lui adjoindre un algorithme d'accélération de la convergence : le Delta-2 de Aitken, par exemple, ou la variante de Steffensen.

L'équation de Kepler se prête particulièrement bien aux algorithmes nécessitant le calcul des dérivées successives élevées, du fait du faible coût en calcul machine nécessaire. En effet :

$f(E)=E-e\cdot \sin E-M$

$f^{\prime }(E)=1-e\cdot \cos E$

$f^{″}(E)=e\cdot \sin E$

$f^{‴}(E)=e\cdot \cos E$

Les dérivées suivantes se déduisant cycliquement des précédentes. Les variantes de la méthode de Newton et de Halley d'ordres plus élevés sont donc très performantes dans ce cas. Il est à noter que ces méthodes peuvent dans certains cas avoir des difficultés à converger (Modèle:Mvar proche de 1 et Modèle:Mvar proche de 0). Il est préférable dans ces zones soit de proposer une valeur de départ moins grossière (germe de Mikkola (Seppo Mikkola) ou de Markley), soit de brider les méthodes itératives pour les forcer à converger (modification de Hamming de la méthode de Newton), ou d'utiliser des méthodes itératives à convergence moins locale (méthode de Laguerre).

Exemple

Lors de son dernier passage en 1986, la comète de Halley a eu la visite de la sonde Giotto. Les données nécessaires pour déterminer la position de la comète lors de cette rencontre sont:

• date du passage au périhélie Modèle:Math : Modèle:Date- à 10:59:55 UTC
• date de la rencontre Modèle:Mvar : Modèle:Date- à 00:03:00 UTC, soit Modèle:Math = 32,54328 jours
• excentricité de l'orbite Modèle:Mvar : 0,96727426
• distance au périhélie Modèle:Mvar : 0,58710224, soit le demi grand axe Modèle:Math et le moyen mouvement Modèle:Math.

L'anomalie moyenne vaut Modèle:Math

L'équation de Kepler à résoudre est :

$E-0,96727426\sin (E)=0,0073673887$

En partant de Modèle:Math et en utilisant la méthode de Newton,

Graphique représentant le nombre d’itérations (de noir=0 à blanc=255 ou plus) nécessaires à la méthode de Newton pour résoudre l'équation de Kepler avec une valeur initiale égale à Modèle:Mvar. Pour chaque couple Modèle:Math, la couleur du point indique le nombre d'itérations

$E_{n+1}=E_n-\frac{E_n-e\cdot \sin E_n-M}{1-e\cdot \cos E_n}$

on trouve successivement :

0,0073673887

0,2249486948

0,1929911041

0,1909186907

0,1909107985

0,1909107984

… (les valeurs suivantes étant identiques) On en déduit l'angle de position de la comète sur son orbite (l'anomalie vraie) Modèle:Math

La distance de la comète au soleil se calcule par Modèle:Math = 0,902374257 UA (légèrement inférieure à la distance terre-soleil)

La vitesse de la comète vaut $v=42.1219\sqrt{1/r-1/(2a)}$ soit Modèle:Unité

Les itérations ne se passent pas toujours aussi bien pour les comètes, comme l'expose le graphique ci contre. Pour des excentricités au-delà de 0.97, la convergence est incertaine avec comme point de départ des itérations Modèle:Math. D'autres points de départ plus précis permettent d'éviter cet écueil.

Dans le cas des comètes, la résolution de la généralisation quasi-parabolique de l'équation de Barker pose deux problèmes :

• le calcul approché de la série qui peut nécessiter un grand nombre de termes, voire être impossible si elle diverge. Il se trouve que cette série se prête particulièrement bien à l'utilisation d'algorithmes d'accélération de la convergence, notamment le Δ² de Aitken ou l'ε-algorithme de Peter Wynn, qui non seulement accélèrent la convergence, mais étendent son domaine de convergence. En pratique, les difficultés interviennent lorsque la comète est très loin de son périastre (elle est alors depuis longtemps invisible), ou son excentricité diffère notablement de 1 (dans ce cas, il est plus judicieux de résoudre l'équation de Kepler elliptique ou hyperbolique).
• La résolution de l'équation proprement dite. Celle-ci peut être effectuée par les méthodes de type Newton avec :

$F(S)=-n(t-t_0)+S+\sum _{p=1}^{\mathrm{\infty }}(-\gamma {)}^{p-1}\cdot (p-(p+1)\cdot \gamma )\cdot \frac{S^{2p+1}}{2p+1}$

en remarquant que la dérivée s'exprime simplement :

$F^{\prime }(S)=\frac{1+S^2}{{(1+\gamma \cdot S^2)}^2}$

les dérivées suivantes s'en déduisent facilement.

On pourra choisir comme valeur initiale de l'itération Modèle:Math, la solution de l'équation cubique obtenue en retenant les premiers termes (différant légèrement de l'équation de Barker), à l'aide de la méthode de Cardan

$n(t-t_0)=S_0+(1-2\cdot \gamma )\cdot \frac{{S_0}^3}{3}.$

Recherches actuelles

Modèle:Section à recycler Les calculs via les intégrateurs symplectiques exigent de rester toujours en butée du nombre de décimales, dans le moindre coût de calcul. Elle dépend beaucoup du doublet Modèle:Math, Modèle:Mvar compris entre 0 et Modèle:Math et de Modèle:Mvar, surtout quand ce dernier paramètre est voisin de 1.

Modèle:Harvsp adopte la méthode de Modèle:Harvsp qui est la méthode de Newton d'ordre 4, en choisissant « adéquatement » le germe Modèle:Math en fonction du doublet Modèle:Math.

Il est clair que dans les calculs numériques, le volume de calculs est essentiel, autant que le nombre de décimales, vu l'instabilité du système solaire évaluée à un coefficient de Liapunov de Modèle:Math. On se heurte à une muraille exponentielle : difficile d'aller plus loin que 25 Myr, même avec un traitement 128 bits.

Ce sont ces calculs (astronomiques... mais informatisés) qui tournent sur les machines de l'IMCCE-Paris. Le calcul de l'ensoleillement terrestre à la latitude 65° Nord, I(65,t) est calculé et on essaie d'en déduire la corrélation avec le climat passé : l'échelle géologique jusqu'au néogène (25 M ans) en est déduite (échelle géologique Gradstein 2004). Prochaine étape prévue : les 65 M ans.

Histoire des sciences

Avant Kepler, l'équation était déjà étudiée pour d'autres motifs :

c'est le problème de la réduction des coordonnées locales aux coordonnées géocentriques : il faut réduire la correction de parallaxe. Habash al Hasib s'y est déjà attaqué.

Avant 1700, il y a déjà beaucoup de tentatives : Kepler naturellement, Curtz (1626), Niele, Boulliau (1645, 1657), Seth Ward (1653), Paganus (1657), Horrebow (1717), Cassini (1669), Newton (1665?), Wren (1658), Wallis (1659), Jeremiah Horrocks (1638)...

Bibliographie

• Modèle:Ouvrage
• Modèle:En John Brinkley, Trans Roy Irish Ac, vol. 7, 1803, Modèle:P.
• Modèle:Ouvrage
• Modèle:Article
• Modèle:Article

Voir aussi

Articles connexes

• Équation du centre
• Mouvement képlérien
• Astronomia nova
• Fonction de Stumpff

Liens externes

• Modèle:Lien web.
• Modèle:En Série d'articles de Danby et Burkardt[1][2]
• Modèle:En Algorithme et programme de calcul par Odell et Gooding [3]
• Modèle:En La méthode de Gauss pour la résolution de l'équation de Kepler
• Modèle:Autorité
• Modèle:Dictionnaires