Bonjour à tous,
J'ai un projet à faire pour la rentré dont l'intitulé est symple : à partir de Ro=4.223e7m et Vo=3071m/s, modeliser la trajectoire d'un satellite. Le but est d'utiliser la méthode "Runge et Kutta" sans l'aide d'un programme MATLAB et la comparer en terme de gain avec ode45.
J'ai fais un programme que vous trouverez ci-dessous.
Mon probleme est que j'ai du mal a le faire marcher mais m'as question est surtout comment modéliser une trajectoire élliptique à partir de cette méthode? Et aussi est-ce la bonne programation pour une équa diff 2de ordre avec X(t) et Y(t).
Merci beaucoup pour votre aide.
h=input('entrez le pas désiré');
j=input('entrez le nombre de jour d étude désiré');
to=0 ; tM =24*3600*j ;
t=to:h:tM ; h2=h/2 ;
N=length(t);
G=6.67e-11 ; MT=5.97e24 ; K=MT*G;
X=zeros(1,N) ; Y=zeros(1,N) ; Xt=zeros(1,N) ; Yt=zeros(1,N);
X(1)=4.223*10^7 ; Xt(1)=0 ; Y(1)=0 ; Yt(1)=3071 ;
global K X Y;
% résolution à l'aide de RK
tic
for k=1:N-1;
t12=t(k)+h2;
j1=h*Xt;
k1=h*feval(f,t(k),X(k));
j2=h*(Xt+k1/2);
k2=h*feval(f,t12,X(k)+0.5*k1);
j3=(Xt+k2/2)*h;
k3=h*feval(f,t12,X(k)+0.5*k2);
j4=(Xt+k3)*h;
k4=h*feval(f,t(k+1),X(k)+k3);
Xt(k+1)=Xt(k)+(j1+2*j2+2*j3+j4 )/6;
X(k+1)=X(k)+(k1+2*k2+2*k3+k4)/6;
end
for p=1:N-1
t12=t(p)+h2;
m1=h*Yt;
p1=h*feval(ff,t(p),Y(p));
m2=h*(Yt+p1/2);
p2=h*feval(ff,t12,Y(p)+0.5*p1) ;
m3=(Yt+p2/2)*h;
p3=h*feval(ff,t12,Y(p)+0.5*p2) ;
m4=(Yt+p3)*h;
p4=h*feval(ff,t(p+1),Y(p)+p3);
Yt(p+1)=Yt(p)+(m1+2*m2+2*m3+m4 )/6;
Y(p+1)=Y(p)+(p1+2*p2+2*p3+p4)/6;
end
toc
plot(X(k),Y(p))
%résolution à l'aide de ODE45
tic
%[t,X]=ode45('f',[to tM], X(1) Xt(1));
%[t,Y]=ode45('ff',[to tM], Y(1) Yt(1));
toc
function Xde=f(t,Y);
global K X Y;
Xde=-K*X./sqrt(X.^2+Y.^2)^(3/2);
function Yde=ff(t,Y);
global K Y X;
Yde=-K*Y./sqrt(X.^2+Y.^2)^(3/2);
-----