romberg.m

来自「matlab中运用Simpson求积时的求积函数」· M 代码 · 共 16 行

function r = romberg(f,a,b,n) 
h = b-a; 
r(1,1) = h*(feval(f,a)+feval(f,b))/2 
for i =2:n 
sum = 0; 
for k =1:2^(i-2) 
sum = sum+feval(f,a+(k-0.5)*h); 
end; 
r(2,1)=(r(1,1)+h*sum)/2; 
for j = 2:i 
r(2,j)=r(2,j-1)+(r(2,j-1)-r(1,j-1))/(4^(j-1)-1); 
end; 
r(2,:) 
h = h/2; 
r(1,:)=r(2,:); 
end;