三次样条插值
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
其中
C0
5, 384
C 1,
1 24
C
2
1 8
.
Cubic Spline
注:提高精度只须增加节点, 而无须提高样条阶数。 在实际应用中,不仅常用S(x) [式(3.1)]计算 f(x) 的
近似值,而且常用 S (x) [式(3.2)]近似计算 f (x) .
稳定性:只要方程组系数矩阵为SDD阵,保证数
x
yi 1
M
i
h2
1 i
6
x
xi hi
x [ xi , xi1], i 0,1,L n 1
这是三次样条插值函数的表达式,当求出Mi 后, S(x)就完全确定.
为了求 Mi ,需要利用S(x)在内结点处一阶导数连续的条件, 由 上式可得
S(x)
Mi
( xi1 2hi
x)2
Mi1
( x xi )2 2hi
三次样条插值函数的构造与Matlab实现
i
hi 1 hi1 hi
,
i
1i
hi hi1 hi
,
i
6
yi1 hi
yi
yi yi1 hi 1
(hi1 hi )
这里有 n+1 个未知数,n1个方程,还需增加 2 个方程。
对于I 类边界条件,S( x0 ) M0 , S( xn ) Mn
2M 0 0 M 1 0
M
n
M
n1
( xi1 6hi
x)3
Mi1
( x xi )3 6hi
C1 x C2
利用插值条件
S( xi ) yi , S( xi1 ) yi1
Cubic Spline
定出积分常数,可以得到
S(x)
Mi
( xi1 6hi
x)3
Mi1
(x
xi )3 6 hi
(3.1)
yi
Mi hi2 6
xi1 hi
§3 三次样条插值 /* Cubic Spline Interpolation */
背景
代数插值
Hermite插值
高次插值出现龙格现象
分段插值 但分段线性插值在节点处不一定光滑
分段Hermite插值 但导数值不容易得到 三次样条插值(先由函数值确定导数值,再由分段 Hermite插值解决问题) 应用最为广泛
➢ 三次样条插值问题
Cubic Spline
定义 设 a x0 x1 ... xn b 。三次样条函数 S( x) C 2[a, b] ,
且在每个[xi , xi1]上为三次多项式 /* cubic polynomial */。若它同 时还满足 S( xi ) f ( xi ), (i 0,1,L , n) ,则称 S(x) 为 f(x) 在结点 xj ( j = 0, 1, …, n) 上的三次样条插值函数 .
共有 3(n-1) 个条件。再加上 n+1 个插值条件,共有4n-2 个
条件因。此,还需要2个条件才能确定S(x)。通常在区间
端点 a = x0 和 b = xn 上各加一个条件(称为边界条件),可 根据实际问题的要求给定。
常用边界条件 /* boundary conditions */
(1)已知两端的二阶导数值,即
0
0
n1
2
n1
M
n1
n
1
0 0 0 0 n 2 M n n
注意其系数矩阵是按行严格对角占优阵,故它有唯一解。
可用追赶法求解。
➢ 三次样条插值函数的误差估计
Cubic Spline
在实际应用中,如果不需要规定内节点处的一阶导 数值,那么使用三次样条插值函数会得到很好的效果。 三次样条插值函数 S(x) 不仅在内节点处的二阶导数是连 续的,而且 S(x) 逼近 f(x) 具有很好的收敛性,也是数值 稳定的。由于误差估计与收敛性定理的证明比较复杂, 下面只给出误差估计的结论。
定理5.5
设函数
f ( x) C 4[a, b],
记
M4
max a xb
f (4)( x) ,
h
max (
0 i n1
xi 1
xi
),
则
x [a,b],
满足 I 类或 II 类边界条件的
三次样条插值函数 S(x) 有估计式
f S C h M (k)
( x)
(k)( x)
4k k
,
4
k 0,1,2
[
x
,
i
x上i1]可以写成
S( x) ai x3 bi x2 ci x di , i 0,1,L , n 1,
共有 4n 个待定参数。S(x) 在[a, b]上二阶导数连续,故在内
结点 xi (i 1, 2,L处,应n -满1)足连续性条件
Sk xi 0 Sk ( xi 0), k 0,1, 2,
2
n n
0 0, 0 2M0 n 0, n 2Mn
Cubic Spline
即得关于Mi (i = 0, 1, …, n)的 n+1 元线性方程组
2 0 0 0 0 0 M0 0
1 2 1 0
0
0
M1
1
0
0
2 2
0O
2
O
0 O
0 0
M2 M
2
M
0
注:三次样条与分段 Hermite 插值的根本区别在于S(x)自 身光滑,不需要知道 f 的导数值(除了在2个端点可能需 要);而Hermite插值依赖于f 在所有插值点的导数值。
f(x)
H(x)
S(x)
ቤተ መጻሕፍቲ ባይዱ
Cubic Spline
3次样条插值函数S( x是) 否存在唯一? 如何计算?误差估计?
三次样条插值函数是分段三次多项式,在每个小区间
Cubic Spline
S( x0 )
f
0
M0 ,
S( xn )
f
n
Mn.
(I 类)
其特殊情况为
S( x0 ) 0, S( xn ) 0,
(自由边界)
对应的样条函数称为自然样条 /* Natural Spline */.
(2)已知两端的一阶导数值,即
S( x0 ) f0 m0 ,
S( xn )
Cubic Spline
➢ 三弯矩法 /* method of bending moment */
三次样条插值函数 S(x) 可以有多种表达式,有时用二阶
导数值 S( xi ) Mi (i 表0示,1,时L ,, n使) 用更方便。Mi 在力学上解释 为细梁在 xi 处的弯矩,并且得到的弯矩与相邻两个弯矩有关, 故称用 Mi 表示 S(x) 的算法称为三弯矩法。
Newton Ln(x),只是形式不同;渐增节点或节点等距 时方便处理.
Hermite: 给出 yi 及 yi . Spline:分段低次, 自身光滑, f 的导数只在边界给出.
yi = interp1(x,Y,xi,method) :用指定的算法计算插值; ’linear’:线性插值(缺省方式),直接完成计算; ’pchip’:分段三次Hermite插值。对于该方法,命令 interp1 调用函数pchip,用于对向量x 与y 执行分段三次 内插值。该方法保留单调性与数据的外形; ’cubic’:与’pchip’操作相同; ’spline’:三次样条函数插值。对于该方法,命令 interp1 调用函数spline、ppval、mkpp、umkpp。这些命 令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值;注意它默认使用的是 ‘not-a-knot’边界条件,也就是第一个点的三次导数和 第二点的三次导数一样;最后一个点的三次导数和倒数 第一个点一样。当Y=[df1,y,df2]时,表示第一点和最后 一个点的一阶导数分别为df1,df2。 pp = csape(x,y,conds): 计算在各种边界条件下的三次样条插值。 >> help csape
hi 1 (i 0,1, 2)
i i
i
i
0
1
6
1 12
12
-3
2 12
12
36
31
-78
得方程组 2 1 0 0 M0 6
1 2
0 0
2 12 0
12 2 1
0 12 2
M1 M2 M3
3 36 78
Cubic Spline
解得
M0
28 3
28 M1 3
106 M2 3
yi1 hi
yi
hi 6
Mi1 Mi
(3.2)
Cubic Spline
由 S(x) 在内结点处一阶导数连续,即
Si( xi 0) Si1( xi 0), i 1, 2,L , n 1
可得关于参数 M j 的方程组,即三弯矩方程的形式为
其中
i Mi1 2Mi i Mi1 i , i 1, 2,L , n 1
由于S(x) 在区间 [xi , xi1] (i 0,1上,L是, n3次1)多项式, 故 S (xi )
在
上是1[ x次i , x多i项1 ] 式, 可表示为
S( x)
Mi
xi1 x hi
Mi1
x xi hi
,
其中 hi = xi+1 - xi . 对S (xi ) 积分两次得
S(x)
Mi
f
n
mn .
(II类)
(3)周期边界条件
S k ( x0 ) S k ( xn ), k 0,1, 2.
(III 类)
此时,对函数值有周期条件 f ( x0 ) f ( xn ).
Cubic Spline 由boundary conditions 唯一确定。
定理 三次样条插值问题的解存在且唯一。
值稳定.
另有三转角法得到样条函数,即设 S (xi) = mi,则易
知[xi, xi+1 ]上的S(x) 就是Hermite函数. 再利用S(x)的 连续性,可导出关于mi 的方程组,加上边界条件即 可解.
Sketch of the Algorithm: Cubic Spline
① 计算 γi , αi , βi ;
例: x = -3:3; y = [-1 -1 -1 0 1 1 1]; t = -3:.01:3; plot(x,y,'o',t,[pchip(x,y,t); spline(x,y,t)]) legend('data','pchip','spline',4)
由于Lagrange插值可能不收敛,所以工程 中很少使用。matlab中没有专门的Lagrange插值 函数。可以自己编写一个Lagrange插值函数: function y=lagrange(x0,y0,x)
② 计算 Mi (追赶法等) ;
③ 找到 x 所在区间 ( 即找到相应的 i ) ;
④ 由该区间上的 S(x) 算出 f(x) 的近似值。
例 由函数表
Cubic Spline
j
0
1
2
3
xj
0
1
1
3
yj
0
2
3
16
求满足边界条件 y(0) 1, y(3) 0的三次样条插值函数。
解:用三弯矩方程(第二种边界条件)计算得
170 M3 3
将此解代入式 (3.1) 即得
x 11x2 14x 3
3
x 0,1
S
(
x)
1 24x3 91x2 108x 35
x 1, 2
3
1
46x3 329x2 732x 525
x 2, 3
3
HW: 习题 #24
插值法小结
Lagrange Ln(x) : 给出 y0 … yn,选基函数 li(x),其次数为 节点数 –1.
举例: ➢ 汽车、船的外形设计,流体力学等要求流线型(光滑)
➢ 木样条的来源
样条是绘图员用于描绘光滑曲线的由一些易弯曲材料制成的窄 条或棒条.在绘制需要通过某点的光滑曲线时,对它在这些点的位置 上“压铁”,它就被强制通过或接近图表上确定的描绘点.“样条函数” 这个术语意在点出这种函数的图象与机械样条画出的曲线很像.