252
7 Equations différentielles ordinaires
Programme 7.10. threebody : second membre pour le système du problème à
trois corps
function f= t h r e ebod y(t ,y )
[n , m ]= size ( y ); f = zeros (n ,m ); Ms =330000; Me =1; Mm =0.1;
D1 = (( y (5) -y (1))^2+( y (7) -y ( 3 ) ) ^2 )^(3 /2) ;
D2 = (y (1)^2+ y ( 3 ) ^ 2) ^(3/ 2);
f (1)= y (2); f (2)=4* pi ^2*( Me / Ms *(y (5) -y (1))/ D1 - y (1)/ D2 );
f (3)= y (4); f (4)=4* pi ^2*( Me / Ms *(y (7) -y (3))/ D1 - y (3)/ D2 );
D2 = (y (5)^2+ y ( 7 ) ^ 2) ^(3/ 2);
f (5)= y (6); f (6)=4* pi ^2*( Mm / Ms *(y (1) -y (5))/ D1 - y (5)/ D2 );
f (7)= y (8); f (8)=4* pi ^2*( Mm / Ms *(y (3) -y (7))/ D1 - y (7)/ D2 );
return
Comparons la méthode de Crank-Nicolson (implicite) et la méthode
adaptative de Runge-Kutta implémentée dans ode23 (explicite). En normalisant la distance Terre-Soleil à 1, la distance Soleil-Mars vaut 1.52 :
on prend donc (1, 0) pour la position initiale de Terre et (1.52, 0) pour
celle de Mars. Supposons de plus que la vitesse horizontale des deux planètes est nulle, que la vitesse verticale de la Terre est égale à −5.1 et que
celle de Mars vaut −4.6 (en unités adimensionnées) : avec ce choix, elles
devraient conserver une orbite relativement stable autour du Soleil. On
choisit 2000 pas de discrétisation pour la méthode de Crank-Nicolson :
[ t23 , u23 ]= ode23 ( @threebody ,[0 10] ,...
[1.52 0 0 -4.6 1 0 0 -5.1]);
[ tcn , ucn ]= cranknic ( @threebody ,[0 10] ,...
[1.52 0 0 -4.6 1 0 0 -5.1] ,2000);
Les graphes de la Figure 7.21 montrent que les deux méthodes reproduisent convenablement les orbites périodiques des deux planètes autour
du Soleil. La méthode ode23 ne nécessite que 543 itérations (avec des
pas non uniformes) pour construire une solution plus précise que celle
−1.5
−1
−0.5
0
0.5
1
1.5
2
−1.5
−1
−0.5
0
0.5
1
1.5
S
−1.5
−1
−0.5
0
0.5
1
1.5
2
−1.5
−1
−0.5
0
0.5
1
1.5
S
Figure 7.21. Les orbites de la Terre (la plus petite) et de Mars autour du
Soleil calculées avec la méthode adaptative ode23 (à gauche) en 543 itérations
et avec la méthode de Crank-Nicolson (à droite) en 2000 itérations
7 Equations différentielles ordinaires
Programme 7.10. threebody : second membre pour le système du problème à
trois corps
function f= t h r e ebod y(t ,y )
[n , m ]= size ( y ); f = zeros (n ,m ); Ms =330000; Me =1; Mm =0.1;
D1 = (( y (5) -y (1))^2+( y (7) -y ( 3 ) ) ^2 )^(3 /2) ;
D2 = (y (1)^2+ y ( 3 ) ^ 2) ^(3/ 2);
f (1)= y (2); f (2)=4* pi ^2*( Me / Ms *(y (5) -y (1))/ D1 - y (1)/ D2 );
f (3)= y (4); f (4)=4* pi ^2*( Me / Ms *(y (7) -y (3))/ D1 - y (3)/ D2 );
D2 = (y (5)^2+ y ( 7 ) ^ 2) ^(3/ 2);
f (5)= y (6); f (6)=4* pi ^2*( Mm / Ms *(y (1) -y (5))/ D1 - y (5)/ D2 );
f (7)= y (8); f (8)=4* pi ^2*( Mm / Ms *(y (3) -y (7))/ D1 - y (7)/ D2 );
return
Comparons la méthode de Crank-Nicolson (implicite) et la méthode
adaptative de Runge-Kutta implémentée dans ode23 (explicite). En normalisant la distance Terre-Soleil à 1, la distance Soleil-Mars vaut 1.52 :
on prend donc (1, 0) pour la position initiale de Terre et (1.52, 0) pour
celle de Mars. Supposons de plus que la vitesse horizontale des deux planètes est nulle, que la vitesse verticale de la Terre est égale à −5.1 et que
celle de Mars vaut −4.6 (en unités adimensionnées) : avec ce choix, elles
devraient conserver une orbite relativement stable autour du Soleil. On
choisit 2000 pas de discrétisation pour la méthode de Crank-Nicolson :
[ t23 , u23 ]= ode23 ( @threebody ,[0 10] ,...
[1.52 0 0 -4.6 1 0 0 -5.1]);
[ tcn , ucn ]= cranknic ( @threebody ,[0 10] ,...
[1.52 0 0 -4.6 1 0 0 -5.1] ,2000);
Les graphes de la Figure 7.21 montrent que les deux méthodes reproduisent convenablement les orbites périodiques des deux planètes autour
du Soleil. La méthode ode23 ne nécessite que 543 itérations (avec des
pas non uniformes) pour construire une solution plus précise que celle
−1.5
−1
−0.5
0
0.5
1
1.5
2
−1.5
−1
−0.5
0
0.5
1
1.5
S
−1.5
−1
−0.5
0
0.5
1
1.5
2
−1.5
−1
−0.5
0
0.5
1
1.5
S
Figure 7.21. Les orbites de la Terre (la plus petite) et de Mars autour du
Soleil calculées avec la méthode adaptative ode23 (à gauche) en 543 itérations
et avec la méthode de Crank-Nicolson (à droite) en 2000 itérations