clear;

t0 = 0;
y0 = 1;
tstop = 2;

exact = exct8(tstop);

i = 0;
for n = 10:10:100
    i = i + 1;
    yeuler = euler8(t0, y0, tstop, n);
    yrk4 = rk4(t0, y0, tstop, n);
    yeulerstop(i) = yeuler(length(yeuler));
    yrk4stop(i) = yrk4(length(yrk4));
    h(i) = (tstop - t0) / n;
    erreuler(i) = abs(exact - yeulerstop(i));
    errrk4(i) = abs(exact - yrk4stop(i));
end

approxeulerslope = (log10(erreuler(i)) - log10(erreuler(1))) / (log10(h(i)) - log10(h(1)));
approxrk4slope = (log10(errrk4(i)) - log10(errrk4(1))) / (log10(h(i)) - log10(h(1)));
slopestreuler = num2str(approxeulerslope);
slopestrrk4 = num2str(approxrk4slope);

figure(1);
plot(log(h), log(erreuler), 'go-');
xlabel('log(h)');
ylabel('log |error|');
streuler = strcat('Euler approx slope =  ', slopestreuler);
title(streuler);

figure(2);
plot(log(h), log(errrk4), 'mo-');
xlabel('log(h)');
ylabel('log |error|');
strrk4 = strcat('RK4 approx slope =  ', slopestrrk4);
title(strrk4);