clear all
close all
odefun = @(t,y) y*(1-y);
y0 = 1/2;
mrif = 50000;
tspan = linspace (0, 1, mrif + 1);
[trif, yrif] = euleroesplicito (odefun, tspan, y0);
plot (trif, yrif, '*')
legend ('numerical solution')
mrange = (10:2:20);
counter = 0;
options.Theta = 1/2;
options.Jacobian = @(t,y) (1-y)-y;
for m = mrange
  counter = counter + 1;
  tspan = linspace (0, 1, m + 1);
  [tout, yout] = thetametodo (odefun, tspan, y0, options);
  err(counter) = abs (yout(m+1) - yrif(mrif+1));
end
figure
loglog (mrange, err, '*', ...
	mrange, err(end) * (mrange/mrange(end)) .^ (-2))
  
