GNU Octave拉格朗日插值偏差模最大值计算性能优化求助

问题背景 #

给定插值节点X，需在区间[X(1), X(end)]内以步长h=(X(end)-X(1))/(100*(size(X,2)-1))遍历所有点，计算插值多项式与原函数exp(x)的偏差模最大值及对应点，但现有代码执行超时。

核心性能瓶颈 #

原代码的最大问题是逐个点调用lagrange函数，每次调用都要执行两层循环计算基函数，且重复计算基函数的分母项，导致时间复杂度陡增（O(M*N²)，M为采样点数，N为插值节点数）。此外，arrayfun的调用开销也会进一步拖慢速度。

优化建议 #

• 向量化计算替代逐点循环：让插值函数直接支持向量输入，一次性完成所有采样点的插值计算，避免重复初始化循环。
• 预计算基函数分母：提前计算拉格朗日基函数的分母部分X(i)-X(j)（i≠j），避免每次计算都重复运算。
• 简化原函数调用：exp函数本身支持向量运算，无需用arrayfun逐个处理。

优化后的代码实现 #

主代码 #

N = size(X, 2) - 1;  % 插值区间数
h = (X(end) - X(1)) / (100 * N);
x_h = X(1):h:X(end);
y_h = exp(x_h);  % 直接用向量运算替代arrayfun
lagr = lagrange_vec(X, exp(X), x_h);  % 用向量化的插值函数
[delta, idx] = delta_max(y_h, lagr);
x_delta_max = x_h(idx);

向量化拉格朗日插值函数（lagrange_vec） #

function y0 = lagrange_vec(X, Y, x0)
    n = length(X);
    m = length(x0);
    % 预计算分母：denom(i,j) = X(i)-X(j)，i≠j
    denom = X' - X;
    denom(logical(eye(n))) = 1;  % 对角线设为1，避免除以0
    
    % 计算每个x0点对应的所有基函数值
    numer = x0' - X;  % m×n矩阵，每行对应一个x0点与所有X的差
    basis = prod(numer ./ denom, 2);  % 逐行乘积得到每个点的基函数值
    
    y0 = basis * Y;  % 加权求和得到所有点的插值结果
endfunction

简化版delta_max函数（向量化实现） #

function [delta, idx] = delta_max(Y1,Y2)
    diffs = abs(Y1 - Y2);
    [delta, idx] = max(diffs);
endfunction

优化效果说明 #

• 时间复杂度从O(MN²)降低到O(MN)，当采样点数M和节点数N较大时，性能提升非常明显。
• 避免了arrayfun的函数调用开销，完全利用Octave/Matlab的向量运算优势。

内容的提问来源于stack exchange，提问作者JustUsNow