用matlab编程romberg算法

下面是用Matlab编程实现Romberg算法的示例代码：

function r = romberg(f, a, b, n)
% f: 被积函数
% a, b: 积分下限和上限
% n: 迭代次数
% r: 最终估计积分值

h = b - a;
R = zeros(n, n);
R(1, 1) = h * (f(a) + f(b)) / 2;

for j = 2:n
    h = h / 2;
    s = 0;
    for i = 1:2^(j-2)
        x = a + (2*i-1)*h;
        s = s + f(x);
    end
    R(j, 1) = R(j-1, 1) / 2 + h * s;
    for k = 2:j
        R(j, k) = (4^(k-1) * R(j, k-1) - R(j-1, k-1)) / (4^(k-1) - 1);
    end
end

r = R(n, n);

使用方法：

假设要求函数 $f(x) = x^2$ 在区间 $[0, 1]$ 上的积分值，迭代次数为 $n=4$。

f = @(x) x.^2;
a = 0;
b = 1;
n = 4;
r = romberg(f, a, b, n);
disp(r);

输出结果为：

r = 0.3333