h=[1/4,1/8,1/16,1/32];
f = @(t,y) -15*y;
ytrue = @(t) exp(-15*t);
y0=1;
ye=zeros(size(h,2),1/h(1,end)*3+1);
ye(:,1)=y0*ones(size(h,2),1);

for j=1:size(h,2)
    for i=1:1/h(1,j)*5
        ye(j,i+1)=ye(j,i)+h(1,j)*f(i*h(1,j),ye(j,i));
    end
end

t=0:h(1,1):3;
ytrue_ev=arrayfun(ytrue,t);
plot(t,ytrue_ev,'b');
hold on
plot(t,ye(1,1:size(t,2)),'r');
title('h=1/4')

figure
t=0:h(1,2):3;
ytrue_ev=arrayfun(ytrue,t);
plot(t,ytrue_ev,'b');
hold on
plot(t,ye(2,1:size(t,2)),'r');
title('h=1/8')

figure
t=0:h(1,3):3;
ytrue_ev=arrayfun(ytrue,t);
plot(t,ytrue_ev,'b');
hold on
plot(t,ye(3,1:size(t,2)),'r');
title('h=1/16')

figure
t=0:h(1,4):3;
ytrue_ev=arrayfun(ytrue,t);
plot(t,ytrue_ev,'b');
hold on
plot(t,ye(4,1:size(t,2)),'r');
title('h=1/32')

