正向差分(forward differences)的插值求解(python,数值积分)

第三十六篇 正向差分的插值求解

差分法

另一种求插值多项式的方法是通过np个数据点(xi, yi)， i = 0,1,2，…，n(其中n = np−1)，首先将插值多项式写成另外一种形式：

其中常数Ci, i = 0,1,2，…，n可以按照要求来确定

重新计算之后得到

计算实例
本文重复使用上篇中给出的例子。使用差分法去推导通过这些点的多项式

然后考虑额外的点去修正多项式

计算x=4.5时的y值
这个方程的第一部分牵涉到3个数据点，因此插值多项式将是二次型的。首先计算方程的常量

然后带入方程得到

这和上篇中的结果相同。
第四个数据点x3 = 5，y3 = 9代入将得到一个三次插值多项式。额外的常量将得到

然后将这个三次项加入之前的二次项中

或

三次项展示在下图然后得到x=4.5和y=8.175。

等间隔的差分法

如果数据以等距的x值提供，使xi−xi−1 = h，则对系数Ci的求导得到了很大程度的简化。
进行求解之前，介绍一个‘从前差分’的概念

和一个‘从后差分’的概念

而且，有一个‘二次差分’，也就是‘差分中的差分’，第二次从前差分能写成

一般形式为

在这次程序中，主要使用’从前差分’的方法
对于之前的方程，系数为

它能展现成一个通用形式

从之前可以得出Δjy0项可以通过递归来求值。为了便于手算，上表已经给出了计算步骤，其中有六个数据点
上表描述的特殊格式称为(牛顿)正向差分，其特征是沿从左上角到右下角向下倾斜的对角线下标保持相同。其他的布局是可能的，如逆向差分和高斯方法，但这些将不在这里讨论。当把表中的数据代入方程，“零”差项，也就是刚开始的差值项，是y值本身，因此

计算实例
给与下面的数据是函数y = cosx 的数据，求cons27
首先计算正向差分的数据，写入表格

一般方程的插值多项式,x0可以选择任意的初始值x。然而,如果我们希望包括所有五个数据点，生成一个四阶插值多项式， x0应设置等于表中的最高值,也就是说,x0 = 20。x值之间的常量间隔为h = 5，系数可以用以表格形式计算如下

保留小数点后五位，方程的插值多项式可以写成

扩展得到

替换需要求解的值，保留到小数点后五位得到

从上面看出，从上面的Q4(x)看出C3和C4项几乎对多项式没有作用。
如果系数变得足够小，截断插值多项式的能力是差分方法的一个有用特征。它不仅节省了计算的工作量，而且还提供了对数据点理论上起源的物理上的观测视角。
在拉格朗日方法中，这样的保存是不可能的，无论是否需要，完整的n阶多项式必须导出。拉格朗日方法确实有简化的优点，但是，是在只提供几个点，相对容易手算的情况下。
程序如下

#正向差分的插值求解
import numpy as np
npt=4;xi=2.8
c=np.zeros((npt))
diff=np.zeros((npt,npt))
x=np.zeros((npt))
x[0:npt]=(0.0,1.0,2.0,3.0)
diff[0:npt,0]=(1.0,1.0,15.0,61.0)
print('正向差分的插值求解')
print('数据点','  x    y')
for i in range(1,npt+1):
    print('{:13.4e}'.format(x[i-1]),end='')
    print('{:13.4e}'.format(diff[i-1,0]))
h=x[1]-x[0]
for i in range(1,npt):
    for j in range(1,npt-i+1):
        diff[j-1,i]=diff[j,i-1]-diff[j-1,i-1]
c[0]=diff[0,0];yi=c[0];term=1.0;factorial=1.0
for i in range(1,npt):
    factorial=factorial*i;c[i]=diff[0,i]/(factorial*h**i)
    term=term*(xi-x[i-1]);yi=yi+term*c[i]
print('多项式系数C')
for i in range(1,npt+1):
    print('{:13.4e}'.format(c[i-1]))
print('插值点','   x   y')
print('{:13.4e}'.format(xi),end='')
print('{:13.4e}'.format(yi))

终端输出结果如下