Analisis numerico

function[ti,xi,yi]=RK4(f1,f2,I,h,alfa) % I=[a,b], alfa[x0,y0] a=I(1);b=I(2);x0=alfa(1);y0=alfa(2)N=(b-a)/h; ti=a:h:b; xi(1)=x0; yi(1)=y0; %extremos a, b, entero N; condiciones iniciales for i=1:Nk(1)=h*f1(ti(i),xi(i),yi(i)); l(1)=h*f2(ti(i),xi(i),yi(i)); k(2)=h*f1(ti(i)+h/2,xi(i)+k(1)/2,yi(i)+l(1)/2);l(2)=h*f2(ti(i)+h/2,xi(i)+k(1)/2,yi(i)+l(1)/2); k(3)=h*f1(ti(i)+h/2,xi(i)+k(2)/2,yi(i)+l(2)/2);l(3)=h*f2(ti(i)+h/2,xi(i)+k(2)/2,yi(i)+l(2)/2); k(4)=h*f1(ti(i)+h/2,xi(i)+k(3),yi(i)+l(3));l(4)=h*f2(ti(i)+h/2,xi(i)+k(3),yi(i)+l(3)); %iteraciones de los k(1),l(1),...l(4)xi(i+1)=xi(i)+1/6*(k(1)+2*k(2)+2*k(3)+k(4)); yi(i+1)=yi(i)+1/6*(l(1)+2*l(2)+2*l(3)+l(4)); %Evaluacion de los xi(i+1) e yi(i+1) EndCorrida de la function con el ejemplo en clases

>> f1=@(t,x,y)-0.5*x+y;f2=@(t,x,y)4-0.1*x-0.3*y;I=[0,2];alfa=[4,6]; >> [ti,xi,yi]=RK4(f1,f2,I,0.5,alfa);[ti;xi;yi]'

y0 =

6

ans =

0 4.00006.0000 0.5000 5.9594 6.7888 1.0000 7.7873 7.3796 1.5000 9.4335 7.8074 2.0000 10.8735 8.1040

>>