function [x, y] = euler(f, a, b, y0, n)
% Euler Verfahren mit Heun Korrektur

h = (b-a)/n;
x = a:h:b;
y = zeros(1, n+1);

y(1) = y0 ;
for i = 2:n+1
    y(i) = y(i-1) + h * f(x(i-1), y(i-1));
    y(i) = y(i-1) + h * (f(x(i-1), y(i-1)) + f(x(i), y(i))) / 2;
end

end