METODO DE RUUGE-KUTTA

clc
disp('METODO DE RUUGE-KUTTA')
syms x
syms y
f=inline(input('ingrese la derivada:','s'));
x=input('ingrese el valor de x:');
y=input('ingrese el valor de y:');
h=input('ingrese el valor de h:');
n=input('ingrese numero de iteraciones:');
clc     
     disp('x(n)     k1      k2         k3          k4  ');
for i=1:n
   k1=h*(feval(f,x,y));
   z=x+(1/2*(h));
   w=y+(1/2*(k1));
   k2=h*(feval(f,z,w));
   w=y+(1/2*(k2));
   k3=h*(feval(f,z,w));
   z=x+h;
   w=y+k3;
   k4=h*(feval(f,z,w));
  fprintf('n%0.1f  %0.4f  %0.4f  %0.4f      %0.4f     %0.4f',x,y,k1,k2,k3,k4);
   x=x+h;
  y=y+((1/6)*(k1+2*k2+2*k3+k4));
end