Metodo de Euler Mejorado

function R = EulerM(a, b, M, y_0)

h = (b - a) / M;

Y = zeros(1, M+1);
Y(1) = y_0;
T = a:h:b;

for i = 1:M
    Y(i + 1) = Y(i) + h * f(T(i), Y(i));
    Y(i + 1) = Y(i) + h / 2 * ( f(T(i), Y(i)) + f(T(i+1), Y(i+1)));
end

R = T' Y';