
function demo_integrazione

%--------------------------------------------------------------------------
% OGGETTO:
%--------------------------------------------------------------------------
% Test sulle formule dei trapezi composte e di Cavalieri-Simpson composte.
%--------------------------------------------------------------------------


% ..................... PARAMETRI DELL'ESPERIMENTO ........................ 

% DEFINIZIONE INTEGRANDA 
% Nota: attenzione ai "punti" per operazioni moltiplicative tra vettori,
% ma non tra "numeri e vettori".
% Alcuni esempi: 
% 1. f=@(x) sin(x)+5*exp(-x.^2)./(1+(cos(x)).^2)-3+log(x.^2+1);
% 2. f=@(x) sin(x+exp(x)); 
% 3. f=@(x) exp(-x.^2);
% 4. f=@(x) abs(x.^4+10*x.^3+19*x.^2-6*x-6).*exp(x);
% 5. f=@(x) exp(-x).*erf(sqrt(1250)*x+1.5);
% 6. f=@(x) exp(-x).*airy(-x); 

f=@(x) sin(x+exp(x)); % Nota: funzione ostica, molto oscillante in [0,4].

% ESTREMI DI INTEGRAZIONE.
a=0; b=4;

% METODO UTILIZZARE
% metodo=1: Trapezi composta, metodo=2: Cavalieri-Simpson composta.
metodo=1;



% .................... ESECUZIONE DELL'ESPERIMENTO ........................  

% INTEGRALE DI RIFERIMENTO DI "f" IN [a,b].
I=integral(f,a,b,'AbsTol',10^(-15),'RelTol',10^(-15));

% APPLICAZIONE DI TRAPEZI COMPOSTA PER VARI "N".
% a) Per ogni "N=2,4,8,...,512,1024", viene registrato ne vettore "errore"
% l'errore assoluto compiuto tra integrale esatto e quello numerico.
% b) Viene stampato a video il numero "N" e l'errore corrispondente, con
%    rispettivamente 2 cifre decimali, 1 cifra prima della virgola e 1 dopo
%    virgola, in notazione esponenziale.
errore=[];

NV=[2 4 8 16 32 64 128 256 512 1024]';

for N=[2 4 8 16 32 64 128 256 512 1024]

    % Vettore dei nodi "x" e vettore dei pesi "w".
    if metodo == 1
        [x,w]=trapezi_composta(N,a,b);
    else
        [x,w]=simpson_composta(N,a,b);
    end

    % Valutazione funzione nei nodi.
    fx=feval(f,x);

    % Calcolo degli integrali.
    I_num=w'*fx;

    % Valutazione dell'errore.
    errore(end+1,1)=abs(I-I_num);

    % Stampa dell'errore.
    fprintf('\n \t N: %2.0f errore: %1.1e',N,errore(end));
end

fprintf('\n');
if metodo == 1
        fprintf('\n \t * Metodo: formula dei trapezi composta.');
    else
        fprintf('\n \t * Metodo: formula di Cavalieri-Simpson composta.');
    end
fprintf('\n \n');

% GRAFICO DELLA INTEGRANDA
t=linspace(a,b,1000);
ft=feval(f,t);
clf; % cancella grafici precedenti
plot(t,ft,'k-','LineWidth',4);
 




function [x,w]=trapezi_composta(N,a,b)

%--------------------------------------------------------------------------
% OGGETTO:
%--------------------------------------------------------------------------
% Formula dei trapezi composta.
%--------------------------------------------------------------------------
% INPUT:
%--------------------------------------------------------------------------
% N: NUMERO SUBINTERVALLI.
% a, b: ESTREMI DI INTEGRAZIONE.
%--------------------------------------------------------------------------
% OUTPUT:
%--------------------------------------------------------------------------
% x: INTEGRAZIONE.
% w: PESI INTEGRAZIONE (INCLUDE IL PASSO!).
%--------------------------------------------------------------------------

h=(b-a)/N;             % PASSO INTEGRAZIONE.
x=a:h:b; x=x';         % NODI INTEGRAZIONE.

w=ones(N+1,1);         % PESI INTEGRAZIONE.
w(1)=0.5; w(N+1)=0.5;
w=w*h;




function [x,w]=simpson_composta(N,a,b)
%--------------------------------------------------------------------------
% OGGETTO:
%--------------------------------------------------------------------------
% Formula di Cavalieri-Simpson composta.
%--------------------------------------------------------------------------
% INPUT:
%--------------------------------------------------------------------------
% N: NUMERO SUBINTERVALLI.
% a, b: ESTREMI DI INTEGRAZIONE.
%--------------------------------------------------------------------------
% OUTPUT:
%--------------------------------------------------------------------------
% x: INTEGRAZIONE.
% w: PESI INTEGRAZIONE (INCLUDE IL PASSO!).
%--------------------------------------------------------------------------

h=(b-a)/N;             % AMPIEZZA INTERVALLO.
x=a:(h/2):b; x=x';     % NODI INTEGRAZIONE.

w=ones(2*N+1,1);           % PESI INTEGRAZIONE.
w(3:2:2*N-1,1)=2*ones(length(3:2:2*N-1),1);
w(2:2:2*N,1)=4*ones(length(2:2:2*N),1);
w=w*h/6;

