使用 Matlab 的拉格朗日多项式

Question

这是我必须解决的问题：

编写一个程序来评估和绘制x 在-5 和5 之间的u(x) = 1/(1+x^2) 的拉格朗日插值Iu(x)。对5,7,9,11,13,15 点插值（5、7、9 等数据点在 -5 和 5 之间，包括 -5 和 5）执行此操作。您的结果应该同时显示函数和插值。

这是我目前想出的代码：

int_pts = [5,7,9,11,13,15]; %various values for no. of point interpolants
n = length(int_pts);

u = @(x) 1./(1+x.^2); %the function we are interested in

numer = 1; 
denom = 1; %set initial values of numerator and denominator 

for k = 1 : n

    L = zeros( 1,int_pts(k) ); %create an array to store Lagrange polynomials 

    x = [-5 : 10/(int_pts(k)-1) : 5]; %create a list of x values

    udata = u(x); %create a list of corresponding values of u(x)

    for i = 1 : int_pts(k)

        if (i ~= k)

            denom = x(k) - x(i);

            syms z;
            numer = symfun( z - x(i) , z);        

            break

        end

        numer = numer * numer;
        denom = denom * denom;

        L(i) = numer / denom;

    end

end

我认为我走在正确的轨道上，但我真的开始纠结于如何提取数据并很好地绘制所有内容。我到处搜索，但大多数代码似乎只想输出某些数字的拉格朗日多项式的值，即 P(2)=...

Answer 1

在大多数情况下，您所开发的内容已引导您走上解决问题的正确道路。

您遇到的一些小问题是：

1) 您缺少for 循环。

在拉格朗日插值中，您应该有 n 项的总和，其中 n 是您拥有的数据点的数量。 并且这些术语由乘积乘以数据点的 y 分量组成，即 yj*∏{(x - xi)/(xj - xi)}，其中如果 i=j 则跳过该术语（这样做是为了确保分母永远不会为零）

2) 你的内部 for 循环和 if 语句不能很好地配合使用。

一旦您的脚本到达 if 语句的末尾，for 循环就会中断，并继续到下一个 k 值。

这应该替换为else，后跟continue。

3) 当您生成 numer 和 denom 变量时，您无意中只是将新值平方并将其存储到变量中。

你有numer = numer * numer这应该是numer = numer * numer_temp

这是解决方案的快速伪代码：

for n equals 'number of data points'

    Generate the n data points, P(x,y)

    for i equals 1 to n

        for j equals 1 to n

            if i is not equal to j
                Calculate Product of Numerators <- Function of x
                Calculate Product of Denominators
            else
                Continue to next iteration of loop
            endIf

        endFor

        Calculate Lagrange Polynomial < Function of x
    endFor

    Plot Lagrange
    Plot Original Function

endFor

由此您可以更新代码以正确确定拉格朗日多项式。

u = @(x) 1./(1+x.^2); % The function we are interested in
number_of_data_Points = [5,7,9,11,13,15];


syms x; % Initializes x to be symbolic

for n = number_of_data_Points

    xn = linspace(-5,5,n);      % Creates x values
    yn = u(xn);                 % Creates y values

    L=0; % Re/Initializes Langrange Polynomial
    for j = 1:n
        numer = 1; % Re/Initializes Numerator
        denom = 1; % Re/Initializes Denominator
        for i = 1:n
            if j~=i
                denom_temp = xn(j) - xn(i);

                numer_temp = x - xn(i);    

                numer =  numer * numer_temp;
                denom = denom * denom_temp;
            else
                continue
            end
        end
        L = L + yn(j) * numer / denom;
    end
    L = matlabFunction(L);
end

我将更改代码以输出您需要的结果给您。