1. 拉格朗日插值
  2. 均差与牛顿插值多项式
  3. 埃尔米特插值
  4. 分段低次插值
  5. 三次样条插值

有些老师如同大山一样压抑着学生。

华中师范大学真是太多好老师了，可惜这辈子，在本科期间是无缘再见了。

潘老师的学习建议

重基础，多练习，勤思考
I hear and I forget, I see and I remember, I do and I understand.

  • 注意掌握各种数值方法的思想和原理，而不是死记硬背
  • 注意数值方法的处理技巧，以及与计算机编程的结合（编程实践）
  • 注重必要的数值计算训练（必要的手算推导能力）
  • 要重视算法的误差分析、收敛性和稳定性

纸上得来终觉浅，绝知此事要躬行。

问题的提出

白话：根据有限数据来推测未知点的值

插值函数：

P(xi)=yi(i=0,1,⋯ ,n)(3.1)P(x_i)=y_i\qquad (i=0,1,\cdots,n) \tag{3.1}

插值多项式：

P(x)=a0+a1x+⋯+anxn(3.2)P(x)=a_0+a_1x+\cdots+a_nx^n\tag{3.2}

确定曲线（插值函数）使其通过插值点，用它近似表示数据的理想函数。

多项式插值

给定n+1n+1个点，求次数不超过nn的多项式。 aia_i为未知数

{a0+a1x0+⋯+anx0n=y0a0+a1x1+⋯+anx1n=y1...a0+a1xn+⋯+anxnn=yn\begin{cases} a_0+a_1x_0+\cdots+a_nx_0^n=y_0\\ a_0+a_1x_1+\cdots+a_nx^n_1=y_1\\ ...\\ a_0+a_1x_n+\cdots+a_nx^n_n=y_n \end{cases}

化成范德蒙矩阵，其行列式值不为零，解唯一。

行列式值与解的关系 待补充

拉格朗日插值

线性插值

只有两个节（数据）点(xk,yk),(xk+1,yk+1)(x_k,y_k),(x_{k+1},y_{k+1})

P(x)=a0+a1x=L1(x)(线性插值)P(x)=a_0+a_1x=L_1(x)\tag{线性插值}

要求插值函数上的值与节点值严格一致，在图像上显示出来就是一条直线,其斜率为：

k=yk+1−ykxk+1−xkk=\frac{y_{k+1}-y_{k}}{x_{k+1}-x_{k}}

转化成：

L1(x)=yklk(x)+yk+1lk+1(x)(3.3)L_1(x)=y_kl_k(x)+y_{k+1}l_{k+1}(x)\tag{3.3}

其中：

lk(x)=x−xk+1xk−xk+1lk+1(x)=x−xkxk+1−xkl_k(x)=\frac{x-x_{k+1}}{x_k-x_{k+1}}\qquad l_{k+1}(x)=\frac{x-x_{k}}{x_{k+1}-x_{k}}

lk(x),lk+1(x)l_k(x),l_{k+1}(x)为线性插值基函数，在节点上满足：

lk(xk)=1lk(xk+1)=0lk+1(xk)=0lk+1(xk+1)=1l_k(x_k)=1\qquad l_k(x_{k+1})=0 \qquad l_{k+1}(x_k)=0 \qquad l_{k+1}(x_{k+1})=1\qquad

习题1： (1,0) (2,4)

l1=x−12−1=x−1l2=x−21−2=2−xl_1=\frac{x-1}{2-1}=x-1\qquad l_2=\frac{x-2}{1-2}=2-x

l1(x)=0+4x−4=4x−4l_1(x)=0+4x-4=4x-4

抛物插值

有三个节点时：

P(x)=a0+a1x+a2x2=L2(x)(抛物插值)P(x)=a_0+a_1x+a_2x^2=L_2(x)\tag{抛物插值}

进而化简得到基函数：下面节点减去所有其他的点，上面用x替换节点

L2(x)=yk−1lk−1(x)+yklk(x)+yk+1lk+1(x)(3.4)L_2(x)=y_{k-1}l_{k-1}(x)+y_{k}l_{k}(x)+y_{k+1}l_{k+1}(x) \tag{3.4}

其中：

lk−1(x)=(x−xk)(x−xk+1)(xk−1−xk)(xk−1−xk+1)l_{k-1}(x)=\frac{(x-x_k)(x-x_{k+1})}{(x_{k-1}-x_{k})(x_{k-1}-x_{k+1})}

基函数序数与x序数的关系：相等为1，不相等为0。

lk−1(xk−1)=1lk−1(xk)=0lk−1(xk+1)=0...l_{k-1}(x_{k-1})=1\qquad l_{k-1}(x_{k})=0\qquad l_{k-1}(x_{k+1})=0\qquad ...

习题2： (0,1) (1,2) (2,3)

l1(x)=(x−1)(x−2)(0−1)(0−2)...l_1(x)=\frac{(x-1)( x-2)}{(0-1)(0-2)}...

L2(x)=l1(x)y1+...L_2(x)=l_1(x)y_1+...

拉格朗日插值多项式

当节点数n比3大时

Ln(xj)=∑k=0nyklk(xj)=yj(j=0,1,⋯ ,n)(3.5)L_n(x_j)=\sum_{k=0}^n y_kl_k(x_j)=y_j\qquad (j=0,1,\cdots,n)\tag{3.5}

n次插值基函数性质：

lj(xk)={1k=j0k≠jl_j(x_k)= \begin{cases} 1 & k=j\\ 0 & k≠j\\ \end{cases}

lk(x)=(x−x0)⋯(x−xk−1)(x−xk+1)⋯(x−xn)(xk−x0)⋯(xk−xk−1)(xk−xk+1)⋯(xk−xn)l_{k}(x)=\frac{\left(x-x_{0}\right) \cdots\left(x-x_{k-1}\right)\left(x-x_{k+1}\right) \cdots\left(x-x_{n}\right)}{\left(x_{k}-x_{0}\right) \cdots\left(x_{k}-x_{k-1}\right)\left(x_{k}-x_{k+1}\right) \cdots\left(x_{k}-x_{n}\right)}

引入记号：

ωn+1(x)=(x−x0)(x−x1)⋯(x−xn)\omega_{n+1}(x)=(x-x_{0})(x-x_{1})\cdots (x-x_{n})

所以：

ωn+1′(xk)=(xk−x0)⋯(xk−xk−1)(xk−xk+1)⋯(xk−xn)\omega^{'}_{n+1}(x_k)=(x_k-x_{0})\cdots (x_k-x_{k-1})(x_k-x_{k+1})\cdots (x_k-x_{n})

于是公式(3.5)可以化简为：

Ln(x)=∑k=0nykωn+1(x)(x−xk)ωn+1′(xk)(3.6)L_{n}(x)=\sum_{k=0}^{n} y_{k} \frac{\omega_{n+1}(x)}{\left(x-x_{k}\right) \omega_{n+1}^{\prime}\left(x_{k}\right)} \tag{3.6}

老师提供的python文件

Lagrange1.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as pylab
from scipy.interpolate import lagrange
from scipy.interpolate import spline
from scipy.optimize import leastsq

#分段线性插值闭包
def xianxingfenduan(xn, yn,x):
    temp = -1 
    #找出x所在的区间
    for i in range(1, len(xn)):
        if x <= xn[i]:
            temp = i-1
            break
        else:
            i += 1
         
    if temp == -1:
        return -100
        #插值
    result = (x-xn[temp+1])*yn[temp]/float((xn[temp]-xn[temp+1])) \
        + (x-xn[temp])*yn[temp+1]/float((xn[temp+1]-xn[temp]))
    return result

#X = [0, 1,4,9,16,25,36,49,64]   #插值点 xi
#Y = [0,1,2,3,4,5,6,7,8]  #插值 yi
#X = np.linspace(-5,5,21)
#Y = list(map(lambda x: 1/(1+x**2), X))
#X = [1.9,2,2.1,2.5,2.7,3.5,4,4.5,4.6,5,5.2,6,6.3,6.5,7.1,8,8.9,9,9.5,10]   #插值点 xi
#Y = [1.4,1.3,1.8,2.5,2.5,3,4,4.2,3.5,5.5,5,5.5,6.4,6,5.3,6.5,7,8,8.1,8.1]  #插值 yi  
X = [i for i in range(1,16)]   #插值点 xi
Y = [4,6.4,8,8.8,9.22,9.5,9.7,9.95,10.2,10.32,10.42,10.5,10.55,10.5,10.6]  #插值 yi  

#x = np.linspace(-5,5,2001)
x = np.linspace(1,15,2001)

#求分段线性插值
L = lagrange(X,Y)
y = list(map(lambda x: L(x), x))

#三次样条
s = list(map(lambda x: spline(X,Y,x), x))

#分段线性
f = list(map(lambda x: xianxingfenduan(X,Y,x), x))

#准确解
#y1 = list(map(lambda x: 1/(1+x**2), x))

# 最小二乘法
def error (p,x,y):                    # 拟合残差
    return np.polyval(p,x)-y 
p0 = [1,1,1,1,1]
para =leastsq(error, p0, args=(X,Y)) # 进行拟合
y_fitted = np.polyval(para[0],x) # 画出拟合后的曲线

#画图
#pylab.ylim(0,10)
pylab.plot(X, Y, 'ro',label='dots')
#pylab.plot(x, y, 'b-')
#pylab.plot(x, s, 'r-')
#pylab.plot(x, f, 'g-')
#pylab.plot(x, y1, 'y-')
#pylab.plot(x, y_fitted)
pylab.show()


pylab.xlim(-1,1)
pylab.ylim(0.8,1)
pylab.plot(X, Y, 'ro')
pylab.plot(x, y, 'b-')
pylab.plot(x, s, 'r-')
pylab.plot(x, f, 'g-')
pylab.plot(x, y1, 'y-')
pylab.show()
Lagrange.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
# -*- coding: utf-8 -*-

import numpy as np
import matplotlib.pyplot as pylab

def Lagr(X,Y,x):
    ans=0.0
    for i in range(len(Y)):
        t=Y[i]
        for j in range(len(Y)):
            if i !=j:
                t*=(x-X[j])/(X[i]-X[j])
        ans +=t
    return ans

X = [0, 1,4,9,16,25,36,49,64]   #插值点 xi
Y = [0,1,2,3,4,5,6,7,8]  #插值 yi

x = 0.3367      #求近似值点

L = Lagr(X, Y, x)
print(L)


import math
err = math.sin(x)-L
print(err)

#求分段线性插值
x = np.linspace(0,64,129)
y = list(map(lambda x: Lagr(X, Y, x), x))

#准确解
y1 = list(map(lambda x: np.sqrt(x), x))

#画图
pylab.plot(X, Y, 'ro')
pylab.plot(x, y, 'b-')
pylab.plot(x, y1, 'y-')
pylab.show()

插值余项与误差估计

罗尔定理

余项定义：

Rn(x)=f(x)−Ln(x)R_n(x)=f(x)-L_n(x)

Rn(x)=f(x)−Ln(x)=f(n+1)(ξ)(n+1)!ωn+1(x)(3.6)R_{n}(x)=f(x)-L_{n}(x)=\frac{f^{(n+1)}(\xi)}{(n+1) !} \omega_{n+1}(x) \tag{3.6}

牛顿插值

为何要用牛顿插值法：

拉格朗日插值与插值节点的依赖性太大，用牛顿插值法可以逐次进行插值。

将拉格朗日线性插值看作零次插值的修正；线性插值看作抛物插值的修正；…

差商（微商的离散形式。）：

f[xi,xi+1,⋯ ,xi+k]=f[xi+1,xi+2,⋯ ,xi+k]−f[xi,xi+1,⋯ ,xi+k−1]xi+k−xif\left[x_{i}, x_{i+1}, \cdots, x_{i+k}\right]=\frac{f\left[x_{i+1}, x_{i+2}, \cdots, x_{i+k}\right]-f\left[x_{i}, x_{i+1}, \cdots, x_{i+k-1}\right]}{x_{i+k}-x_{i}}

差商的性质

  1. k阶差商

    f[x0,x1,⋅⋅⋅,xk]=∑j=0kf(xj)ωk+1′(xj)f[x_0, x_1, · · · , x_k] =\sum_{j=0}^k\frac{f(x_j)}{\omega^{'}_{k+1}(x_j)}

  2. 差商与其所含节点的排列次序无关，即

    f[xi,xi+1]=f[xi+1,xi]f[x_i , x_{i+1}] = f[x_{i+1}, x_i]

  3. 设f(x)f(x)在包含互异节点x0,x1,⋅⋅⋅,xnx_0, x_1, · · · , x_n的闭区间[a,b][a, b]上有n阶导数，则

    f[x0,x1,⋅⋅⋅,xn]=fn(ξ)n!ξ∈(a,b)f[x_0, x_1, · · · , x_n] =\frac{f^{n}(\xi)}{n!}\quad \xi \in(a,b)

差商表

牛顿插值

2021-10-16 12:31:19

不学了！！！我要去收拾东西，晚上再坐上32个小时的火车，去向不可知的远方。楼下母亲快把饭做好了，饿！饿！饿！

回到学校再来过！！！

Pn(x)=f(x0)+f[x0,x1](x−x0)+f[x0,x1,x2](x−x0)(x−x1)+⋯+f[x0,x1,⋯ ,xn](x−x0)(x−x1)⋯(x−xn−1)R~n(x)=f[x,x0,x1,⋯ ,xn](x−x0)(x−x1)⋯(x−xn)\begin{aligned} P_{n}(x)=& f\left(x_{0}\right)+f\left[x_{0}, x_{1}\right]\left(x-x_{0}\right)+f\left[x_{0}, x_{1}, x_{2}\right]\left(x-x_{0}\right)\left(x-x_{1}\right)+\cdots \\ &+f\left[x_{0}, x_{1}, \cdots, x_{n}\right]\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n-1}\right) \\ \tilde{R}_{n}(x)=& f\left[x, x_{0}, x_{1}, \cdots, x_{n}\right]\left(x-x_{0}\right)\left(x-x_{1}\right) \cdots\left(x-x_{n}\right) \end{aligned}

称Nn(x)N_n(x)为nn次牛顿插值多项式,R~n(x)R̃_n (x)为牛顿型插值余项。

  • 插值多项式存在且唯一：Pn(x)≡Ln(x)P_n(x)≡ L_n(x)

差分形式

步长：等距节点xk=x0+kh (k=0,1,…,n)x_k=x_0+kh\ (k=0,1,\dots,n)

称Δfk=fk+1−fk\Delta f_k=f_{k+1}-f_k为一阶差分

Δfk=fk+1−fk\Delta f_{k}=f_{k+1}-f_{k}

  • 差分与差商的关系：

f[x0,x1,⋯ ,xk]=Δky0k!hk,k=1,2,⋯ ,nf[xn,xn−1,⋯ ,xn−k]=∇kynk!hk,k=1,2,⋯ ,n\begin{aligned} &f\left[x_{0}, x_{1}, \cdots, x_{k}\right]=\frac{\Delta^{k} y_{0}}{k ! h^{k}}, \quad k=1,2, \cdots, n \\ &f\left[x_{n}, x_{n-1}, \cdots, x_{n-k}\right]=\frac{\nabla^{k} y_{n}}{k ! h^{k}}, \quad k=1,2, \cdots, n \end{aligned}

  • 差分与导数的关系

    Δnyo=hnf(n)(ξ) ξ∈(x0,xn)\Delta^n y_o=h^nf^{(n)}(\xi)\ \xi \in (x_0,x_n)

  • 牛顿向前插值多项式

    Pn(x0+th)=f(x0)+tΔy0+t(t−1)2!Δ2y0+⋯+t(t−1)⋯(t−n+1)n!Δny0P_{n}\left(x_{0}+t h\right)=f\left(x_{0}\right)+t \Delta y_{0}+\frac{t(t-1)}{2 !} \Delta^{2} y_{0}+\cdots+\frac{t(t-1) \cdots(t-n+1)}{n !} \Delta^{n} y_{0}

    其余项为：

    Rn(x0+th)=t(t−1)⋯(t−n+1)n!hn+1f(n+1)(ξ)ξ∈(x0,xn)R_n(x_0+th)=\frac{t(t-1) \cdots(t-n+1)}{n !}h^{n+1}f^{(n+1)}(\xi)\quad\xi \in (x_0,x_n)