改进Euler格式 | matlab

function [x, y]=mod_euler(f, a, b, y0, h)
% 改进的欧拉格式
% f是带求函数的一阶导形式
% a，b分别是积分上下限
% y0 是初始条件y（0）
% h是步长

s = (b - a) / h; % 求步数
X = zeros(1, s+1);
Y = zeros(1, s+1);
X = a:h:b;
Y(1) = y0;
for k = 1:s
    k1 = f(X(k), Y(k));
    k2 = f(X(k+1), Y(k) + h * k1);
    Y(k+1) = Y(k) + (h/2) * (k1 + k2)
end
x = X';
y = Y';

[x, y] = mod_euler(@f_test, 0, 0.5, 1, 0.1);
figure(1);
plot(x, y);
title('改进的欧拉格式')

%=> y = 1.0000    1.1105    1.2483    1.4248    1.6587    1.9833