% numerical evidence of the order of accuracy
clear all;
k = 2; F = [1 2 3];
nPoints = 100;
err_vec = zeros(1,nPoints);
tau_vec = pi./(2*(1:nPoints));

for i = 1:nPoints

    x1 = psi_Extrap(0,2,tau_vec(i),k,F);

    xTemp = psi_Extrap(0,2,tau_vec(i)/2,k,F);
    x2 = psi_Extrap(tau_vec(i)/2,xTemp,tau_vec(i)/2,k,F);

    err_vec(i) = abs(x1-x2);

end

plot(log(tau_vec),log(err_vec),'*');axis equal;

f = fit(log(tau_vec)', log(err_vec)', 'poly1' );
display(f);

f = 

     Linear model Poly1:
     f(x) = p1*x + p2
     Coefficients (with 95% confidence bounds):
       p1 =       6.833  (6.738, 6.929)
       p2 =     -0.3618  (-0.678, -0.04554)

Published with MATLAB® R2014a