% ex1rk,  y'=sin(ty), y(0)=6  med Runge-Kuttas metod 
%  function f=fsin(t,y)
%     f=sin(t*y);
% --------------------------------------------------
  clear, clf
  tslut=2; 
  disp('     h        yslut ')                 
  t=0; y=6;  T=t; Y=y; h=0.2;         % stega med h=0.2 
  while t<tslut-h/2
    f1=fsin(t,y);
    f2=fsin(t+h/2,y+h*f1/2);
    f3=fsin(t+h/2,y+h*f2/2);
    f4=fsin(t+h,  y+h*f3);
    y=y+h/6*(f1+2*f2+2*f3+f4); t=t+h;  
    T=[T; t]; Y=[Y; y]; 
  end
  disp([h y]), plot(T,Y,'r', T,Y,'o')
  axis([0 tslut  5.8  6.4]), hold on
  t=0; y=6;  T=t; Y=y; h=0.1;         % räkna om med halverat steg
  while t<tslut-h/2
    f1=fsin(t,y);
    f2=fsin(t+h/2,y+h*f1/2);
    f3=fsin(t+h/2,y+h*f2/2);
    f4=fsin(t+h,  y+h*f3);
    y=y+h/6*(f1+2*f2+2*f3+f4); t=t+h;  
    T=[T; t]; Y=[Y; y]; 
  end
  disp([h y]), plot(T,Y,'g--', T,Y,'x')
  text(0.4,6.36,'Runge-Kuttas metod på  dy/dt=sin(ty),  y(0)=6')



