
function I=romberg(a,b,n,ftype)

for i=1:(n+1)
    [x_tpz,w_tpz]=simpson_composta(2^(i-1),a,b);
    fx_tpz=f(x_tpz,ftype);                 
    A(i,1)=w_tpz'*fx_tpz;
end

for j=2:(n+1)
    for i=j:(n+1)
        A(i,j)=(4^(j-1)*A(i,j-1)-A(i-1,j-1))/(4^(j-1)-1);
    end
end

I=A(n+1,n+1);