Aidé moi en Matlab ??

[invite479ee65f]

Salut, notre professeur nous propose un projet sous Matlab MAIS LE problème que j'ai affronté c'est que la fonction ne s'affiche pas.
Svp est ce que vous pouvez le courigé si c'etais possible:

% boucle i de balayage des valeurs de dθ/dt(t=0)
 for i=-2:2
 % boucle m de balayage des valeurs de θ(τ0)
 for m=1:17;
 y=zeros(100,6);
 final=0;
 y(1,5)=m*5/180*pi; % angle θ initial (de 5° à 85 ° par pas de 5°)
 options = odeset('RelTol',1e-4);
 [t,y] = ode45('equations',[0 10],[0 25*cos(y(1,5)) 2 25*sin(y(1,5)) y(1,5) i],options);
 e=size(y);
 % boucle k enregistre à quel pas de temps k le javelot atteint le sol
 for k=1:e(1);
 if (y(k,3)>0) & (y(k+1,3)<0)
 final=k;
 break;
 end
 end
 % fin boucle k
 % Le vecteur à 3 dimensions range contiendra les résultats finaux
 range(1,m,i+3)=m*5; % angle θ initial en degrés
 range(2,m,i+3)=-y(final,3)*(y(final+1,1)-y(final,1))/(y(final+1,3)-y(final,3))+y(final,1);
 range(3,m,i+3)=i; % vitesse angulaire dθ/dt(t=0) initiale
 axis([0 150 0 50]);
 title(['Trajectoire pour une vitesse angulaire \theta/dt=',int2str(i),'rad/sec']);
 xlabel(['Portée [m]= ',int2str(range(2,m,i+3))]);
 ylabel('Hauteur [m]');
 comet(y(1:final+1,1),y(1:final+1,3));
 plot(y(1:final+1,1),y(1:final+1,3));
 hold on % force chacune des trajectoires à s'imprimer sur la même figure
 end
 % fin boucle m
 hold off
 end
 % fin boucle i
 for i=-2:2
 plot(range(1,:,i+3),range(2,:,i+3),'LineWidth',2);
 hold on;
 end
 title('Diagramme Portée en fonction de Angle Attaque \theta ');
 xlabel('Angle \theta [°]');
 ylabel('Portée[m]');
 %Sauvegarde dans le fichier javlot.txt de la portée du javelot en fonction de θ(t=0) pour
 chacune des 5 valeurs initiales de dθ/dt(t=0)
 fid = fopen('c:\javelot.txt','w');
 fprintf(fid,'angle deg; portée ; dteta/dt(t=0)\n');
 fprintf(fid,'% 6.2f ; %6.2f ; %6.2f \n',range);
 fclose(fid);
 hold off;
 Le fichier equation.m dont le module de résolution ode45.m a besoin est donné par:
 function z = equation(t,y)
 z = zeros(6,1);
 Delta=0.255;
 Speed = sqrtm(y(2)*y(2) + y(4)*y(4));
 Phi = atan(y(4)/y(2));
 Alpha = sin(y(5) - Phi);
 Drag = 0.00024*exp(5.157*abs(Alpha));
 Lift = 0.0;
 if Alpha > 0.001
 Lift = 0.0127*exp(1.34*log(abs(Alpha)));
 end
 z(1) = y(2);
 z(2) = - (Drag*y(2)*Speed + Lift*y(4)*Speed)/0.80625;
 z(3) = y(4);
 z(4) = - (Drag*y(4)*Speed - Lift*y(2)*Speed)/0.80625 - 9.81;
 z(5) = y(6);
 z(6) = - Delta*Speed*Speed*(Drag*sin(Alpha) + Lift*cos(Alpha))/0.42;

J'att vos répense
-----