%improved euler

% y0 = 1;
% Tmax = 0.1;
% h = 0.05;
% 
% yy = improvedeuler(y0, Tmax, h)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function yy = improvedeuler(y0, Tmax, h) 
%Tmax    %final time
%h: step size, make sure you choose step size such that T/h is an integer
%y0: initial value at t=0

%yy output vector of the numerical solution

N = Tmax/h

t = 0;   %initial time
y = y0;   %initial value at t=0

yy = zeros(N+1,1); yy(1) = y0; %vector, saved solution

for i=1:N
    k1 = f(t,y);
    k2 = f(t+h,y+h*k1);
    y = y+h/2*(k1+k2);
    t = t+h;
    
    yy(i+1) = y;
end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%
function f=f(t,y)
    f = 2*y-3*t;
end

