function Tn=mil(f,n,a,b)
h=(b-a)/n;
Tn=0
for i=0:(n-1)
Tn=Tn+f(a+(i+1/2)*h)
end
Tn=h*Tn;
endfunction
//
deff('y=f(x)','y=x*cos(x)^3');
a=0;b=%pi;
valexacte=-14/9;
err=[];valN=[];
for i=1:5
    err=[err,mil(f,10^i,a,b)-valexacte];   
    valN=[valN,10^i];
end
plot2d(valN,abs(err),logflag = 'll')
plot2d([10^2,10^3],[10^(-2),10^(-4)],2,logflag = 'll') // exemple de pente égale à -2
