Python使用cubic方法的时间插值方法 python插值计算

第2章-函数插值方法及其Python实现

• 1、前言
• 2、编程环境介绍及环境搭建
• 3、插值问题的提出

• 科学背景
• 插值和拟合的区别

• 4、插值问题的数学知识

• 插值问题
• 多项式插值的存在唯一性定理

• 5、常见插值公式

• 5.1 Lagrange插值公式

• 5.1.1 插值基函数与插值函数
• 5.1.2 Lagrange插值的Python实现

• 5.2 Newton插值公式

• 5.2.1 均差及其性质
• 5.2.2 Newton插值公式
• Newton插值的Python实现

• 5.3 Hermite插值公式

1、前言

这个学期，上了一门叫做《应用数值分析(华南理工大学出版社)》(如下图所示)的数学课，为了巩固一下所学知识，同时也为了给其他同学一个参考，特撰写《数值分析》系列博客，主要包括：

【数值分析】第2章-函数插值方法及其Python实现
【数值分析】第3章-曲线拟合/函数逼近及其Python实现
【数值分析】第4章-数值微分/数值积分及其Python实现
【数值分析】第5章-线性代数方程组数值解法——直接法及其Python实现
【数值分析】第6章-线性代数方程组数值解法——迭代法及其Python实现
【数值分析】第7章-非线性代数方程组数值解法及其Python实现
【数值分析】第9章-常微分方程数值解法及其Python实现

本篇主要讲解函数插值方法然后使用Python实现所讲解的算法。

2、编程环境介绍及环境搭建

Python由于其拥有大量的开源包(例如numpy、scipy、matplotlib)，最近几年，颇受各领域的科研工作者追捧。本系列文章都是使用python编程实现的。

3、插值问题的提出

科学背景

在科学研究或者工程问题当中，我们通常会遇到下列情况：

• 已知函数太复杂，难以输入到计算机，或者利用计算机运行一次，需要进行大量计算；
• 已知的函数只是一张离散的表格，也就是一些$(x_i, y_i)$点(由实验或观测得到)；
• 已知的数据往往会存在不完整的情况，希望通过某种方法补全缺失的数据。

插值和拟合的区别

这两种情况下，我们通常希望寻求简单的函数来逼近原函数。那么，有的同学可能就会问：这和拟合有什么区别？如果不知道什么是拟合，可以跳过这一小节。

拟合长这个样子：

插值长这个样子：

更直白地说，拟合通常不要求拟合函数完全穿过所给的数据点，但是插值，则通常要求插值函数经过所给的数据点。

4、插值问题的数学知识

插值问题

定义：设$f$是定义在区间$[a,b]$上的实值函数，并已知在$[a,b]$上$n+1$个互异节点$x_i$及相应的函数值$f_i=f(x_i)(i=0,1,\cdots,n)$，要求建造一个简单的、便于计算的满足下列条件的函数$p$：
$p(x_i)=f_i(i=0,1,\cdots,n)$

从而以$p$作为$f$的近似函数，其中
$f$：被插函数；
$p$：插值函数；
$x_i$：插值节点；
$[a,b]$：插值区间。

多项式插值的存在唯一性定理

设已知$[a,b]$上的函数$f$在$n+1$个互异节点$x_i \in [a,b]$上的值$f_i=f(x_i)(i=0,1,\cdots,n)$，则存在唯一的次数$\le n$的多项式$p(x) \in P_n$满足：
$p_n(x_i)=f_i(i=0,1,\cdots,n)$

5、常见插值公式

5.1 Lagrange插值公式

5.1.1 插值基函数与插值函数

根据多项式插值的存在唯一性定理，已知$[a,b]$上的函数$f$在$n+1$个互异节点$x_i \in [a,b]$，则必存在唯一的次数$\le n$的多项式$L_n(x)$满足：
$L_n(x_i)=f_i(i=0,1,\cdots,n)$
定义插值基函数$l_i(x)$：
$l_i(x)=\frac{(x-x_0)\cdots(x-x_{i-1})\cdots(x-x_n)}{(x_i-x_0)\cdots(x_i-x_{i-1})\cdots(x_i-x_n)}=\prod^n_{j=0,j\not=i}{\frac{x-x_j}{x_i-x_j}}(i=0,1,\cdots,n)$
则插值多项式可以写为：
$L_n(x)=\sum_{i=0}^n{l_i(x) \cdot f_i}$
相对应的插值余项(用于误差估计)：
$R_n(X)=\frac{1}{(n+1)!}f^{(n+1)}(\xi)\omega_{n+1}(x)$
其中，$\omega_{n+1}(x)=\prod^n_{j=0}{x-x_j}$

5.1.2 Lagrange插值的Python实现

1、定义插值基函数

class DimException(Exception):
    def __init__(self, X,Y):
        err = 'X和Y的维度不相等,请检查数据输入是否正确!'
        Exception.__init__(self, err)
        self.X = X
        self.Y = Y


def lagrange_inter_base(X,Y,x):
    '''
    n次Lagrange插值的基函数数值计算
    X:已知的n+1个点的x_i
    Y:已知的n+1个点的y_i
    X和Y的维度必须一致，否则报错
    x:待通过插值函数L(x)计算具体函数值的x
    Return:
    l:插值基函数计算x后返回的结果，它是一个列表，维度和X、Y一致
    '''
    if len(X)!=len(Y):
        raise DimException(X,Y)
    
    n = len(Y)
    l = []
    for i in range(n):
        li = 1 #基函数初始值设置为1，方便进行循环求乘
        # 计算第i个插值基函数的值
        for j in range(n):
            if j!=i:
                li=li*(x-X[j])/(X[i]-X[j])
        l.append(li)
    return l

2、定义插值函数

def lagrange_inter(X,Y,x):
    '''
    Lagrange插值
    n次Lagrange插值的数值计算
    X:已知的n+1个点的x_i
    Y:已知的n+1个点的y_i
    x:一维变量
    '''
    if len(X)!=len(Y):
        raise DimException(X,Y)
    
    n = len(Y)
    result=0.0
    l = lagrange_inter_base(X,Y,x) # 计算基函数值
    for i in range(n):
        result = result + l[i]*Y[i]  #通过插值多项式计算L(x)的值
    return result

3、实例求解
课本50页 例2.2.1 设已知$f(x)=e^{-x}$的4个函数值如下表所示：

x	0	1	2	3
$f(x)=e^{-x}$	1	0.3679	0.1352	0.0498

试用Lagrange公式的线性插值求$e^{-1.4}$的近似值，用二次插值求$e^{-2.1}$。

# 用Lagrange公式的线性插值求$e^{-1.4}$的近似值
print("(1) 用Lagrange公式的线性插值求$e^{-1.4}$的近似值")
X = [1,2]
Y = [0.3679,0.1353]
x = 1.4
result = lagrange_inter(X,Y,x)
print("e^{-1.4}的Lagrange插值计算结果为：",result)

返回结果如下：
(1) 用Lagrange公式的线性插值求$e^{-1.4}$的近似值
e^{-1.4}的Lagrange插值计算结果为： 0.27486000000000005

# 用二次插值求$e^{-2.1}$
print("(2) 用Lagrange公式的线性插值求$e^{-1.4}$的近似值")
X = [1,2,3]
Y = [0.3679,0.1353,0.0498]
x = 2.1
result = lagrange_inter(X,Y,x)
print("e^{-2.1}的Lagrange插值计算结果为：",result)

返回结果如下：
(2) 用Lagrange公式的线性插值求$e^{-1.4}$的近似值
e^{-2.1}的Lagrange插值计算结果为： 0.1201305

课本50页 例2.2.2 试证：由$n+1$个互异节点$x_i \in [a,b]$构成的$n+1$个插值基函数$l_i(x)(i=0,1,\cdots,n)$具有性质：
$\sum_{i=0}^n{l_i(x)=1}, \forall x \in [a,b]$

# 这里采用数值方法简单验证
import numpy as np
X = [0,1,2,3]
Y = [1,0.3679,0.1353,0.0498]
x_list = np.linspace(0.1,2.9,100)
for x in x_list:
    li = lagrange_inter_base(X,Y,x)
    sum_li = sum(li)
    print(sum_li)

结果略！

5.2 Newton插值公式

5.2.1 均差及其性质

均差（也称为差商）是数值方法中的一个重要概念，它可以反映出列表函数的性质，并能对Lagrange插值公式给出新的表达形式，也就是Newton插值公式 。

• 零阶均差
  $f[x_i]=f(x_i)$
• 1阶均差
  $f[x_i,x_j]=\frac{f(x_i)-f(x_j)}{x_i-x_j}$
• 2阶均差
  $f[x_i,x_j,x_k]=\frac{f[x_i,x_j]-f[x_j,x_k]}{x_i-x_k}$
• 递推地，n阶均差
  $f[x_0,x_1,\cdots,x_k]=\frac{f[x_0,x_1,\cdots,x_{k-1}]-f[x_1,x_2,\cdots,x_k]}{x_1-x_k}$

均差具有如下特性：

• k阶均差可以表示成k+1个函数值的线性组合；
• 均差对其节点是对称的；

5.2.2 Newton插值公式

未完待续

Newton插值的Python实现

5.3 Hermite插值公式

未完待续