ctrapz.m

来自「数值分析各种公式的解法」· M 代码 · 共 19 行

function ctrapz(a,b,eps)
m=1
t0=0;
t2=(f(a)+f(b))/4+f((a+b)/2)
while abs(t2-t0)>eps
    m=m+1
    p=t2;
    t1=0;
    n=2^m;
    h=(b-a)/n;
    for k=0:(n-1)
        t1=t1+h*((f(a+k*h)+f(a+(k+1)*h))/4+f(a+(k+1/2)*h)/2);
    end
    t2=t1
    t0=p;
end
I=t2
function f=f(x)
f=exp(sin(x));