romberg.m

来自「非线性方程牛顿迭代法和数值积分的龙贝格解法」· M 代码 · 共 38 行

 Romberg.m
function I=Romberg(fun,a,b,e)
m=1;
h=b-a;
I=h/2*(feval(fun,a)+feval(fun,b));
T(1,1)=I;
while 1
    N=2^(m-1);
    h=h/2;
    I=I/2;
    
    
    for i=1:N
        I=I+h*feval(fun,a+(2*i-1)*h);
    end
    T(m+1,1)=I;
    M=2*N;k=1;
    
    
    
    
    while(M>1)
        T(m+1,k+1)=(4^k*T(m+1,k)-T(m,k))/(4^k-1);
        M=M/2;
        k=k+1;
    end




    if (abs(T(k,k)-T(k-1,k-1))<e)
        break;
    end
    m=m+1;
end
I=T(k,k);
T