三次樣條插值算法詳解課件

1樣條函數(shù)的定義定義4.1

設(shè)區(qū)間[a,b]上給定一個節(jié)點劃分a=x0<x1<……<xn-1<xn=b 如果存在正整數(shù)k使得[a,b]上的分段函數(shù)s(x)滿足如下兩條：(1)在[a,b]上有直到k-1階連續(xù)導(dǎo)數(shù)。(2)在每個小區(qū)間[xi,xi+1]上是次數(shù)不大于k的多項式。則稱分段函數(shù)s(x)是以(2.6)為節(jié)點集的k次樣條函數(shù)。1樣條函數(shù)的定義定義4.1設(shè)區(qū)間[a,b]上給定一個節(jié)點劃2三次樣條插值函數(shù)的定義并且關(guān)于這個節(jié)點集的三次樣條函數(shù)s(x)滿足插值條件：則稱這個三次樣條函數(shù)s(x)為三次樣條插值函數(shù)。2三次樣條插值函數(shù)的定義并且關(guān)于這個節(jié)點集的三次樣條函數(shù)s(3三次樣條插值函數(shù)的邊界條件插值條件：連續(xù)性條件：一階導(dǎo)數(shù)連續(xù)條件：二階導(dǎo)數(shù)連續(xù)條件：3三次樣條插值函數(shù)的邊界條件插值條件：連續(xù)性條件：一階導(dǎo)數(shù)連4（1）因為s(x)在每個小區(qū)間上是一個次小于三次的多項式，故有四個未知系數(shù)；（2）因為s(x)有n分段，從而共有4n個未知系數(shù)！（3）但插值條件與樣條條件僅給出4n-2個條件，無法定出4n個未知系數(shù)，還差2個條件！這2個條件我們用邊界條件給出！

4（1）因為s(x)在每個小區(qū)間上是一個次小于三次的多項式，5通常我們對插值多項式在兩端點的狀態(tài)加以要求也就是所謂的邊界條件:第一邊界條件：由區(qū)間端點處的一階導(dǎo)數(shù)給出即 5通常我們對插值多項式在兩端點的狀態(tài)加以要求也就是第一邊界條6第二邊界條件：由區(qū)間端點處的二階導(dǎo)數(shù)給出即

特殊情況為自然邊界條件：由區(qū)間端點處的二階導(dǎo)數(shù)恒為0給出即

6第二邊界條件：由區(qū)間端點處的二階導(dǎo)數(shù)給出即 特殊情況為自7這樣三次樣條插值問題就分成三類！其實不止這三類!第三類又稱周期邊界條件：由區(qū)間端點處的函數(shù)值或?qū)?shù)值滿足周期條件給出

7這樣三次樣條插值問題就分成三類！其實不止這三類!第三類又稱8樣條函數(shù)的例子容易驗證：是滿足如下數(shù)據(jù)的第一類邊界樣條插值問題解：x0123y0000y’108樣條函數(shù)的例子容易驗證：是滿足如下數(shù)據(jù)的第一類邊界樣條插值9樣條函數(shù)的例子9樣條函數(shù)的例子10通常有三轉(zhuǎn)角法、三彎矩法、B樣條基函數(shù)法。三次樣條插值函數(shù)的求法這三種方法的基本思想是類似的，都是通過待定某些參數(shù)來確定插值函數(shù)，但肯定不是待定4n個參數(shù)。而是利用已知條件將待定參數(shù)減小到最少。比如：待定一階導(dǎo)數(shù)、待定二階導(dǎo)數(shù)、采用基函數(shù)方法來確定插值函數(shù)。10通常有三轉(zhuǎn)角法、三彎矩法、B樣條基函數(shù)法。三次樣條插值函11三轉(zhuǎn)角法:待定一階數(shù)為了確定三次樣條插值函數(shù)的表達式S(x)，我們采用待定系數(shù)法來求解，我們待定什么系數(shù)呢？考慮到帶一階導(dǎo)數(shù)的分段三次Hermite插值多項式11三轉(zhuǎn)角法:待定一階數(shù)為了確定三次樣條插值函數(shù)的表達式S12我們采用待定一階導(dǎo)數(shù)的方法即設(shè)因為分段三次Hermite插值多項式已經(jīng)至少是一階連續(xù)可導(dǎo)了，為了讓它成為三次樣條函數(shù)只需確定節(jié)點處的一階導(dǎo)數(shù)使這些節(jié)點處的二階導(dǎo)數(shù)連續(xù)即可！12我們采用待定一階導(dǎo)數(shù)的方法即設(shè)因為分段三次Hermite131314由于在內(nèi)部節(jié)點處二階導(dǎo)數(shù)連續(xù)條件：整理化簡后得：14由于在內(nèi)部節(jié)點處二階導(dǎo)數(shù)連續(xù)條件：整理化簡后得：15稱為三轉(zhuǎn)角法基本方程組以上推導(dǎo)還沒有考慮邊界條件！針對不同類型的三次樣條問題，就可以導(dǎo)出不同的方程組！15稱為三轉(zhuǎn)角法基本方程組以上推導(dǎo)還沒有考慮邊界條件！針對不16第一類三次樣條插值問題方程組基本方程組化為n-1階方程組由于已知：化為矩陣形式16第一類三次樣條插值問題方程組基本方程組化為n-1階方程組17這是一個嚴(yán)格對角占優(yōu)的三對角方程組，用追趕法可以求解！17這是一個嚴(yán)格對角占優(yōu)的三對角方程組，18第二類三次樣條插值問題的方程組由于已知：故得：18第二類三次樣條插值問題的方程組由于已知：故得：19稍加整理得聯(lián)合基本方程組得一個n+1階三對角方程組，化成矩陣形式為：仍然是嚴(yán)格對角占優(yōu)19稍加整理得聯(lián)合基本方程組得一個n+1階三對角方程組，20第三類樣條插值問題的方程組立即可得下式：由于：20第三類樣條插值問題的方程組立即可得下式：由于：21其中：聯(lián)合基本方程得一個廣義三對角或周期三對角方程組：這個方程組的系數(shù)矩陣仍然是嚴(yán)格對角占優(yōu)陣！21其中：聯(lián)合基本方程得一個廣義三對角或周期三對角方程組：這22求解這些不同類型的樣條插值問題的方程組，我們可得所要待定的一階導(dǎo)數(shù)：稱為三次樣條插值問題三轉(zhuǎn)角公式！再代入S(x)的每一段表達式，就求得三次樣條函數(shù)的表達式!22求解這些不同類型的樣條插值問題的方程組，我們可得稱為三次23例1.對于給定的節(jié)點及函數(shù)值解:這是自然邊界條件下的樣條問題。23例1.對于給定的節(jié)點及函數(shù)值解:這是自然邊界條件下的樣24我們可以將上述計算列于表中：k0123xk1245yk1342mk????Mk0??0hk121*λk*2/31/3*μk*1/32/3*gk69/2-7/2-624我們可以將上述計算列于表中：k0123xk1245yk125由些得如下方程組：利用三轉(zhuǎn)角公式：25由些得如下方程組：利用三轉(zhuǎn)角公式：26同樣可以求得第三段表達式！26同樣可以求得第三段表達式！2727282829三彎矩法:待定二階導(dǎo)數(shù)選擇二階導(dǎo)數(shù)作為待定參數(shù):由于三次樣條S(x)是三次多項式,故它的二階導(dǎo)數(shù)是一次多項式,從而思考:(1)的原因?29三彎矩法:待定二階導(dǎo)數(shù)選擇二階導(dǎo)數(shù)作為待定參數(shù):由于三次303031從而推導(dǎo)出了三次樣條S(x)在第k個小區(qū)間[xk,xk+1]上的表達式為:它的系數(shù)都是用二階導(dǎo)數(shù)與函數(shù)值表示!31從而推導(dǎo)出了三次樣條S(x)在第k個小區(qū)間[xk,xk+32對所有中間節(jié)點xk,k=1,2,…,n-1,左邊小區(qū)間與右邊小區(qū)間上的三次多項式的一階導(dǎo)數(shù)應(yīng)當(dāng)連續(xù)!確定二階導(dǎo)數(shù)32對所有中間節(jié)點xk,k=1,2,…,n-1,左邊小區(qū)間與33三彎矩法基本方程注意到這個基本方程只包括了n-1個方程!但卻有n個二階導(dǎo)數(shù)需要待定,這是一個欠定方程組,還需要根據(jù)邊界條件再確定兩個方程!33三彎矩法基本方程注意到這個基本方程只包括了n-1個方程!34曲率調(diào)整樣條這種樣條的邊界條件是已知兩端點的二階導(dǎo)數(shù)值!這樣從三彎矩基本方程可以導(dǎo)數(shù)確定其它n-2個待定參數(shù)的方程組:34曲率調(diào)整樣條這種樣條的邊界條件是已知兩端點的二階導(dǎo)數(shù)值!35自然樣條這種樣條的邊界條件是:已知兩端點的二階導(dǎo)數(shù)值為0!這樣從三彎矩基本方程可以導(dǎo)數(shù)確定其它n-2個待定參數(shù)的方程組:35自然樣條這種樣條的邊界條件是:已知兩端點的二階導(dǎo)數(shù)值為036固支樣條這種樣條的邊界條件是:已知兩端點的一階導(dǎo)數(shù)值!根據(jù)前面推導(dǎo)過程中得到的樣條函數(shù)S(x)的一階導(dǎo)數(shù)的表達式(2.11),得方程36固支樣條這種樣條的邊界條件是:已知兩端點的一階導(dǎo)數(shù)值!根37固支樣條這樣從三彎矩基本方程可以導(dǎo)數(shù)確定n個待定參數(shù)的方程組:37固支樣條這樣從三彎矩基本方程可以導(dǎo)數(shù)確定n個待定參數(shù)的方38非扭結(jié)樣條這種樣條的邊界條件是:要求樣條S(x)在開始的兩個小區(qū)間[x0,x1],[x1,x2]上的三階導(dǎo)數(shù)相同,在最后兩個小區(qū)間[xn-2,xn-1],[xn-1,xn]上的三階導(dǎo)數(shù)相同.對表達式(2.9)再求一次導(dǎo)數(shù)得方程38非扭結(jié)樣條這種樣條的邊界條件是:要求樣條S(x)在開始的39非扭結(jié)樣條再由三彎矩基本方程,可得39非扭結(jié)樣條再由三彎矩基本方程,可得40周期樣條這種樣條的邊界條件是:要求樣條S(x)及其導(dǎo)數(shù)是以區(qū)間長度xn-x0為周期的函數(shù)即這些條件可以確定如下兩個方程:40周期樣條這種樣條的邊界條件是:要求樣條S(x)及其導(dǎo)數(shù)是41再由三彎矩基本方程,可得周期樣條41再由三彎矩基本方程,可得周期樣條42在Matlab中數(shù)據(jù)點稱之為斷點。如果三次樣條插值沒有邊界條件，最常用的方法，就是采用非扭結(jié)（not-a-knot）條件。這個條件強迫第1個和第2個三次多項式的三階導(dǎo)數(shù)相等，對最后一個和倒數(shù)第2個三次多項式也做同樣地處理。Matlab中三次樣條插值也有現(xiàn)成的函數(shù)：

y=interp1(x0,y0,x,'spline')；

y=spline(x0,y0,x)；

pp=csape(x0,y0,conds), pp=csape(x0,y0,conds,valconds)，

y=ppval(pp,x)。 其中x0,y0是已知數(shù)據(jù)點，x是插值點，y是插值點的函數(shù)值。對于三次樣條插值，我們提倡使用函數(shù)csape，csape的返回值是pp形式，要求插值點的近似函數(shù)值，必須調(diào)用函數(shù)ppval。MATLAB中三次樣條函數(shù)法實現(xiàn)42在Matlab中數(shù)據(jù)點稱之為斷點。如果三次樣條插值沒有邊43pp=csape(x0,y0,conds,valconds) conds指定插值的邊界條件，其值可為：

'complete' 邊界為一階導(dǎo)數(shù)，一階導(dǎo)數(shù)的值在valconds參數(shù)中給出。

'not-a-knot' 非扭結(jié)條件

'periodic' 周期條件

'second' 邊界為二階導(dǎo)數(shù)，二階導(dǎo)數(shù)的值在valconds參數(shù)中給出，若忽略valconds參數(shù)，二階導(dǎo)數(shù)的缺省值為[0,0]。MATLAB中三次樣條函數(shù)法實現(xiàn)43pp=csape(x0,y0,conds,valcond44例2：第一邊界條件的例題x1245y1342y’17/8-19/8clear;x=[1,2,4,5];y=[1,3,4,2];pp=csape(x,y,'complete',[17/8,-19/8]);pp.coefsMATLAB代碼44例2：第一邊界條件的例題x1245y1342y’17/845第一邊界條件的例題pp=form:'pp'breaks:[1245]coefs:[3x4double]pieces:3order:4dim:1pp.coefs-0.12500 2.12501.0000-0.1250-0.37501.75003.00000.3750-1.1250-1.25004.0000分段多項式結(jié)構(gòu)與系數(shù)矩陣45第一邊界條件的例題pp=分段多項式結(jié)構(gòu)與系數(shù)矩陣46如下代碼求解上述樣條問題：x-4-3-2-101234y00.151.122.362.361.460.490.060y’00例3：自然邊界條件例題46如下代碼求解上述樣條問題：x-4-3-2-101234y47x=[-4 -3 -2 -1 0 1 2 3 4];y=[0 0.15 1.12 2.36 2.36 1.46 0.49 0.06 0];pp=csape(x,y,'second');xx=-4:0.01:4;yy=ppval(pp,xx);holdon;plot(x,y,'ok');plot(xx,yy,'k-');holdoff;MATLAB程序47x=[-4 -3 -2 -1 0 1 2 3 4];MA48pp=form:'pp'breaks:[-4-3-2-101234]coefs:[8x4double]pieces:8order:4dim:1分段多項式結(jié)構(gòu)：48pp=分段多項式結(jié)構(gòu)：49

pp.coefs=0.180856038291610.00000000000000-0.030856038291610-0.084280191458030.542568114874820.511712076583210.15000000000000-0.393735272459500.289727540500741.344007731958761.120000000000000.14922128129602-0.891478276877760.742256995581742.360000000000000.13685014727540-0.44381443298969-0