% set the max. number of refinements
nrefmax = 4;

% set the order of the Lagrange elements used
p = 2;

% preallocate vectors:
Elem   = zeros(nrefmax+1,1);
Npts   = zeros(nrefmax+1,1);
DoF    = zeros(nrefmax+1,1);
normsq = zeros(nrefmax+1,1);
err    = zeros(nrefmax+1,1);
Rr     = zeros(nrefmax+1,1);
Qr     = zeros(nrefmax+1,1);


% obtain square of the norm of the FEM error
% on the refinement level r:
for r=0:nrefmax
    [e,nElem, nVertex, nDofs]=getmodel(r,p);
    normsq(r+1)=e;
    Elem(r+1) = nElem;
    Npts(r+1) = nVertex;
    DoF(r+1)  = nDofs;
end

for r=0:nrefmax
    err(r+1) = sqrt(normsq(r+1));
    if r>=1
        Rr(r+1)=err(r)/err(r+1);
        Qr(r+1)=log(Rr(r+1))/log(2);
    end
end
fprintf('Lagrange Elements with order p = %2d and nrefmax = %3d \n',p, nrefmax)
fprintf('    r     N_e    N_v    DoF  enorminfsq      enorminf        Rr       Qr\n')
for r = 0:nrefmax
    fprintf('%5d  %5d  %5d  %5d %11.3e %15.3e %9.2f %9.2f\n',r,Elem(r+1),Npts(r+1), DoF(r+1),normsq(r+1),err(r+1),Rr(r+1),Qr(r+1))
end