% Funzione integranda
f = @(x) x./(x.^2+1);
% Estremi di integrazione
a = 0; b = 5;

% Valore esatto
Int_exact = log(26)/2;

% Massima suddivisione dell'intervallo
Nmax = 20;
% Inizializzo i vettori errore
errPT = zeros(1,Nmax); errT = zeros(1,Nmax); errS = zeros(1,Nmax);
% Ciclo su N
for N = 1:Nmax
    % Genero i punti dei vari metodi
    [xPT,wPT] = ptomed_composta(a,b,N);
    [xT,wT] = trapezi_composta(a,b,N);
    [xS,wS] = simpson_composta(a,b,N);
    % Calcolo integrale approssimato con le formule
    IntPT = sum(wPT.*f(xPT));
    IntT = sum(wT.*f(xT));
    IntS = sum(wS.*f(xS));
    % Errori relativi utilizzando la suddivisione in N sottointervalli 
    % nella posizione N-esima del vettore 
    errPT(N) = abs(IntPT - Int_exact)/abs(Int_exact);
    errT(N) = abs(IntT - Int_exact)/abs(Int_exact);
    errS(N) = abs(IntS - Int_exact)/abs(Int_exact);
end
% Grafico sovrapposto in scala semilogaritmica
figure(1)
semilogy(1:Nmax, errPT,'bo-');
hold on;
semilogy(1:Nmax,errT,'mo-');
semilogy(1:Nmax,errS,'ko-');
legend('Pto Medio', 'Trapezi','Simpson');