Astronef en rotation entre terre et lune

(Cliquez-ici pour accéder à la version originale de cette discussion avec couleurs et images)





Posted by: vinch

Bonjour
je programme en ce moment un petit programme censé simuler grossièrement la trajectoire d'un astronef qu'on aurait lancé de la terre vers la lune.
Mon programme est créé, il tourne mais malheuresement le résultat diverge, j'ai très vite des résultats énormes ...
Je ne crois pas qu'il y ait d'erreurs dans mon code mais je soupçonne qu'il y en ait dans mes équations.
Si quelqu'un à le temps de jeter un coup d'oeil et de me dire ce qu'il en pense je lui serai grandement reconnaissant !
J'utilise le formalisme lagrangien et hamiltonien .
allons y :
Nous nous plaçons dans le référentiel de la terre, la trajectoire de la lune est suposée circulaire, toutes les coordonnées sont données dans le système polaire.
Nous notons :
mt masse terre
ml masse lune
rl rayon de l'orbite de la lune
r la distance astronef-terre
rp ou toute autre lettre données avec "p" veut dire sa dérivée par rapport au temps (point)
te = teta l'angle de la position de l'astro avec une origine (l'horizontale par exemple)
tel= la meme chose pour la lune

L'énergie cinétique de l'astronef est T=m/2*(rp^2+r^2*tep^2)
les énergies potentielles dues à la lune et à la terre :
Vl=-G*m*ml/|AL| et Vt= -Gm*mt/r
avec AL la distance lune astronef que je trouve égale à :
|AL|=racine(r^2+rl^2-2*r*rl*cos(te-tel))

le lagrangien de ce système s'écrit donc :
L=T-V=m/2*(rp^2+r^2*tep^2)+G*m*(ml/|AL|+mt/r)

j'en déduis deux moments conjugués :
Pr=m*rp et Pte= m*r^2*tep^2

Je trouve l'hamiltonien du système :
H= Pr^2/(2*m)+Pte^2/(2*m*r^2)-Gm(mt/r^2+ml/|AL|)

et grâces aux équations canoniques de hamilton j'en déduis :

Prp= Pte^2/(m*r) -Gm*(mt/r^2 + ml*(r-rl*cos(te-tel))/AL^3)

et

Ptep=-G*m*ml*(r*rl*sin(tel-te))/|AL|^3

Qui sont les deux dérivées par rapport au temps de mes moments conjugués.
Grace à ces deux équations et la méthode d'euler (ou runge kutta) = y=y+yp*dt je peux trouver les valeurs de Pr et Pte d'un instant t à l'instant t+dt

Dans mon programme c'est l'expression de Prp qui prends des valeurs énormes ... est ce normal ? la valeur que j'ai trouvé est elle correcte ?


Posted by: Dominique Lefebvre

Bonsoir,

Ta question mérite réflexion... Tu essais de programmer un problème à trois corps, qui n'a pas de solution analytique.
Sur le plan numérique, les méthodes d'intégration sont délicates, car souvent on observe une divergence du schéma, ce qui semble être ton cas. A première vue, je dirais que ta solution pêche à la fois sur le plan algorithmique et numérique. Il est certain qu'Euler ne convient pas et j'ai des doutes sur RK4.
Je regarderai cela de plus près quand j'aurai un moment.


Posted by: vinch

ça devrait marcher au moins pour les des temps petits et proches de l'origine (je ne me fais pas d'illusion, apres 1000 intégrations numériques les données obtenues n'auront plus grand chose à voir avec "la réalité")

il est possible d'avoir un bon petit modèle, mon prof l'a fait (mais je tiens à ne pas regarder son code !)

pourquoi tu penses que RK4 n'est pas utilisable dans un cas comme ça ? pour intégrer des equations différencielles d'ordre 1 ça me parait une méthode puissante non ?


Posted by: Dominique Lefebvre

Citation:
Posté par vinch
ça devrait marcher au moins pour les des temps petits et proches de l'origine (je ne me fais pas d'illusion, apres 1000 intégrations numériques les données obtenues n'auront plus grand chose à voir avec "la réalité")

il est possible d'avoir un bon petit modèle, mon prof l'a fait (mais je tiens à ne pas regarder son code !)

pourquoi tu penses que RK4 n'est pas utilisable dans un cas comme ça ? pour intégrer des equations différencielles d'ordre 1 ça me parait une méthode puissante non ?

Bonjour,
Tout dépend de la précision que tu souhaites dans ta simulation.
Si tu fais trop d'approximation au niveau du modèle, ce dernier n'aura plus aucun sens physique.
Si tu prends en compte le modèle à trois corps (celui habituellement utilisé dans ce cas), alors le schéma RK4 présente des défauts de convergence qui rendront tes résultats instables.

Ceci dit, tu dis vouloir le faire marcher pour des temps petits, proches de l'origine. Dans ce cas, je ne vois pas comment tu vas pouvoir simuler la trajectoire de ton vaisseau entre la Terre et la Lune?


Posted by: vinch

Dans mes equations le temps est donné en heures, des temps petits ça peut vouloir dire environ 10h (si tu fais 1000 itérations avec un pas de 0.01h)
normalement les équations ne divergent pas trop, de toute façon près de la terre les équations ressemblent fortement au problème à deux corps (l'astronef fait environ 1 tone l'attraction de la lune sur lui est assez faible au début)

Mais mes équations sont juste d'apres toi ?
Parce que là je viens de recopier mon problème dans scilab (alors que mon programme est en C++) pour m'affranchir de toute possibilité d'erreur de programmation et pouvoir tester plein de choses ... et j'ai tjrs le meme résultat, Prp la dérivée de la qtt de mouvement suivant r est trop grande !


Posted by: Dominique Lefebvre

Je viens de relire ton problème et je me pose un certain nombre de questions:

1 - l'hypothèse de base est qu'on lance l'astronef avec une vitesse v0 supérieure à la vitesse de libération et que cet astronef n'est pas propulsé. Ces deux points ne ressortent pas clairement. Sont-ils exacts?

2 - tu supposes le mouvement de la lune autour de la terre circulaire, avec une vitesse de révolution de la lune autour de la terre constante ? Est-ce exact?

3 - tu sembles travailler dans le plan orbital de la lune, le mouvement ayant lieu dans ce plan, sans que cela soit bien clair. Est-ce exact?

4 - l'influence des autres masses (le soleil, jupiter, saturne) est négligée.

Si la réponse est oui à toutes ces questions, le mouvement de l'astronef est défini par deux équations de deuxième ordre qui sont :
Ax = AxT + AxL et Ay = AyT + AyL A étant l'accélaration de l'astronef, l'indice T désigne la composante de l'accélaration terrestre et l'indice L désigne la composante de l'accélaration lunaire.

L'intégration numérique de des deux équations (par RK4) te donne la trajectoire dans le plan orbital lunaire de l'astronef. Il s'agit d'une trajectoire d'interception dans la mesure où l'astronef n'est pas propulsé.


Posted by: vinch

Les réponses à tes 4 hypothèses sont oui.

Par contre si j'utilise le formalisme lagrangien ou même l'équation de Newton, je trouve 2 équations différentielles d'ordre 2, je suis d'accord.
Mais l'idée qu'on avait c'était de passer par les équations canoniques de Hamilton, à ce moment là j'ai non plus 2 équations mais 4 (rp, tep, Prp, Ptep) mais seulement d'ordre 1.
J'aurai moins de mal à les intégrer grace à RK4.


Posted by: Dominique Lefebvre

Citation:
Posté par vinch
Les réponses à tes 4 hypothèses sont oui.

Par contre si j'utilise le formalisme lagrangien ou même l'équation de Newton, je trouve 2 équations différentielles d'ordre 2, je suis d'accord.
Mais l'idée qu'on avait c'était de passer par les équations canoniques de Hamilton, à ce moment là j'ai non plus 2 équations mais 4 (rp, tep, Prp, Ptep) mais seulement d'ordre 1.
J'aurai moins de mal à les intégrer grace à RK4.


cela revient strictement au même... Si tu veux intégrer une EDO d'ordre 2 avec RK4, tu la décompose en 2 EDO d'ordre 1!
Et je trouve personnelement que le formalisme Hamiltonien est un peu lourd pour ce problème.
Par contre, je ne vois pas apparaître la vitesse initiale de ton astronef dans tes équations, j'ai mal lu?


Posted by: vinch

c'est que j'ai écrit que 2 des 4 équations, je les mets toutes ici :

rp=Pr/m
tep=Pte/(m*r^2)

Ptep=Gm*ml*(r-rl*sin(tel-te))
Prp=Pte^2/(m*r^2)-G*m*(mt/r^2+ml*r-rl*cos(te-tel)/|AL|^3)

donc les dérivées de r et teta dépendent de leurs quantités de mouvement donc si ils existe une vitesse initiale, il y a une quantité de mouvement initiale, c'est à cet endroi qu'elle est prise en compte.


Posted by: Dominique Lefebvre

Citation:
Posté par vinch
c'est que j'ai écrit que 2 des 4 équations, je les mets toutes ici :

rp=Pr/m
tep=Pte/(m*r^2)

Ptep=Gm*ml*(r-rl*sin(tel-te))
Prp=Pte^2/(m*r^2)-G*m*(mt/r^2+ml*r-rl*cos(te-tel)/|AL|^3)

donc les dérivées de r et teta dépendent de leurs quantités de mouvement donc si ils existe une vitesse initiale, il y a une quantité de mouvement initiale, c'est à cet endroi qu'elle est prise en compte.


OK, bien que ce ne soit pas très lisible...


Posted by: vinch

et donc elles ont l'air correcte ces equations ?
je vois vraiment pas ce qui peut clocher ...










-