曲线拟合 5 点

Question

我正在尝试在 C 中曲线拟合 5 个点。我使用上一篇文章 (Can sombody simplify this equation for me?) 中的这段代码来做 4 个点，但现在我需要添加另一个点。

// Input data: arrays x[] and y[]
// x[1],x[2],x[3],x[4] - X values
// y[1],y[2],y[3],y[4] - Y values

// Calculations
A = 0
B = 0
C = 0
D = 0
S1 = x[1] + x[2] + x[3] + x[4]
S2 = x[1]*x[2] + x[1]*x[3] + x[1]*x[4] + x[2]*x[3] + x[2]*x[4] + x[3]*x[4]
S3 = x[1]*x[2]*x[3] + x[1]*x[2]*x[4] + x[1]*x[3]*x[4] + x[2]*x[3]*x[4]
for i = 1 to 4 loop
   C0 = y[i]/(((4*x[i]-3*S1)*x[i]+2*S2)*x[i]-S3)
   C1 = C0*(S1 - x[i])
   C2 = S2*C0 - C1*x[i]
   C3 = S3*C0 - C2*x[i]
   A = A + C0
   B = B - C1
   C = C + C2
   D = D - C3
end-loop

// Result: A, B, C, D

我一直试图将其转换为 5 点曲线拟合，但无法弄清楚循环内部的内容：

// Input data: arrays x[] and y[]
// x[1],x[2],x[3],x[4],x[5] - X values
// y[1],y[2],y[3],y[4],y[5] - Y values

// Calculations
A = 0
B = 0
C = 0
D = 0
E = 0
S1 = x[1] + x[2] + x[3] + x[4]
S2 = x[1]*x[2] + x[1]*x[3] + x[1]*x[4] + x[2]*x[3] + x[2]*x[4] + x[3]*x[4]
S3 = x[1]*x[2]*x[3] + x[1]*x[2]*x[4] + x[1]*x[3]*x[4] + x[2]*x[3]*x[4]
S4 = x[1]*x[2]*x[3]*x[4] + x[1]*x[2]*x[3]*[5] + x[1]*x[2]*x[4]*[5] + x[1]*x[3]*x[4]*[5] + x[2]*x[3]*x[4]*[5]

for i = 1 to 4 loop
   C0 = ??
   C1 = ??
   C2 = ??
   C3 = ??
   C4 = ??
   A = A + C0
   B = B - C1
   C = C + C2
   D = D - C3
   E = E + C4
end-loop

// Result: A, B, C, D, E

任何填写 C0...C4 的帮助将不胜感激。我知道这与矩阵有关，但我无法弄清楚。带有伪代码或真实代码的示例将是最有帮助的。

谢谢

Answer 1

我拒绝错过这个概括的机会。 :)

相反，我们将稍微了解一下Lagrange polynomials 和Newton Divided Difference Method 的计算。

拉格朗日多项式

给定n+1 数据点，插值多项式为

l_j(i) 在哪里

.

这意味着我们可以通过对这些多项式求和来找到逼近n+1 点的多项式，而不管间距等。但是，这有点痛苦，我不想在 C 中这样做。让我们看看Newton Polynomials。

牛顿多项式

同样的开始，给定n+1 数据点，近似多项式将是

其中每个n(x) 是

，系数为

，为分差。

最终的形式看起来像

.

如您所见，考虑到除差值，公式非常简单。到目前为止，您只需执行每个新的划分差异并乘以每个点。应该注意的是，您最终会得到来自n+1 点的次数为n 的多项式。

分差

剩下的就是定义这两张图片最能解释的划分差异：

和

.

有了这些信息，C 实现应该是合理的。我希望这会有所帮助，我希望你学到了一些东西！ :)

Answer 2

如果x 值与x2-x1=h、x3-x2=h、x4-x3=h 和x5-x4=h 等距，则

C0 = y1;
C1 = -(25*y1-48*y2+36*y3-16*y4+3*y5)/(12*h);
C2 = (35*y1-104*y2+114*y3-56*y4+11*y5)/(24*h*h);
C3 = -(5*y1-18*y2+24*y3-14*y4+3*y5)/(12*h*h*h);
C4 = (y1-4*y2+6*y3-4*y4+y5)/(24*h*h*h*h);

y(x) = C0+C1*(x-x1)+C2*(x-x1)^2+C3*(x-x1)^3+C4*(x-x1)^4
// where `^` denotes exponentiation (and not XOR).