% Euleresim1virhe.m
% Tutkitaan virhekäytöstä:

clear; close all       
n=2; h=1/n; y1=[];  

for laskuri=1:5
  n=2*n; h=1/n;
  t=linspace(0,1,n);
  y(1)=3;   % Indeksointi alkaa 1:stä Matlabissa.
  
  for k=1:n-1
      y(k+1)=y(k)+h*sin(t(k)*y(k));
  end;
  
  y1=[y1,y(n)]      % Kerätään y1-vektoriin Euler-approksimaatioita pisteessä t=1.
  TY{laskuri}=[t;y]   % Matlabin "uusi piirre". ty-arvojen taulukko.
  figure(laskuri)
  plot(t,y); hold on
  if n <=16, plot(t,y,'o'), end
  title(['h=',num2str(h)])
  hold off
  %ylim([2.8 10])
end;

% Lasketaan Matlabin valmiilla pitkälle viilatulla ode45-ratkaisijalla.
%  Menetelmänä Runge-Kutta_Fehlberg.
sinty=inline('sin(t*y)','t','y')
[T,Y]=ode45(sinty,[0,1],3);

figure

plot(T,Y)
title('ode45:lla laskettu')
ode45y1=Y(end)




% Pidetään ode45-laskua tarkkana arvona ja verrataan Euler-laskelmia siihen:

virheet=Y(end)-y1
suhteet=virheet(2:end)./virheet(1:end-1)  % Kunkin seuraavan virheen suhde edelliseen.

%virheet =
%    0.1071    0.0358    0.0149    0.0068    0.0033
%suhteet =
%    0.3341    0.4160    0.4573    0.4785

% Nähdään Eulerille tyypillinen käytös: Kun askel puolittuu, niin virhe puolittuu, eli
% virhekäytös on luokkaa O(h).
