%-------------------------------------------------------------------------------
% QUESTO CODICE SEGUE [ATKINSON, AN INTRODUCTION TO NUMERICAL ANALYSIS, p.414].
% NECESSITA DELLA FUNCTION: "g".
%-------------------------------------------------------------------------------

theta=1; mvect=[2 4 8 16]; sigma=10; tfin=5;       % SETTINGS.

% ESEMPIO.
G=inline('(-0.1+pi^2)*(exp((-0.1)*t).*sin(pi*x))',...
    'x','t');
d0=inline('zeros(size(t))','t');
d1=inline('exp((-0.1)*t).*sin(pi)','t');
f=inline('sin(pi*x)','x');
solution=inline('exp((-0.1)*t).*sin(pi*x)','x','t');

err_hist_prev_m=[];
fprintf('\n \t [THETA]:%3.3f [TFIN]:%3.3f [SIGMA]:%1.1e',theta,tfin,sigma);
for mindex=1:length(mvect)
    err_hist=[]; m=mvect(mindex); % m IN UNA LISTA.
    [V_hist,x_mid,t_mid,hx,ht]=cnheat(theta,tfin,m,f,d0,d1,G,sigma);
    [X, T]=meshgrid(x_mid,t_mid);  % VALUTAZIONE ERRORI.
    U=feval(solution,X,T); % SOL. ESATTA NELLA GRIGLIA.
    % ERRORE COMPIUTO PER "m" FISSATO SU TUTTA LA GRIGLIA (X,T).
    err=norm(U(size(U,1),:)-V_hist(size(V_hist,1),:),inf);
    err_hist=[err_hist; err];
    fprintf('\n \t [m]: %3.0f [ERROR]: %2.2e [hx]: %2.2e [ht]: %2.2e',m,err,hx,ht);
    if length(err_hist_prev_m) > 0
        fprintf(' [RATIO]: %2.2f',err_hist_prev_m(size(err_hist_prev_m,1))/err);
    end
    err_hist_prev_m=err_hist;
end
fprintf('\n \n');