分段⼆次插值例题_计算⽅法（2）——插值法（附Python程
序）
给定⼀些数据，⽣成函数的⽅式有两种：插值，回归。
插值⽽得到的函数通过数据点，回归得到的函数不⼀定通过数据点。
Hermite插值的程序，
拉格朗⽇插值，⽜顿插值和Hermite插值
下⾯给出拉格朗⽇插值，⽜顿插值
具体原理可参考课本，不再赘述。
拉格朗⽇插值法
线性插值 ⼀次精度 需要2个节点
⼆次插值 ⼆次精度 需要3个节点
n次插值 n次精度 需要n+1个节点
拉格朗⽇插值代码段（根据传⼊数据⾃动判断次数）：
# 返回多项式
def p(x,a):
"""
p(x,a)是x的函数，a是各幂次的系数
"""
s = 0
for i in range(len(a)):
s += a[i]*x**i
return s
# n次拉格朗⽇插值
def lagrange_interpolate(x_list,y_list,x):
"""
x_list 待插值的x元素列表
y_list 待插值的y元素列表
插值以后整个lagrange_interpolate是x的函数
"""
if len(x_list) != len(y_list):
raise ValueError("list x and list y is not of equal length!")
# 系数矩阵
A = []
for i in range(len(x_list)):
A.append([])
for j in range(len(x_list)):
A[i].append(pow(x_list[i],j))linspace numpy
b = []
for i in range(len(x_list)):
b.append([y_list[i]])
# 求得各阶次的系数
a = lu_solve(A, b) # ⽤LU分解法解线性⽅程组，可以使⽤numpy的类似函数
a = transpose(a)[0] # change col vec a into 1 dimension
val = p(x,a)
print(x,val)
return val
其中lu_solve(A,b)是⾃⼰写的轮⼦，可以⽤numpy的numpy.linalg.sovle(A,b)来代替
到后⾯会有⼀期讲矩阵⽅程直接法，会有讲到如何写lu_solve()
看⼀看插值的效果如何
import numpy as np
import matplotlib.pyplot as plt
# 待插值的元素值
x_points = [0,1,2,3,4,5]
y_points = [1,5,4,8,7,12]
# 拉格朗⽇插值
x = np.linspace(0,5)
y = list(map(lambda t: lagrange_interpolate(x_points,y_points,t),x)) # 画图
plt.scatter(x_points,y_points,color = "orange")
plt.plot(x,y)
plt.legend(["lagrange interpolation","scattered points"])
plt.show()
拉格朗⽇插值
可以看到，插值之后的函数可以较好地反映原数据点的情况。
⽜顿插值法
涉及到的两个表，差分表和差商表：
差分表
差商表
递归打印差分表
def difference_list(dlist): # Newton
if len(dlist)>0:
print(dlist)
prev,curr = 0,0
n = []
for i in dlist:
curr = i
n.append(curr - prev)
prev = i
n.pop(0)
difference_list(n)
打印差商表，并以列表的形式返回图中蓝⾊圈起来的部分
def difference_quotient_list(y_list,x_list = []):
if x_list == []:
x_list = [i for i in range(len(y_list))]
print(y_list)
prev_list = y_list
dq_list = []
dq_list.append(prev_list[0])
for t in range(1,len(y_list)):
prev,curr = 0,0
m = []
k = -1
for i in prev_list:
curr = i
m.append((curr - prev)/(x_list[k+t]-x_list[k]))
prev = i
k+=1
m.pop(0)
prev_list = m
dq_list.append(prev_list[0])
print(m)
return dq_list
⽜顿插值代码段（调⽤了之前求差商表的代码）
def newton_interpolate(x_list,y_list,x):
coef = difference_quotient_list(y_list,x_list)
p = coef[0]
for i in range(1,len(coef)):
product = 1
for j in range(i):
product *= (x - x_list[j] )
p += coef[i]*product
return p
if __name__ == '__main__':
import numpy as np
import matplotlib.pyplot as plt
# 待插值的元素值
x_points = [0,1,2,3,4,5]
y_points = [1,5,4,8,7,12]
# ⽜顿插值
x = np.linspace(0,5)
y = list(map(lambda t: newton_interpolate(x_points,y_points,t),x))
# 画图
plt.scatter(x_points,y_points,color = "orange")
plt.plot(x,y)
plt.legend(["newton interpolation","scattered points"])
plt.show()
插值效果图
⽜顿插值
可以看到，相同的插值节点，使⽤拉格朗⽇插值和⽜顿插值的结果是相同的。
Hermite 插值法
三次Hermite插值不但要求在插值节点上插值多项式与函数值相等，还要求插值多项式的导数与函数导数相等。进⾏Hermite插值需要六个信息：
这个⽐较好理解，通过
两点的插值有⽆限种⽅式。⽐如：
但是，通过限制住两个端点的导数值，就可以确定下来⼤致的插值形状。（对于Hermite插值，六个条件能确定出唯⼀的插值结果）例如，可以编程求出

版权声明：本站内容均来自互联网，仅供演示用，请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198，我们将在24小时内删除。