2

我正在尝试编写一些执行多项式插值的代码，但我无法让它工作。我是一名学生，我正在遵循此视频中的逻辑https://www.youtube.com/watch?v=VpI-wC94RKw以便在 Matlab 中以代码形式重新创建它，但每当我尝试创建我的视频中显示的矩阵的自己版本我得到了一个几乎完全用零填充的矩阵（一个元素除外）。我不明白为什么会这样。

我的代码：

x=[150, 200, 300, 500, 1000, 2000, 99999]';
y=[2, 3, 4, 5, 6, 7, 8]';

function f = interPoly(x,y)
    % Skapar en matris A var varje rad är [x_1^6, x_1^5,..., 1] till [x_n^6, x_n^5,..., 1]
    A = [x.^6 x.^5 x.^4 x.^3 x.^2 x ones(numel(x),1) y];
    % Gaussar matris A
    R = rref(A);
    % Plockar sista kolumnen ur R
    c = R(:,end);
    f = c(1)*x.^6+c(2)*x.^5+c(3)*x.^4+c(4)*x.^3+c(5)*x.^2+c(6)*x+c(7);
end

（矩阵'A'是这里有问题的矩阵。我最后得到的函数也只是用零作为值填充。也很抱歉评论是瑞典语）

我在 x 和 y 中有 7 个值，因此是 6 阶多项式，但我真的不知道倒数第二列中的常数应该是什么，所以我只是在那里放了一堆（我是新手，所以我我有点不确定逻辑）。

无论如何，我已经尝试将相同的功能与其他一些输入数据一起使用，并且效果很好。

替代输入数据：

x=[0, 0.5, 1, 1.5, 2, 2.99, 3]';
y=[0, 0.52, 1.09, 1.75, 2.45, 3.5, 4]';

它是否因为元素溢出而给我零（例如 99999^6 是一个非常高的数字）？我真的不明白这里发生了什么以及为什么它与一组不同的输入数据工作得很好。帮助？

谢谢！

编辑：这项任务的全部要点（由我的学校给出）是将“最小二乘法”方法（我也编写了代码但未发布）与多项式插值方法（上面代码中的方法）进行比较。上面'x'中的最后一个值应该是无穷大（f（inf）= 8），所以我只是用一个非常高的数字替换它，因此它不是“均匀”分布的。有一个更好的方法吗？

4

1 回答 1

4

在您的示例中，您得到以下矩阵R：

   1.00000   0.00000  -0.00000  -0.00000  -0.00000  -0.00000  -0.00000  -0.00000   1.6925e-05
   0.00000   1.00000   0.00051   0.00000   0.00000   0.00000   0.00000   0.00000   7.1286e-05
   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   5.4078e-04
   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   6.9406e-03
   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   2.2098e-01
   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   7.0000e+00
   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   8.0000e+00

如果您定义的方程组有一个唯一解，那么您应该有一个单位矩阵（右侧有这个附加列）。在你的情况下，你有（除了前两行）方程

0 == 5.4078e-04
0 == 6.9406e-03
0 == 2.2098e-01
0 == 7.0000e+00
0 == 8.0000e+00

这表明您的系统不承认解决方案。这是由于您选择的输入理论上可能有解决方案，但实际上数值属性非常糟糕。事实上，如果我们计算条件数cond(A)，我们会得到一个4.3691e+06对于这么小的矩阵来说非常糟糕的值。这很可能是由于您的x价值观分布不均。（您正在尝试解决Vandermonde 系统。）如果您以更“均匀”的分布方式，事情看起来会好很多。

除此之外，您的代码看起来不错，但您可能希望在用于定义它们的值之外的值中评估插值函数。此外，已知这种“精确”方法在数值上不是很理想，因此可能值得使用比分解更稳定的算法来求解rref系统qr。我建议进行以下更改：

function f = interPoly(x,y,xnew)
    A = [x.^6 x.^5 x.^4 x.^3 x.^2 x ones(numel(x),1)];
    [q,r] = qr(A);
    c = r\(q' * y);
    disp(cond(A));
    f = c(1)*xnew.^6+c(2)*xnew.^5+c(3)*xnew.^4+c(4)*xnew.^3+c(5)*xnew.^2+c(6)*xnew+c(7);
end

x=[150, 200, 300, 500, 1000, 2000, 99999]';
y=[2, 3, 4, 5, 6, 7, 8]';


disp(interPoly(x,y,x+eps));

在线尝试！

于 2019-01-23T15:52:52.317 回答