1.函数式
import numpy as np
'''
(X,Y):需要插值的一系列点,array数据类型
point:需要计算插值的点
'''
def _Lagrange(X,Y,point):
'''挑战一行代码实现Lagrange插值'''
return np.sum([(((np.prod(point-X)/(point-X[i]))/np.prod(np.delete((X[i]-X),i)))*Y[i]) for i in range(len(X)) if point != X[i]])
def Lagrange(X,Y,point):
'''正常写法'''
result = 0.0
for i in range(len(X)):
if point != X[i]:
Li = (np.prod(point-X)/(point-X[i])) / np.prod(np.delete((X[i]-X),i))
result = result + Li*Y[i]
else:
result = Y[i]
return result
def Difference(X,Y):
'''差分函数'''
return np.sum([Y[i]/np.prod(np.delete((X[i]-X),i)) for i in range(len(X))])
def Newton(X,Y,point):
'''牛顿插值'''
result = Y[0]
for i in range(1,len(X)):
result = result + Difference(X[:i+1],Y[:i+1])*np.prod(point-X[:i])
return result
X = np.arange(10)
Y = X**2
point = 3.5
print('真实值 ',point**2)
print('拉格朗日插值1:',_Lagrange(X,Y,point))
print('拉格朗日插值2:',Lagrange(X,Y,point))
print('牛顿插值',Newton(X,Y,point))
2.面向对象式
import numpy as np
class Interpolation(object):
def __init__(self,X,Y,point):
'''
@author:zengwei
(X,Y):给定的用于插值的点
point:需要用插值函数计算值的点
'''
self.X = np.array(X)
self.Y = Y
self.point = point
def _Lagrange(self):
'''挑战一行代码实现Lagrange插值'''
return np.sum([(((np.prod(self.point-self.X)/(self.point-self.X[i]))/np.prod(np.delete((self.X[i]-self.X),i)))*self.Y[i]) for i in range(len(self.X)) if self.point!=self.X[i]])
def Lagrange(self):
'''Lagrange插值的正常写法'''
result = 0.0
for i in range(len(self.X)):
if self.point != self.X[i]:
Li = (np.prod(self.point-self.X)/(self.point-self.X[i]))/np.prod(np.delete((self.X[i]-self.X),i))
result = result + Li*self.Y[i]
else:
result = self.Y[i]
return result
def Difference(X,Y):
return np.sum([Y[i]/np.prod(np.delete((X[i]-X),i)) for i in range(len(X))])
def Newton(self):
'''牛顿插值'''
result = self.Y[0]
for i in range(1,len(self.X)):
result = result + Interpolation.Difference(self.X[:i+1],self.Y[:i+1])*np.prod(self.point-self.X[:i])
return result
def Hermite(self):
'''Hermite插值'''
pass
def CubicSpline(self):
'''三次样条插值'''
pass
X = np.arange(10)
Y = X**2
point = 3.5
myop = Interpolation(X,Y,point)
print('真实值 ',point**2)
print('拉格朗日插值1:',myop._Lagrange())
print('拉格朗日插值2:',myop.Lagrange())
print('牛顿插值',myop.Newton())
print('计算差分',Interpolation.Difference(X,Y))