%
% RK4 integration of IVP y' = y, y(0) = 1, for 0 <= t <= 1
%
f = @(t,y) y;

a = 0;
b = 1;

n = 2^5;
h = (b - a)/(n - 1);
yy = zeros(n, 1);
tt = zeros(n, 1);

t = a;
y = 1;

tt(1) = a;
yy(1) = y;
for k = 2:n
    y = rk4step(y, t, h, f);
    t = t + h;
    yy(k) = y;
    tt(k) = t;
end

hold on
plot(tt, exp(tt), 'ro')
plot(tt, yy, 'k')
hold off
