Python|Python数值插值和曲线拟合

通常参与技术计算。 数据来源可以是实验观察或数值计算。 插值和曲线拟合之间是有区别的。 在插值中,我们通过数据点构造一条曲线。 这样做时,我们隐含了一个假设,即数据点是准确且不同的。 相反,曲线拟合应用于通常由测量误差引起的包含散射(噪声)的数据。 在这里,我们要找到一条平滑的曲线,可以从某种意义上近似该数据。 因此,曲线不一定会击中数据点。 插值和曲线拟合之间的差异如图1所示。
图1
多项式插值 拉格朗日法
插值的最简单形式是多项式。 始终可以构造一个唯一的次数为 n n n的多项式,该多项式通过 n + 1 n+1 n+1个不同的数据点。 获得该多项式的一种方法是拉格朗日公式,
P n ( x ) = ∑ i = 0 n y i ? i ( x ) , ( 1 a ) P_{n}(x)=\sum_{i=0}^{n} y_{i} \ell_{i}(x),\qquad(1a) Pn?(x)=i=0∑n?yi??i?(x),(1a)
其中下标 n n n表示多项式的次数,而
? i ( x ) = x ? x 0 x i ? x 0 ? x ? x 1 x i ? x 1 ? ? x ? x i ? 1 x i ? x i ? 1 ? x ? x i + 1 x i ? x i + 1 ? x ? x n x i ? x n = ∏ j = 0 j ≠ i n x ? x i x i ? x j , i = 0 , 1 , … , n , ( 1 b ) \begin{aligned}\ell_{i}(x) &=\frac{x-x_{0}}{x_{i}-x_{0}} \cdot \frac{x-x_{1}}{x_{i}-x_{1}} \cdots \cdots \frac{x-x_{i-1}}{x_{i}-x_{i-1}} \cdot \frac{x-x_{i+1}}{x_{i}-x_{i+1}} \cdots \frac{x-x_{n}}{x_{i}-x_{n}} \\&=\prod_{j=0 \atop j \neq i}^{n} \frac{x-x_{i}}{x_{i}-x_{j}}, \quad i=0,1, \ldots, n\end{aligned},\qquad(1b) ?i?(x)?=xi??x0?x?x0???xi??x1?x?x1????xi??xi?1?x?xi?1???xi??xi+1?x?xi+1???xi??xn?x?xn??=j?=ij=0?∏n?xi??xj?x?xi??,i=0,1,…,n?,(1b)
被称为基数函数。
例如,如果 n + 1 n+1 n+1,则插值是直线 P 1 ( x ) = y 0 ? 0 ( x ) + y 1 ? 1 ( x ) P_{1}(x)=y_{0} \ell_{0}(x)+y_{1} \ell_{1}(x) P1?(x)=y0??0?(x)+y1??1?(x),其中
? 0 ( x ) = x ? x 1 x 0 ? x 1 ? 1 ( x ) = x ? x 0 x 1 ? x 0 \ell_{0}(x)=\frac{x-x_{1}}{x_{0}-x_{1}} \quad \ell_{1}(x)=\frac{x-x_{0}}{x_{1}-x_{0}} ?0?(x)=x0??x1?x?x1???1?(x)=x1??x0?x?x0??
当 n = 2 n=2 n=2时,插值是抛物线形: P 2 ( x ) = y 0 ? 0 ( x ) + y 1 ? 1 ( x ) + y 2 ? 2 ( x ) P_{2}(x)=y_{0} \ell_{0}(x)+y_{1} \ell_{1}(x)+y_{2} \ell_{2}(x) P2?(x)=y0??0?(x)+y1??1?(x)+y2??2?(x),其中
? 0 ( x ) = ( x ? x 1 ) ( x ? x 2 ) ( x 0 ? x 1 ) ( x 0 ? x 2 ) ? 1 ( x ) = ( x ? x 0 ) ( x ? x 2 ) ( x 1 ? x 0 ) ( x 1 ? x 2 ) ? 2 ( x ) = ( x ? x 0 ) ( x ? x 1 ) ( x 2 ? x 0 ) ( x 2 ? x 1 ) \begin{array}{l}\ell_{0}(x)=\frac{\left(x-x_{1}\right)\left(x-x_{2}\right)}{\left(x_{0}-x_{1}\right)\left(x_{0}-x_{2}\right)} \\\ell_{1}(x)=\frac{\left(x-x_{0}\right)\left(x-x_{2}\right)}{\left(x_{1}-x_{0}\right)\left(x_{1}-x_{2}\right)} \\\ell_{2}(x)=\frac{\left(x-x_{0}\right)\left(x-x_{1}\right)}{\left(x_{2}-x_{0}\right)\left(x_{2}-x_{1}\right)}\end{array} ?0?(x)=(x0??x1?)(x0??x2?)(x?x1?)(x?x2?)??1?(x)=(x1??x0?)(x1??x2?)(x?x0?)(x?x2?)??2?(x)=(x2??x0?)(x2??x1?)(x?x0?)(x?x1?)??
基数函数是次数为n的多项式,并具有以下性质
? i ( x j ) = { 0ifi ≠ j 1ifi = j } = δ i j , ( 2 ) \ell_{i}\left(x_{j}\right)=\left\{\begin{array}{l}0 \text { if } i \neq j \\1 \text { if } i=j\end{array}\right\}=\delta_{i j},\qquad(2) ?i?(xj?)={0 if i?=j1 if i=j?}=δij?,(2)
其中 δ i j \delta_{i j} δij?是Kronecker增量。 图2中针对 x 0 = 0 , x 1 = 2 x_{0}=0, x_{1}=2 x0?=0,x1?=2和 x 2 = 3 x_{2}=3 x2?=3的三点插值( n = 2 n=2 n=2)演示了此属性。
图2
为了证明插值多项式通过数据点,我们将 x = x j x=x_{j} x=xj?代入等式(1a),然后使用等式(2)。 结果是
P n ( x j ) = ∑ i = 0 n y i ? i ( x j ) = ∑ i = 0 n y i δ i j = y j P_{n}\left(x_{j}\right)=\sum_{i=0}^{n} y_{i} \ell_{i}\left(x_{j}\right)=\sum_{i=0}^{n} y_{i} \delta_{i j}=y_{j} Pn?(xj?)=i=0∑n?yi??i?(xj?)=i=0∑n?yi?δij?=yj?
可以证明多项式插值的误差为
f ( x ) ? P n ( x ) = ( x ? x 0 ) ( x ? x 1 ) ? ( x ? x n ) ( n + 1 ) ! f ( n + 1 ) ( ξ ) , ( 3 ) f(x)-P_{n}(x)=\frac{\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n}\right)}{(n+1) !} f^{(n+1)}(\xi),\qquad(3) f(x)?Pn?(x)=(n+1)!(x?x0?)(x?x1?)?(x?xn?)?f(n+1)(ξ),(3)
ξ \xi ξ位于区间 ( x 0 , x n ) \left(x_{0}, x_{n}\right) (x0?,xn?)的某处;其值是未知的。 请注意,数据点离 x x x越远,则它对 x x x处的误差的贡献就越大。
牛顿法
尽管拉格朗日法从概念上讲很简单,但它并不适合有效的算法。 使用牛顿法可以获得更好的计算过程,其中插值多项式以以下形式编写:
P n ( x ) = a 0 + ( x ? x 0 ) a 1 + ( x ? x 0 ) ( x ? x 1 ) a 2 + ? + ( x ? x 0 ) ( x ? x 1 ) ? ( x ? x n ? 1 ) a n P_{n}(x)=a_{0}+\left(x-x_{0}\right) a_{1}+\left(x-x_{0}\right)\left(x-x_{1}\right) a_{2}+\cdots+\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n-1}\right) a_{n} Pn?(x)=a0?+(x?x0?)a1?+(x?x0?)(x?x1?)a2?+?+(x?x0?)(x?x1?)?(x?xn?1?)an?
该多项式有助于进行有效的评估程序。 例如,考虑四个数据点( n = 3 n=3 n=3)。 这里的插值多项式是
P 3 ( x ) = a 0 + ( x ? x 0 ) a 1 + ( x ? x 0 ) ( x ? x 1 ) a 2 + ( x ? x 0 ) ( x ? x 1 ) ( x ? x 2 ) a 3 = a 0 + ( x ? x 0 ) { a 1 + ( x ? x 1 ) [ a 2 + ( x ? x 2 ) a 3 ] } \begin{aligned}P_{3}(x) &=a_{0}+\left(x-x_{0}\right) a_{1}+\left(x-x_{0}\right)\left(x-x_{1}\right) a_{2}+\left(x-x_{0}\right)\left(x-x_{1}\right)\left(x-x_{2}\right) a_{3} \\&=a_{0}+\left(x-x_{0}\right)\left\{a_{1}+\left(x-x_{1}\right)\left[a_{2}+\left(x-x_{2}\right) a_{3}\right]\right\}\end{aligned} P3?(x)?=a0?+(x?x0?)a1?+(x?x0?)(x?x1?)a2?+(x?x0?)(x?x1?)(x?x2?)a3?=a0?+(x?x0?){a1?+(x?x1?)[a2?+(x?x2?)a3?]}?
可以通过以下递归关系向后求值:
P 0 ( x ) = a 3 P 1 ( x ) = a 2 + ( x ? x 2 ) P 0 ( x ) P 2 ( x ) = a 1 + ( x ? x 1 ) P 1 ( x ) P 3 ( x ) = a 0 + ( x ? x 0 ) P 2 ( x ) \begin{array}{l}P_{0}(x)=a_{3} \\P_{1}(x)=a_{2}+\left(x-x_{2}\right) P_{0}(x) \\P_{2}(x)=a_{1}+\left(x-x_{1}\right) P_{1}(x) \\P_{3}(x)=a_{0}+\left(x-x_{0}\right) P_{2}(x)\end{array} P0?(x)=a3?P1?(x)=a2?+(x?x2?)P0?(x)P2?(x)=a1?+(x?x1?)P1?(x)P3?(x)=a0?+(x?x0?)P2?(x)?
对于任意 n n n我们有
P 0 ( x ) = a n P k ( x ) = a n ? k + ( x ? x n ? k ) P k ? 1 ( x ) , k = 1 , 2 , … , n , ( 4 ) P_{0}(x)=a_{n} \quad P_{k}(x)=a_{n-k}+\left(x-x_{n-k}\right) P_{k-1}(x), \quad k=1,2, \ldots, n,\qquad(4) P0?(x)=an?Pk?(x)=an?k?+(x?xn?k?)Pk?1?(x),k=1,2,…,n,(4)
用xData表示数据点的 x x x坐标数组,用 n n n表示多项式的次数,我们有以下算法来计算 P n ( x ) P_{n}(x) Pn?(x):

p = a[n] for k in range(1,n+1): p = a[n-k] + (x - xData[n-k])*p

通过强制多项式通过每个数据点来确定 P n P_{n} Pn?的系数: y i = P n ( x i ) , i = 0 , 1 , … , n y_{i}=P_{n}\left(x_{i}\right), \quad i=0,1, \ldots, n yi?=Pn?(xi?),i=0,1,…,n。 这产生了这些联立方程:
y 0 = a 0 y 1 = a 0 + ( x 1 ? x 0 ) a 1 y 2 = a 0 + ( x 2 ? x 0 ) a 1 + ( x 2 ? x 0 ) ( x 2 ? x 1 ) a 2 ? y n = a 0 + ( x n ? x 0 ) a 1 + ? + ( x n ? x 0 ) ( x n ? x 1 ) ? ( x n ? x n ? 1 ) a n , ( a ) \begin{aligned}&\begin{array}{l}y_{0}=a_{0} \\y_{1}=a_{0}+\left(x_{1}-x_{0}\right) a_{1} \\y_{2}=a_{0}+\left(x_{2}-x_{0}\right) a_{1}+\left(x_{2}-x_{0}\right)\left(x_{2}-x_{1}\right) a_{2}\end{array}\\& \vdots\\& y_{n}=a_{0}+\left(x_{n}-x_{0}\right) a_{1}+\cdots+\left(x_{n}-x_{0}\right)\left(x_{n}-x_{1}\right) \cdots\left(x_{n}-x_{n-1}\right) a_{n}\end{aligned},\qquad(a) ?y0?=a0?y1?=a0?+(x1??x0?)a1?y2?=a0?+(x2??x0?)a1?+(x2??x0?)(x2??x1?)a2???yn?=a0?+(xn??x0?)a1?+?+(xn??x0?)(xn??x1?)?(xn??xn?1?)an??,(a)
引入除数差
? y i = y i ? y 0 x i ? x 0 , i = 1 , 2 , … , n ? 2 y i = ? y i ? ? y 1 x i ? x 1 , i = 2 , 3 , … , n ? 3 y i = ? 2 y i ? ? 2 y 2 x i ? x 2 , i = 3 , 4 , … n ? ? n y n = ? n ? 1 y n ? ? n ? 1 y n ? 1 x n ? x n ? 1 , ( 5 ) \begin{aligned}&\begin{aligned}\nabla y_{i} &=\frac{y_{i}-y_{0}}{x_{i}-x_{0}}, \quad i=1,2, \ldots, n \\\nabla^{2} y_{i} &=\frac{\nabla y_{i}-\nabla y_{1}}{x_{i}-x_{1}}, \quad i=2,3, \ldots, n \\\nabla^{3} y_{i} &=\frac{\nabla^{2} y_{i}-\nabla^{2} y_{2}}{x_{i}-x_{2}}, \quad i=3,4, \ldots n\end{aligned}\\&\qquad\qquad\vdots\\&\nabla^{n} y_{n}=\frac{\nabla^{n-1} y_{n}-\nabla^{n-1} y_{n-1}}{x_{n}-x_{n-1}}\end{aligned},\qquad(5) ??yi??2yi??3yi??=xi??x0?yi??y0??,i=1,2,…,n=xi??x1??yi???y1??,i=2,3,…,n=xi??x2??2yi???2y2??,i=3,4,…n???nyn?=xn??xn?1??n?1yn???n?1yn?1???,(5)
等式(a)的解为
a 0 = y 0 a 1 = ? y 1 a 2 = ? 2 y 2 ? a n = ? n y n , ( 6 ) a_{0}=y_{0} \quad a_{1}=\nabla y_{1} \quad a_{2}=\nabla^{2} y_{2} \quad \cdots \quad a_{n}=\nabla^{n} y_{n},\qquad(6) a0?=y0?a1?=?y1?a2?=?2y2??an?=?nyn?,(6)
如果系数是手工计算的,则使用下表1中的格式( n = 4 n=4 n=4所示)很方便。
表1
表中的对角项( y 0 , ? y 1 , ? 2 y 2 , ? 3 y 3 y_{0}, \nabla y_{1}, \nabla^{2} y_{2}, \nabla^{3} y_{3} y0?,?y1?,?2y2?,?3y3?和 ? 4 y 4 \nabla^{4} y_{4} ?4y4?)是多项式的系数。 如果以不同的顺序列出数据点,则表中的条目将发生变化,但所得多项式将是相同的-请记住, n n n次插值 n + 1 n+1 n+1个不同数据点的多项式是唯一的。
可以使用以下算法在一维数组 a a a中进行机器计算(我们使用符号 m = n + 1 = m=n+1= m=n+1=数据点数):
a = yData.copy() for k in range(1,m): for i in range(k,m): a[i] = (a[i] - a[k-1])/(xData[i] - xData[k-1])

最初,a包含数据的y坐标,因此它与上表1中的第二列相同。 每次通过外部循环都会在下一列中生成条目,这些条目将覆盖a的相应元素。 因此,最终包含表1的对角项(即多项式的系数)。
牛顿插值算法代码
newtonPoly模块包含牛顿法进行插值所需的两个函数。 给定数据点数组xData和yData,函数coeffts返回系数数组a。 找到系数后,可以使用函数evalPoly在x的任何值上评估内插 P n ( x ) P_{n}(x) Pn?(x)。

内维尔法
牛顿的插值方法包括两个步骤:系数的计算,然后是多项式的求值。 如果使用相同的多项式在 x x x的不同值处重复执行插值,则效果很好。 如果仅要插值一点,则使用内维尔算法一步计算插值的方法是更方便的选择。
令 P k [ x i , x i + 1 , … , x i + k ] P_{k}\left[x_{i}, x_{i+1}, \ldots, x_{i+k}\right] Pk?[xi?,xi+1?,…,xi+k?]表示 k k k次多项式通过 k + 1 k+1 k+1个数据点 ( x i , y i ) , ( x i + 1 , y i + 1 ) , … , ( x i + k , y i + k ) \left(x_{i}, y_{i}\right),\left(x_{i+1}, y_{i+1}\right), \ldots,\left(x_{i+k}, y_{i+k}\right) (xi?,yi?),(xi+1?,yi+1?),…,(xi+k?,yi+k?)。 对于单个数据点,我们有
P 0 [ x i ] = y i , ( 7 ) P_{0}\left[x_{i}\right]=y_{i},\qquad(7) P0?[xi?]=yi?,(7)
基于两个数据点的插值是
P 1 [ x i , x i + 1 ] = ( x ? x i + 1 ) P 0 [ x i ] + ( x i ? x ) P 0 [ x i + 1 ] x i ? x i + 1 P_{1}\left[x_{i}, x_{i+1}\right]=\frac{\left(x-x_{i+1}\right) P_{0}\left[x_{i}\right]+\left(x_{i}-x\right) P_{0}\left[x_{i+1}\right]}{x_{i}-x_{i+1}} P1?[xi?,xi+1?]=xi??xi+1?(x?xi+1?)P0?[xi?]+(xi??x)P0?[xi+1?]?
很容易验证 P 1 [ x i , x i + 1 ] P_{1}\left[x_{i}, x_{i+1}\right] P1?[xi?,xi+1?]通过两个数据点; 也就是说,当 x = x i x=x_{i} x=xi?时, P 1 [ x i , x i + 1 ] = y i P_{1}\left[x_{i}, x_{i+1}\right]=y_{i} P1?[xi?,xi+1?]=yi?,当 x = x i + 1 x=x_{i+1} x=xi+1?时, P 1 [ x i , x i + 1 ] = y i + 1 P_{1}\left[x_{i}, x_{i+1}\right]=y_{i+1} P1?[xi?,xi+1?]=yi+1?。
三点插值是
P 2 [ x i , x i + 1 , x i + 2 ] = ( x ? x i + 2 ) P 1 [ x i , x i + 1 ] + ( x i ? x ) P 1 [ x i + 1 , x i + 2 ] x i ? x i + 2 P_{2}\left[x_{i}, x_{i+1}, x_{i+2}\right]=\frac{\left(x-x_{i+2}\right) P_{1}\left[x_{i}, x_{i+1}\right]+\left(x_{i}-x\right) P_{1}\left[x_{i+1}, x_{i+2}\right]}{x_{i}-x_{i+2}} P2?[xi?,xi+1?,xi+2?]=xi??xi+2?(x?xi+2?)P1?[xi?,xi+1?]+(xi??x)P1?[xi+1?,xi+2?]?
为了证明该插值确实与数据点相交,我们首先将 x = x i x=x_{i} x=xi?代入,得到
P 2 [ x i , x i + 1 , x i + 2 ] = P 1 [ x i , x i + 1 ] = y i P_{2}\left[x_{i}, x_{i+1}, x_{i+2}\right]=P_{1}\left[x_{i}, x_{i+1}\right]=y_{i} P2?[xi?,xi+1?,xi+2?]=P1?[xi?,xi+1?]=yi?
同样, x = x i + 2 x=x_{i+2} x=xi+2?得出
P 2 [ x i , x i + 1 , x i + 2 ] = P 2 [ x i + 1 , x i + 2 ] = y i + 2 P_{2}\left[x_{i}, x_{i+1}, x_{i+2}\right]=P_{2}\left[x_{i+1}, x_{i+2}\right]=y_{i+2} P2?[xi?,xi+1?,xi+2?]=P2?[xi+1?,xi+2?]=yi+2?
最后,当 x = x i + 1 x=x_{i+1} x=xi+1?时,我们有
P 1 [ x i , x i + 1 ] = P 1 [ x i + 1 , x i + 2 ] = y i + 1 P_{1}\left[x_{i}, x_{i+1}\right]=P_{1}\left[x_{i+1}, x_{i+2}\right]=y_{i+1} P1?[xi?,xi+1?]=P1?[xi+1?,xi+2?]=yi+1?
于是
P 2 [ x i , x i + 1 , x i + 2 ] = ( x i + 1 ? x i + 2 ) y i + 1 + ( x i ? x i + 1 ) y i + 1 x i ? x i + 2 = y i + 1 P_{2}\left[x_{i}, x_{i+1}, x_{i+2}\right]=\frac{\left(x_{i+1}-x_{i+2}\right) y_{i+1}+\left(x_{i}-x_{i+1}\right) y_{i+1}}{x_{i}-x_{i+2}}=y_{i+1} P2?[xi?,xi+1?,xi+2?]=xi??xi+2?(xi+1??xi+2?)yi+1?+(xi??xi+1?)yi+1??=yi+1?
建立了模式之后,我们现在可以推导出通用的递归公式:
P k [ x i , x i + 1 , … , x i + k ] = ( x ? x i + k ) P k ? 1 [ x i , x i + 1 , … , x i + k ? 1 ] + ( x i ? x ) P k ? 1 [ x i + 1 , x i + 2 , … , x i + k ] x i ? x i + k , ( 8 ) \begin{aligned}& P_{k}\left[x_{i}, x_{i+1}, \ldots, x_{i+k}\right]=\\& \frac{\left(x-x_{i+k}\right) P_{k-1}\left[x_{i}, x_{i+1}, \ldots, x_{i+k-1}\right]+\left(x_{i}-x\right) P_{k-1}\left[x_{i+1}, x_{i+2}, \ldots, x_{i+k}\right]}{x_{i}-x_{i+k}}\end{aligned},\qquad(8) ?Pk?[xi?,xi+1?,…,xi+k?]=xi??xi+k?(x?xi+k?)Pk?1?[xi?,xi+1?,…,xi+k?1?]+(xi??x)Pk?1?[xi+1?,xi+2?,…,xi+k?]??,(8)
给定 x x x的值,可以按照以下表格格式(针对四个数据点显示)进行计算:
表2
用 m m m表示数据点的数量,计算表元素的算法为

该算法适用于一维数组y,该数组最初包含数据的y值(表2中的第二列)。 每次通过外部循环,都会计算下一列中的y元素,这些元素将覆盖先前的条目。 在过程的最后,y包含表格的对角线项。 穿过所有数据点的插值的值(在x处评估)是y的第一个元素
内维尔插值算法代码

多项式插值的局限性
多项式插值应使用最少可行的数据点数进行。 如果数据点之间的距离很近,则使用最接近的两个点进行线性插值通常就足够了。 在大多数情况下,三到六个最近的邻居点会产生良好的结果。 相交超过六个点的插值器必须令人怀疑。 原因是远离兴趣点的数据点不会影响插值的准确性。 实际上,它们可能是有害的。
使用过多点的危险如图3所示。 圆圈表示11个等距的数据点。 实线是内插值,它是所有点相交的次数为10的多项式。 如图所示,如此高的多项式有在数据点之间过度振荡的趋势。 通过使用跨越四个最近邻点的三次插值,可以获得更加平滑的结果。
图3
多项式外推(在数据点范围之外进行内插)是危险的。 作为示例,请考虑图4。 六个数据点显示为圆圈。 第五级插值多项式由实线表示。 插值在数据点的范围内看起来很好,但是当x> 12时,明显偏离明显的趋势。例如,在这种情况下,将x推算为14时,y是不合理的。
图4
如果无法避免外推,则可以使用以下措施:
  • 绘制数据并直观地验证外推值是否有意义。
  • 使用基于最近邻居数据点的低阶多项式。 例如,对于图4中的数据,线性或二次插值将产生 y ( 14 ) y(14) y(14)的合理估计。
  • 处理对数x与对数y的关系图,该图通常比x-y曲线平滑得多,因此推断起来更安全。 通常,该图几乎是一条直线。 这在图5中进行了说明,该图表示图4中数据的对数图。
图5
示例1
给定数据点
使用拉格朗日法确定 x = 1 x=1 x=1处的 y y y
? 0 = ( x ? x 1 ) ( x ? x 2 ) ( x 0 ? x 1 ) ( x 0 ? x 2 ) = ( 1 ? 2 ) ( 1 ? 3 ) ( 0 ? 2 ) ( 0 ? 3 ) = 1 3 ? 1 = ( x ? x 0 ) ( x ? x 2 ) ( x 1 ? x 0 ) ( x 1 ? x 2 ) = ( 1 ? 0 ) ( 1 ? 3 ) ( 2 ? 0 ) ( 2 ? 3 ) = 1 ? 2 = ( x ? x 0 ) ( x ? x 1 ) ( x 2 ? x 0 ) ( x 2 ? x 1 ) = ( 1 ? 0 ) ( 1 ? 2 ) ( 3 ? 0 ) ( 3 ? 2 ) = ? 1 3 y = y 0 ? 0 + y 1 ? 1 + y 2 ? 2 = 7 3 + 11 ? 28 3 = 4 \begin{array}{l}\ell_{0}=\frac{\left(x-x_{1}\right)\left(x-x_{2}\right)}{\left(x_{0}-x_{1}\right)\left(x_{0}-x_{2}\right)}=\frac{(1-2)(1-3)}{(0-2)(0-3)}=\frac{1}{3} \\\\\ell_{1}=\frac{\left(x-x_{0}\right)\left(x-x_{2}\right)}{\left(x_{1}-x_{0}\right)\left(x_{1}-x_{2}\right)}=\frac{(1-0)(1-3)}{(2-0)(2-3)}=1 \\\\\ell_{2}=\frac{\left(x-x_{0}\right)\left(x-x_{1}\right)}{\left(x_{2}-x_{0}\right)\left(x_{2}-x_{1}\right)}=\frac{(1-0)(1-2)}{(3-0)(3-2)}=-\frac{1}{3} \\\\y=y_{0} \ell_{0}+y_{1} \ell_{1}+y_{2} \ell_{2}=\frac{7}{3}+11-\frac{28}{3}=4\end{array} ?0?=(x0??x1?)(x0??x2?)(x?x1?)(x?x2?)?=(0?2)(0?3)(1?2)(1?3)?=31??1?=(x1??x0?)(x1??x2?)(x?x0?)(x?x2?)?=(2?0)(2?3)(1?0)(1?3)?=1?2?=(x2??x0?)(x2??x1?)(x?x0?)(x?x1?)?=(3?0)(3?2)(1?0)(1?2)?=?31?y=y0??0?+y1??1?+y2??2?=37?+11?328?=4?
示例2
数据点
在多项式上。 通过构造一个类似于表1的除差表来确定该多项式的次数。
这是用于得出表中数字的一些示例计算:
? y 2 = y 2 ? y 0 x 2 ? x 0 = 59 ? ( ? 1 ) 4 ? ( ? 2 ) = 10 ? 2 y 2 = ? y 2 ? ? y 1 x 2 ? x 1 = 10 ? 1 4 ? 1 = 3 ? 3 y 5 = ? 2 y 5 ? ? 2 y 2 x 5 ? x 2 = ? 5 ? 3 ? 4 ? 4 = 1 \begin{aligned}\nabla y_{2} &=\frac{y_{2}-y_{0}}{x_{2}-x_{0}}=\frac{59-(-1)}{4-(-2)}=10 \\\\\nabla^{2} y_{2} &=\frac{\nabla y_{2}-\nabla y_{1}}{x_{2}-x_{1}}=\frac{10-1}{4-1}=3 \\\\\nabla^{3} y_{5} &=\frac{\nabla^{2} y_{5}-\nabla^{2} y_{2}}{x_{5}-x_{2}}=\frac{-5-3}{-4-4}=1\end{aligned} ?y2??2y2??3y5??=x2??x0?y2??y0??=4?(?2)59?(?1)?=10=x2??x1??y2???y1??=4?110?1?=3=x5??x2??2y5???2y2??=?4?4?5?3?=1?
从表中我们可以看到,牛顿多项式的最后一个非零系数(最后一个非零对角线项)为 ? 3 y 3 \nabla^{3} y_{3} ?3y3?,这是三次项的系数。 因此,多项式是三次的。
示例3
给定数据点
用内维尔法确定 y ( x ) = 0 y(x)=0 y(x)=0的根。
这是逆插值的示例,其中 x x x和 y y y的角色互换。 我们不是在给定的 x x x上计算 y y y,而是找到与给定 y y y对应的 x x x(在这种情况下, y = 0 y=0 y=0)。 使用表2的格式(当然 x x x和 y y y互换了),我们得到
以下是表中使用的示例计算:
表格中的所有 P P P都是根的估计值,根是由涉及不同数据点的不同插值顺序得出的。 例如, P 1 [ y 0 , y 1 ] P_{1}\left[y_{0}, y_{1}\right] P1?[y0?,y1?]是根据前两个点通过线性插值获得的根,而 P 2 [ y 1 , y 2 , y 3 ] P_{2}\left[y_{1}, y_{2}, y_{3}\right] P2?[y1?,y2?,y3?]是使用后三个点进行二次插值的结果。 从所有四个数据点的三次插值获得的根是 x = P 3 [ y 0 , y 1 , y 2 , y 3 ] = 3.8317 x=P_{3}\left[y_{0}, y_{1}, y_{2}, y_{3}\right]=3.8317 x=P3?[y0?,y1?,y2?,y3?]=3.8317。
示例4
表中的数据点位于 f ( x ) = 4.8 cos ? π x 20 f(x)=4.8 \cos \frac{\pi x}{20} f(x)=4.8cos20πx?的图上。 通过牛顿法在 x = 0 , 0.5 , 1.0 , … , 8.0 x=0,0.5,1.0, \ldots, 8.0 x=0,0.5,1.0,…,8.0处插值此数据,然后将结果与“精确”值 y i = f ( x i ) y_{i}=f\left(x_{i}\right) yi?=f(xi?)进行比较。
【Python|Python数值插值和曲线拟合】详情参阅 - 亚图跨际

    推荐阅读