function Q=Romberg(f,a,b,kmax)
%用romberg公式求f(x)在[a,b]上的積分;
%最多計算kmax步;
Q(1,1)=trap(f,a,b,1);
Q(2,1)=trap(f,a,b,2);
Q(2,2)=(4*Q(2,1)-Q(1,1))/3;
for k=2:kmax
    n=2^k
    Q(k+1,1)=trap(f,a,b,n);
    %梯形值序列
    for j=2:k+1;
        c=4^(j-1);
        Q(k-j+2,j)=(c*Q(k-j+3,j-1)-Q(k-j+2,j-1))/(c-1);
    end
end