拉格朗日插值多项式

import numpy as np
import sympy # 符号运算库

class LagrangeInterpolation:
    '''
    拉格朗日插值
    '''
    def __init__(self,x,y):
        '''
        拉格朗日必要参数的初始化，及各健壮性的检测；健壮性：系统在不正常输入情况下仍能表现正常的程度
        :param x: 已知数据x的坐标点
        :param y: 已知数据y的坐标点
        '''
        #构造一下作为类的属性，防止后续进行数组的运算,做一个类型转换
        self.x = np.asarray(x, dtype=np.float)
        self.y = np.asarray(y, dtype=np.float) #类型转换，数据结构采用array

        if len(self.x)>1 and len(self.x)==len(self.y): #如果x，y是一个数据点，就不用做了，且个数应相同
            self.n=len(self.x) #有多少个数据点，即已知离散数据点的个数
        else:
            raise ValueError("插值数据（x，y）维度不匹配")
        #插值最终结果是一个多项式的形式，可作为一个类的属性由用户进行返回
        self.polynomial = None #最终的插值多项式，符号表示
        self.poly_coefficient = None # 最终插值多项式的系数向量，幂次从高到低
        self.coefficient_order = None # 对应多项式系数的阶次（所有项中最高次的幂）
        self.y0 = None # 所求插值点的值，单个值或向量

    def fit_interp(self): # 拟合
        '''
        核心算法：生成拉格朗日插值多项式
        :return:
        '''
        t = sympy.Symbol("t") # 定义符号变量，因为参数已经用过x
        self.polynomial = 0.0 # 插值多项式实例化
        for i in range(self.n): # 由多少个样本点，n个数据生成n-1次
            # 根据每个数据带点，构造插值基函数
            basis_fun = self.y[i]  # 插值基函数
            # 计算每个 yi*li（x）
            for j in range(i):
                basis_fun *= (t-self.x[j]) / (self.x[i]-self.x[j])
                
            for j in range(i + 1 , self.n):
                basis_fun *= (t - self.x[j]) / (self.x[i] - self.x[j])
            
            # 下面开始累加
            self.polynomial +=basis_fun # 插值多项式累加
        self.polynomial = sympy.expand(self.polynomial) #多项式的展开
        print(self.polynomial)


    def cal_interp_x0(self, x0): # 计算插值的值
        '''
        计算所给定的插值点的值，即插值
        :param x0: 所求插值的x坐标
        :return:
        '''
        pass

测试

from lagelangri.lagrange_interpolation import  LagrangeInterpolation
import numpy as np

if __name__ =='__main__':
    x = np.linspace(0, 2*np.pi, 10 ,endpoint=True) #0到2*pi，等分十个，最后一个值包含
    y=np.sin(x)
    """print(x)
    print(y)"""
    lag_interp = LagrangeInterpolation(x=x,y=y)
    lag_interp.fit_interp()

未简化的结果

0.920725428958529*t*(1.125 - 0.179049310978382*t)*(1.14285714285714 - 0.204627783975294*t)*(1.16666666666667 - 0.238732414637843*t)*(1.2 - 0.286478897565412*t)*(1.25 - 0.358098621956765*t)*(1.33333333333333 - 0.477464829275686*t)*(1.5 - 0.716197243913529*t)*(2.0 - 1.43239448782706*t) + 0.705316598492019*t*(1.28571428571429 - 0.204627783975294*t)*(1.33333333333333 - 0.238732414637843*t)*(1.4 - 0.286478897565412*t)*(1.5 - 0.358098621956765*t)*(1.66666666666667 - 0.477464829275686*t)*(2.0 - 0.716197243913529*t)*(3.0 - 1.43239448782706*t)*(1.43239448782706*t - 1.0) + 0.413496671566344*t*(1.5 - 0.238732414637843*t)*(1.6 - 0.286478897565412*t)*(1.75 - 0.358098621956765*t)*(2.0 - 0.477464829275686*t)*(2.5 - 0.716197243913529*t)*(4.0 - 1.43239448782706*t)*(0.716197243913529*t - 0.5)*(1.43239448782706*t - 2.0) + 0.122476942006377*t*(1.8 - 0.286478897565412*t)*(2.0 - 0.358098621956765*t)*(2.33333333333333 - 0.477464829275686*t)*(3.0 - 0.716197243913529*t)*(5.0 - 1.43239448782706*t)*(0.477464829275686*t - 0.333333333333333)*(0.716197243913529*t - 1.0)*(1.43239448782706*t - 3.0) - 0.0979815536051016*t*(2.25 - 0.358098621956765*t)*(2.66666666666667 - 0.477464829275686*t)*(3.5 - 0.716197243913529*t)*(6.0 - 1.43239448782706*t)*(0.358098621956765*t - 0.25)*(0.477464829275686*t - 0.666666666666667)*(0.716197243913529*t - 1.5)*(1.43239448782706*t - 4.0) - 0.206748335783172*t*(3.0 - 0.477464829275686*t)*(4.0 - 0.716197243913529*t)*(7.0 - 1.43239448782706*t)*(0.286478897565412*t - 0.2)*(0.358098621956765*t - 0.5)*(0.477464829275686*t - 1.0)*(0.716197243913529*t - 2.0)*(1.43239448782706*t - 5.0) - 0.201519028140577*t*(4.5 - 0.716197243913529*t)*(8.0 - 1.43239448782706*t)*(0.238732414637843*t - 0.166666666666667)*(0.286478897565412*t - 0.4)*(0.358098621956765*t - 0.75)*(0.477464829275686*t - 1.33333333333333)*(0.716197243913529*t - 2.5)*(1.43239448782706*t - 6.0) - 0.115090678619816*t*(9.0 - 1.43239448782706*t)*(0.204627783975294*t - 0.142857142857143)*(0.238732414637843*t - 0.333333333333333)*(0.286478897565412*t - 0.6)*(0.358098621956765*t - 1.0)*(0.477464829275686*t - 1.66666666666667)*(0.716197243913529*t - 3.0)*(1.43239448782706*t - 7.0) - 3.89817183251938e-17*t*(0.179049310978382*t - 0.125)*(0.204627783975294*t - 0.285714285714286)*(0.238732414637843*t - 0.5)*(0.286478897565412*t - 0.8)*(0.358098621956765*t - 1.25)*(0.477464829275686*t - 2.0)*(0.716197243913529*t - 3.5)*(1.43239448782706*t - 8.0)

self.polynomial = sympy.expand(self.polynomial)多项式展开后

-2.29369601281319e-6*t**9 + 6.48527268907878e-5*t**8 - 0.00061946678451133*t**7 + 0.00167479859899278*t**6 + 0.0040407632285091*t**5 + 0.00707566335693954*t**4 - 0.173843738744177*t**3 + 0.00401585715826069*t**2 + 0.999072579744674*t