复合高斯求积matlab代码

以下是一个简单的 MATLAB 代码实现复合高斯求积：

function I = composite_gauss_quad(f, a, b, n)
% f: 被积函数
% a: 积分下限
% b: 积分上限
% n: 子区间数量
% I: 积分结果

% 预先计算节点和权重
[x, w] = lgwt(n, a, b);

% 复合高斯求积
I = 0;
for i = 1:n
    x1 = x(i);
    x2 = x(i+1);
    I = I + w(i) * f((x1 + x2)/2) * (x2 - x1)/2;
end

end

function [x, w] = lgwt(N, a, b)
% N: 节点数量
% a: 积分下限
% b: 积分上限
% x: 节点
% w: 权重

% 找到初始节点和权重
x = linspace(-1, 1, N)';
x_old = ones(N, 1) * inf;
while max(abs(x - x_old)) > eps
    x_old = x;
    Pk = legendre_poly(N-1, x);
    dPk = diff(Pk, 1, 2);
    x = x_old - Pk(:,end) ./ dPk(:,end);
end
w = 2 ./ ((1 - x.^2) .* (diff(Pk, 1, 2).^2));

% 转换节点和权重到[a,b]区间
x = (b-a)/2 * x + (b+a)/2;
w = (b-a)/2 * w;
end

function P = legendre_poly(n, x)
% n: 多项式次数
% x: 自变量
% P: 勒让德多项式

if n == 0
    P = ones(size(x));
elseif n ==