數(shù)值分析第二章 插值

1、第二章 插值方法/* Interpolation */Interpolation_introduction1 問題提出問題提出函數(shù)逼近函數(shù)逼近 / /* *problem formulation-function approximationproblem formulation-function approximation* */ /用用Interpolation_introduction函數(shù)逼近的方法有很多，例如函數(shù)逼近的方法有很多，例如Taylor級(jí)數(shù)，級(jí)數(shù)，F(xiàn)ourier級(jí)數(shù)，有限元方法、邊界元方法，小級(jí)數(shù)，有限元方法、邊界元方法，小波分析等，大學(xué)科叫波分析等，大學(xué)科叫逼近論逼近論。本書

2、討論連續(xù)函數(shù)的逼近，主要介紹本書討論連續(xù)函數(shù)的逼近，主要介紹插值法插值法(chapter 2)和和最佳一直逼近、最小平方逼近最佳一直逼近、最小平方逼近離散數(shù)據(jù)擬合離散數(shù)據(jù)擬合(chapter 3)Interpolation_introduction插值節(jié)點(diǎn)插值節(jié)點(diǎn)插值條件插值條件-插值問題插值問題多項(xiàng)式插值是數(shù)值分析的基本工具，常用來(lái)計(jì)算被插函數(shù)多項(xiàng)式插值是數(shù)值分析的基本工具，常用來(lái)計(jì)算被插函數(shù)的近似的近似函數(shù)值函數(shù)值，零、極點(diǎn)零、極點(diǎn)，導(dǎo)數(shù)、積分導(dǎo)數(shù)、積分（第四章（第四章 數(shù)值積數(shù)值積分和數(shù)值微分），分和數(shù)值微分），解微分方程解微分方程（第五章）、（第五章）、積分方程積分方程插值插值多項(xiàng)式插

3、值多項(xiàng)式插值-polynomial interpolationProblem I. 給定給定y=f(x)的函數(shù)表的函數(shù)表, xi a,bniyxPiin,., 0,)(= = =求求 次數(shù)不超過次數(shù)不超過 n 的多項(xiàng)式的多項(xiàng)式 使得使得nnnxaxaaxP = =10)(條件：條件：無(wú)重合節(jié)點(diǎn)無(wú)重合節(jié)點(diǎn)，即，即jixx ji Interpolation intervalInterpolation conditionInterpolation polynomialInterpolation pointsInterpolation polynomial (2.1)(2.2)x0 x1x2x3x4x

4、Pn(x) f(x)多項(xiàng)式插值的幾何意義多項(xiàng)式插值的幾何意義Interpolation polynomial 求求插值多項(xiàng)式的唯一性插值多項(xiàng)式的唯一性 提問：提問：Problem I 中的中的Pn(x)是否存在？是否存在？ 若存在，是否唯一？如何求？若存在，是否唯一？如何求？Interpolation polynomial Interpolation polynomial 如何求？解線性方程如何求？解線性方程組（組（2.3）-待定系待定系數(shù)法數(shù)法Interpolation polynomial 2 拉格朗日多項(xiàng)式拉格朗日多項(xiàng)式 /* Lagrange Polynomial */niyxPiin

5、,., 0,)(= = =求求 n 次多項(xiàng)式次多項(xiàng)式 使得使得nnnxaxaaxP = =10)(條件：條件：無(wú)重合節(jié)點(diǎn)，即無(wú)重合節(jié)點(diǎn)，即jixx ji n = 1已知已知 x0 , x1 ; y0 , y1 ，求，求xaaxP101)( = =使得使得111001)(,)(yxPyxP= = =可見可見 P1(x) 是過是過 ( x0 , y0 ) 和和 ( x1, y1 ) 兩點(diǎn)的直線。兩點(diǎn)的直線。)()(0010101xxxxyyyxP- - - - = =101xxxx- - -010 xxxx- - -= y0 + y1l0(x)l1(x) = = =10)(iiiyxl稱為稱為拉氏

6、基函數(shù)拉氏基函數(shù) /* Lagrange Basis */，滿足條件滿足條件 li(xj)= ij /* Kronecker Delta */2 Lagrange Polynomialn 1希望找到希望找到li(x)，i = 0, , n 使得使得 li(xj)= ij ；然后令；然后令，則顯然有，則顯然有Pn(xi) = yi 。li(x)每個(gè)每個(gè) li 有有 n 個(gè)根個(gè)根 x0 xi xn= = =jiC0 = =- -nj i jxx)(- - - -inxxixxxxC0).().(ixl)( - -= =j i jxixiC)(1= =iixl1)(= = - - -= =njijj

7、ijixxxxxl0)()()( = = =niiinyxlxL0)()(Lagrange Polynomial與節(jié)點(diǎn)有關(guān)，而與與節(jié)點(diǎn)有關(guān)，而與 f無(wú)關(guān)無(wú)關(guān) = = =niinxlxP0)()(yi 基函數(shù)法（基函數(shù)法（n=1情形的推廣情形的推廣）2 Lagrange Polynomial定理定理 (唯一性唯一性) 滿足滿足 的的 n 階插值多階插值多項(xiàng)式是唯一存在的。項(xiàng)式是唯一存在的。niyxPii,., 0,)(= = =證明證明： ( 前面已利用前面已利用Vandermonde 行列式行列式論證論證)反證：若不唯一，則除了反證：若不唯一，則除了Ln(x) 外還有另一外還有另一 n 階多項(xiàng)

8、階多項(xiàng)式式 Pn(x) 滿足滿足 Pn(xi) = yi ?？疾炜疾?則則 Qn 的階數(shù)的階數(shù), )()()(xLxPxQnnn- -= = n而而 Qn 有有 個(gè)不同的根個(gè)不同的根n + 1x0 xn注注：若不將多項(xiàng)式次數(shù)限制為若不將多項(xiàng)式次數(shù)限制為 n ，則插值多項(xiàng)式，則插值多項(xiàng)式不唯一不唯一。例如例如 也是一個(gè)插值也是一個(gè)插值多項(xiàng)式，其中多項(xiàng)式，其中 可以是任意多項(xiàng)式?？梢允侨我舛囗?xiàng)式。= =- - = =niinxxxpxLxP0)()()()()(xp2 Lagrange Polynomial 插值余項(xiàng)插值余項(xiàng) /* Remainder */設(shè)節(jié)點(diǎn)設(shè)節(jié)點(diǎn))1( nf在在a , b內(nèi)存

9、在內(nèi)存在, 考察截?cái)嗾`差考察截?cái)嗾`差)()()(xLxfxRnn- -= =, baCfn bxxxan 10，且，且 f 滿足條件滿足條件 ,Rolles Theorem: 若若 充分光滑，充分光滑， ，則，則存在存在 使得使得 。)(x 0)()(10= = =xx ),(10 xx 0)(= = 推廣：推廣：若若0)()()(210= = = =xxx ),(),(211100 xxxx 使得使得0)()(10= = = = ),(10 使得使得0)(= = 0)()(0= = = =nxx 存在存在),(ba 使得使得0)()(= = nRn(x) 至少有至少有 個(gè)根個(gè)根n+1 = =

10、- -= =niinxxxKxR0)()()(任意固定任意固定 x xi (i = 0, , n), 考察考察 = =- - -= =niixtxKtRnt0)()()()( (x)有有 n+2 個(gè)不同的根個(gè)不同的根 x0 xn x),(, 0)()1(baxxn = = !)1()()()1(-nxKRxnn 注意這里是對(duì)注意這里是對(duì) t 求導(dǎo)求導(dǎo)= = - - - !)1)()()()1()1(nxKLfxnnxn !)1()()()1( = = nfxKxn = = - - = =niixnnxxnfxR0)1()(! ) 1()()( 2 Lagrange Polynomial1 La

11、grange Polynomial注：注： 通常不能確定通常不能確定 x , 而是估計(jì)而是估計(jì) , x (a,b) 將將 作為誤差估計(jì)上限。作為誤差估計(jì)上限。1)1()( nnMxf= = - - niinxxnM01|)!1(當(dāng)當(dāng) f(x) 為任一個(gè)次數(shù)為任一個(gè)次數(shù) n 的的多項(xiàng)式多項(xiàng)式時(shí)，時(shí)， , 可知可知 ，即插值多項(xiàng)式對(duì)于次數(shù)，即插值多項(xiàng)式對(duì)于次數(shù) n 的的多項(xiàng)多項(xiàng)式是式是精確精確的。的。0)()1( xfn0)( xRnQuiz: 給定給定 xi = i +1, i = 0, 1, 2, 3, 4, 5. 下面哪個(gè)是下面哪個(gè)是 l2(x)的圖像？的圖像？ y 0 - - - 1 0.

12、5 -0.5 1 2 3 4 5 6 x y 0 - - - 1 0.5 -0.5 1 2 3 4 5 6 x y 0 - - - 1 0.5 -0.5 1 2 3 4 5 6 x ABC 1 Lagrange Polynomial例：例：已知已知233sin,214sin,216sin= = = = 分別利用分別利用 sin x 的的1次、次、2次次 Lagrange 插值計(jì)算插值計(jì)算 sin 50 并估計(jì)誤差。并估計(jì)誤差。 解：解：0 x1x2x185500 =n = 1分別利用分別利用x0, x1 以及以及 x1, x2 計(jì)算計(jì)算4,610 =xx利用利用216/4/6/214/6/4/

13、)(1 - - - - - -= = xxxL這里這里)3,6(,sin)(,sin)()2( - -= = =xxxfxxf而而)4)(6(!2)()(,23sin21)2(1 - - -= = xxfxRxx00762. 0)185(01319. 01- - - - Rsin 50 = 0.7660444)185(50sin10 L0.77614外推外推 /* extrapolation */ 的實(shí)際誤差的實(shí)際誤差 - -0.010010.010013,421 = = =xx利用利用sin 50 0.76008, 00660. 018500538. 01 R內(nèi)插內(nèi)插 /* interpol

14、ation */ 的實(shí)際誤差的實(shí)際誤差 0.005960.00596內(nèi)插通常優(yōu)于外推。選擇內(nèi)插通常優(yōu)于外推。選擇要計(jì)算的要計(jì)算的 x 所在的區(qū)間的所在的區(qū)間的端點(diǎn)，插值效果較好。端點(diǎn)，插值效果較好。1 Lagrange Polynomialn = 223)()(21)()(21)()()(4363463464363646342 - - - - - - - - - - - - - - -= = xxxxxxxL)185(50sin20 L0.7654323cos21;)3)(4)(6(!3cos)(2 - - - - -= =xxxxxxR 00077. 018500044. 02 Rsin 5

15、0 = 0.76604442次插值的實(shí)際誤差次插值的實(shí)際誤差 0.000610.00061高次插值通常優(yōu)于高次插值通常優(yōu)于低次插值低次插值但絕對(duì)不是次數(shù)越但絕對(duì)不是次數(shù)越高就越好，嘿高就越好，嘿嘿嘿 When you start writing the program, you will find how easy it is to calculate the Lagrange polynomial.Oh yeah? What if I find the current interpolation not accurate enough? Then you might want to take

16、 more interpolating points into account.Right. Then all the Lagrange basis, li(x), will have to be re-calculated. Excellent point !We will come to discuss this problemnext time.3 逐次線性插值逐次線性插值 /* Lagrange Polynomial */ 0,1010011100,100100,202200,200200,20,10,1 20,1121:1Lagrange,;:1LagrangeIxx xxf xx

17、f xf xf xIxf xxxxxIxx xf xf xIxf xxxxxIxIxIxIxxxxx-=-=-=-，以 ，為節(jié)點(diǎn)的 次插值公式，實(shí)際上是過，點(diǎn)的直線，采用點(diǎn)斜式以 ，為節(jié)點(diǎn)的 次插值公式,同理有仿上，令易證 0,200,100,1 200,10010,100210,210,110,1 210,11110,111210,220,120,1 220,12210,222210,1 20122IxIxIxIxxxIxf xxxIxIxIxIxxxIxf xxxIxIxIxIxxxIxf xxxIxx x x-=-=-=-=-=-=-，由插值公式的唯一性知，是以 ，, 為節(jié)點(diǎn)的 次Lag

18、range插值多項(xiàng)式。以此 0,12,0,110,10,111101kkkkkkkkkIxIxIxIxxxxxx xxk-=-，類推，是以 ，, , 為節(jié)點(diǎn)的 次Lagrange插值多項(xiàng)式。 10,10,110,12,11kkkkkkkkkkx xx xIxIxIxxxxx-=-，實(shí)際上實(shí)際上，是對(duì)兩個(gè)低次插值的線性插值，這種通過低次插值再作是對(duì)兩個(gè)低次插值的線性插值，這種通過低次插值再作線性插值生成高次插值的方法稱為線性插值生成高次插值的方法稱為逐次線性插值逐次線性插值。 Aitken法法(按下表計(jì)算按下表計(jì)算) 0,12,0,110,10,1111kkkkkkkkIxIxIxIxxxxx-

19、=-，線性插值基函數(shù)線性插值基函數(shù) 0,0,1,0,1,2,0,1,2,3,00110,10220,20,1,2 kkkkkkxf xIxIxIxIxxf xxf xIxx xxf xIxIx- 1330,30,1,30,1,2,32440,40,1,40,1,2,40,1,2,3,43 x xxf xIxIxIxx xxf xIxIxIxIxx x-0123 x xx xx xx x-增加增加如果精度不構(gòu)，增加節(jié)點(diǎn)如果精度不構(gòu)，增加節(jié)點(diǎn)x4, ,同時(shí)表中增加一行，三同時(shí)表中增加一行，三角形斜邊上即為所要求的各次插值多項(xiàng)式。角形斜邊上即為所要求的各次插值多項(xiàng)式。k1k0k2k3k4 Nevil

20、le法法(按下表計(jì)算按下表計(jì)算) 1,1, ,11, ,11, ,1,200110,10221,20,1,2 kkkkkk kkk kkk kkxf xIxIxIxIxxf xxf xIxx xxf xIxIx- 0332,31,2,30,1,2,30443,42,3,40,1,2,40,1 x xxf xIxIxIxx xxf xIxIxIxI- ,2,3,401111 kkkkxx xx xx xx xx x-增加增加如果精度不構(gòu)，增加節(jié)點(diǎn)如果精度不構(gòu)，增加節(jié)點(diǎn)x4, ,同時(shí)表中增加一行，三角形斜邊上同時(shí)表中增加一行，三角形斜邊上即為所要求的各次插值多項(xiàng)式。即為所要求的各次插值多項(xiàng)式。k1

21、k0k1k1k1 11,0,110,10,1100kkkkkkIxIxIxIxxxxx-=-，HW: 用類似于前面的方法用類似于前面的方法構(gòu)造構(gòu)造Neville計(jì)算公式計(jì)算公式注：注：Atkin方法和方法和Neville方法與方法與Lagrange公式相比，當(dāng)公式相比，當(dāng)需要增加節(jié)點(diǎn)時(shí)，很容易由低次插值構(gòu)造高次插值，而需要增加節(jié)點(diǎn)時(shí)，很容易由低次插值構(gòu)造高次插值，而Lagrange插值公式中，每個(gè)基函數(shù)都需要作適當(dāng)變化。插值公式中，每個(gè)基函數(shù)都需要作適當(dāng)變化。誤差估計(jì)：誤差估計(jì)：由插值多項(xiàng)式的存在唯一性知，仍有由插值多項(xiàng)式的存在唯一性知，仍有= = - - = =niixnnxxnfxR0)1

22、()(! ) 1()()( 但這里可采用一種更簡(jiǎn)便的方法。當(dāng)?shù)@里可采用一種更簡(jiǎn)便的方法。當(dāng)f (n+1)(x)在插在插值區(qū)間變化不大時(shí)，設(shè)值區(qū)間變化不大時(shí)，設(shè)f (n+1)(x)L，則有則有 0,110-10,12,0-2-!-!kkkkkkLf xIxx xx xkLf xIxx xx xx xk-， 10,110,12,0,1110,12,0,11 if kkkkkkkkkkxxf xIxIxIxxxIxIx-，可認(rèn)為可認(rèn)為 滿足精度要求。滿足精度要求。 0,11kf xIx-， 0,11-10,12,-kkkkkf xIxx xf xIxx x-，根據(jù)前面的計(jì)算結(jié)果估計(jì)當(dāng)根據(jù)前面的計(jì)算

23、結(jié)果估計(jì)當(dāng)前的誤差：事后誤差估計(jì)前的誤差：事后誤差估計(jì)（實(shí)用），前面給出的誤差（實(shí)用），前面給出的誤差估計(jì)（事先誤差估計(jì)）不實(shí)估計(jì)（事先誤差估計(jì)）不實(shí)用用HW: p.58-59#1-84 牛頓插值牛頓插值 /* Newtons Interpolation */Lagrange 插值雖然易算，但若要增加一個(gè)節(jié)點(diǎn)時(shí)，插值雖然易算，但若要增加一個(gè)節(jié)點(diǎn)時(shí)，全部基函數(shù)全部基函數(shù) li(x) 都需重新算過。公式不具有繼承性，都需重新算過。公式不具有繼承性，不利于編程。不利于編程。將將 Ln(x) 改寫成改寫成.)()(102010 - - - - - xxxxaxxaa).(10- - - - nnxxx

24、xa的形式，希望每加一個(gè)節(jié)點(diǎn)時(shí)，的形式，希望每加一個(gè)節(jié)點(diǎn)時(shí)，只附加一項(xiàng)只附加一項(xiàng)上去即可。上去即可。? 差商差商( (亦稱均差亦稱均差) ) /* divided difference */),()()(,jijijijixxjixxxfxfxxf - - -= =1階差商階差商 /* the 1st divided difference of f w.r.t. xi and xj */)(,kixxxxfxxfxxxfkikjjikji - - -= =2階差商階差商f(x0)1010()()f xf xxx-2010201021()()()()f xf xf xf xxxxxxx-1階差商

25、的幾何階差商的幾何意義：弦截線意義：弦截線的斜率的斜率4 Newtons Interpolation4 Newtons Interpolation11101010111010,.,.,.,.,., - - - - - -= =- - -= =kkkkkkkkkkkxxxxxfxxxfxxxxxfxxxfxxf(k+1)階差商階差商： = = = =kiikikxxfxxf010)()(,., 事實(shí)上事實(shí)上其中其中,)()(01= = - -= =kiikxxx = = - -= = kijjjiikxxx01)()( Warning: my head is explodingWhat is t

26、he point of this formula?差商的值與差商的值與 xi 的順序無(wú)關(guān)！的順序無(wú)關(guān)！1.1.線性線性：1201020( )( )( ) , , ,kkkf xaf xbf xf xxaf xxbf xx=2.2.差商差商可以表示為可以表示為函數(shù)值的線性組合函數(shù)值的線性組合：0000-111( )( ),.,( )kkiikiiiiiiiikkif xf xf xxxxxxxxxxx=-,)()(01= = - -= =kiikxxx = = - -= = kijjjiikxxx01)()( 3.3. 對(duì)稱性：對(duì)稱性：由由2 2知，知，差商的值與節(jié)點(diǎn)的順序無(wú)關(guān)！差商的值與節(jié)點(diǎn)的

27、順序無(wú)關(guān)！4.4. 差商的另一種定義：由差商的另一種定義：由2，3及均差定義可得及均差定義可得10100 ,.,.,.,kkkkf xxf xxf xxxx-=-4 Newtons Interpolation6.5.5. 差商與導(dǎo)數(shù)的關(guān)系：差商與導(dǎo)數(shù)的關(guān)系：f(x)C ka,b,則則 a,b, s.t. 0,.,!kkff xxk=HW: 證明差商的性質(zhì)證明差商的性質(zhì)2，44 Newtons Interpolation4 Newtons Interpolation 牛頓插值牛頓插值 /* Newtons Interpolation */,)()()(000 xxfxxxfxf- - = =,)

28、(,101100 xxxfxxxxfxxf- - = =,.,)(,.,.,0010nnnnxxxfxxxxfxxxf- - = =- -).(.)()()(10102010- - - - - - - - - = =nnnxxxxaxxxxaxxaaxN .)(,)(,)()(102100100 - - - - - = =xxxxxxxfxxxxfxfxf).(,.,100- - - - nnxxxxxxf)().(,.,100nnnxxxxxxxxxf- - - - - -Nn(x)n次多次多項(xiàng)式，滿足：項(xiàng)式，滿足： Nn(xi)= f(xi)Rn(x)插值余項(xiàng)，插值余項(xiàng)，滿足滿足Rn(xi

29、)=0，i=0,n ai = f x0, , xi (1)(2)(n)(n)(n-1) (2) (1)4 Newtons Interpolation注：注： 由由唯一性可知唯一性可知 Nn(x) Ln(x)， 只是算法不同，表達(dá)只是算法不同，表達(dá)形式不同，故其余項(xiàng)也相同，即形式不同，故其余項(xiàng)也相同，即(1)011() ,.,( )( )(1)!nxnnnff x xxxxn=),(,!)(,.,maxmin)(0 xxkfxxfkk = = 實(shí)際計(jì)算過程為實(shí)際計(jì)算過程為f (x0)f (x1)f (x2)f (xn- -1)f (xn)f x0, x1f x1, x2 f xn- -1, xn

30、f x0, x1 , x2 f xn- -2, xn- -1, xnf x0, , xn f (xn+1) f xn, xn+1 f xn- -1, xn, xn+1 f x1, , xn+1 f x0, , xn+1增加增加如果精度不構(gòu)，增加節(jié)點(diǎn)如果精度不構(gòu)，增加節(jié)點(diǎn)xn+1, ,同時(shí)表中增加一行，三角形斜邊同時(shí)表中增加一行，三角形斜邊上即為所要求的各次項(xiàng)系數(shù)。上即為所要求的各次項(xiàng)系數(shù)。4 Newtons Interpolation 等距節(jié)點(diǎn)公式等距節(jié)點(diǎn)公式 /* Formulae with Equal Spacing */向前差分向前差分 /* forward difference */i

31、iifff- -= = 1ikikikikffff1111)(- - - - - - - = = = = 向后差分向后差分 /* backward difference */111- - - - - - = = ikikikfffi- -1iifff- -= = 中心差分中心差分 /* centered difference */212111- - - - - -= =ikikikfff 其中其中)(221hiixff = = 當(dāng)節(jié)點(diǎn)當(dāng)節(jié)點(diǎn)等距等距分布時(shí)分布時(shí):),.,0(0nihixxi= = = =fi= f(xi) 差分計(jì)算可通過構(gòu)造差分表得到差分計(jì)算可通過構(gòu)造差分表得到 2342345

32、0000000234111111232222233 = kkkkkkkx f xfffffxffffffxfffffxffffxf23344455 ffxffxf增加增加 23450011122222333 = kkkkkkkkxf xffffffxfxffxfffxff233323444444423455555555 ffxfffffxffffff1112211112211122identity calculus translation calculus: -, kkkkkkkkkkIffEffE ffE ffEffE II EEE-=-：， 差分的重要性質(zhì)：差分的重要性質(zhì)： 線性：例如線性

33、：例如gbfaxgbxfa = = )()( 各階差分可用函數(shù)值表示：各階差分可用函數(shù)值表示：!)1).(1(jjnnnjn - - -= = 其中其中/* binomial coefficients */4 Newtons Interpolation0011nnnjjnnjkkkn kjjjnnfEIfEffjj- -=-=-=-1()0011nnnnjnjnnjkkkkj njjnnfIEfEffjj- -=-=-=- 函數(shù)值可用各階差分表示：函數(shù)值可用各階差分表示：0e.g. Ennnjn kkkkjnffIffj= = 差商與差分的關(guān)系：差商與差分的關(guān)系：111121211222211

34、-1,-,-,1,-22, 1,!11,!kkkkkkkkkkkkkkkkkkkkkmkkkmkmmmkkkmkmkmmffff xxfxxxhf xxf xxfff xxxfxxhhgenerally we havef xxxfm hf xxxffm hm h-=- = 若若 f (x)是是 n 次多項(xiàng)式，則次多項(xiàng)式，則 是是 次多次多項(xiàng)式，而項(xiàng)式，而 ( ) (0)mf xmnnm-( )0 ()mf xmn= 差分與導(dǎo)數(shù)的關(guān)系差分與導(dǎo)數(shù)的關(guān)系( (由差分與差商、差商與導(dǎo)數(shù)的關(guān)系得）：由差分與差商、差商與導(dǎo)數(shù)的關(guān)系得）： ()1()1,!mmkkkmkmmmkmff xxxfmm hffh

35、=4 Newtons Interpolation等距節(jié)點(diǎn)牛頓公式等距節(jié)點(diǎn)牛頓公式 牛頓前差公式牛頓前差公式 /* Newtons forward-difference formula */001100010010-12000001( )( )( ),01, , ,( )( -)1( )(),( -),( -)( -)1111( -1)2!( ),nnjjnnnjjnnnnnf xNxRxxxthtxxjhxxtj hxx xht ttnNxf xf xxx xf xxxx xx xfftft tft ttnnRxf x xxx=-=-=-= -=01(1)( -)( -)1( )(1)!nnn

36、nx xx xt ttnhfn-=注：注：一般當(dāng)一般當(dāng) x 靠近靠近 x0 時(shí)用前插，靠近時(shí)用前插，靠近 xn 時(shí)用后插，故兩時(shí)用后插，故兩種公式亦稱為種公式亦稱為表初公式表初公式和和表末公式表末公式。 牛頓后差公式牛頓后差公式 /* Newtons forward-difference formula */nn110nnn 1nnn 10n12nnnnnnn 10, 10, , , ( )( -)1( )(),( -),( -)( -)1111( -1)2!( ),( -jjnnnjjnnxxthtxxjhxxtj hxx xht ttnNxf xf x xx xf x xxx xx xff

37、tft tft ttnnR xf x x xxx x=-=- =-= =節(jié)點(diǎn)倒排，n01(1)( -)1( )(1)!nnx xt ttnhfn=HW: p. 59#9-16 5 厄米插值厄米插值 /* Hermite Interpolation */不僅要求函數(shù)值重合，而且要求若干階不僅要求函數(shù)值重合，而且要求若干階導(dǎo)數(shù)導(dǎo)數(shù)也重合。也重合。即：要求插值函數(shù)即：要求插值函數(shù) (x) 滿足滿足 (xi) = f (xi), (xi) = f (xi), (mi) (xi) = f (mi) (xi).注：注： N 個(gè)條件可以確定個(gè)條件可以確定 階多項(xiàng)式。階多項(xiàng)式。N - - 1要求在要求在1個(gè)節(jié)

38、點(diǎn)個(gè)節(jié)點(diǎn) x0 處直到處直到m0 階導(dǎo)數(shù)都重合的插階導(dǎo)數(shù)都重合的插值多項(xiàng)式即為值多項(xiàng)式即為Taylor多項(xiàng)式多項(xiàng)式00)(!)(.)()()(000)(000mmxxmxfxxxfxfx- - - - = = 其余項(xiàng)為其余項(xiàng)為)1(00)1(00)()!1()()()()(-=-=mmxxmfxxfxR 一般只考慮一般只考慮 f 與與f 的值。的值。厄米插值厄米插值3 Hermite Interpolation例：例：設(shè)設(shè) x0 x1 x2, 已知已知 f(x0)、 f(x1)、 f(x2) 和和 f (x1), 求多項(xiàng)式求多項(xiàng)式 P(x) 滿足滿足 P(xi) = f (xi)，i = 0,

39、 1, 2，且且 P(x1) = f (x1), 并估計(jì)誤差。并估計(jì)誤差。模仿模仿 Lagrange 多項(xiàng)式的思想，設(shè)多項(xiàng)式的思想，設(shè)解：解：首先，首先，P 的階數(shù)的階數(shù) = 3=213)()()()()(=0iiixhx1f xhxfxP h0(x)有根有根x1, x2，且且 h0(x1) = 0 x1 是重根。是重根。)()()(22100 xxxxCxh- - -= =又又: h0(x0) = 1 C0 )()()()()(202102210 xxxxxxxxxh- - - - -= =h2(x)h1(x)有根有根 x0, x2 )()()(201xxxxBAxxh- - - = =由余

40、下條件由余下條件 h1(x1) = 1 和和 h1(x1) = 0 可解。可解。與與h0(x) 完全類似。完全類似。 (x) h1有根有根 x0, x1, x2 h1)()()(2101xxxxxxCx- - - -= = h1又又: (x1) = 1 C1 可解?？山?。其中其中 hi(xj) = ij , hi(x1) = 0, (xi) = 0, (x1) = 1 h1 h1),()()()()()(221033xxxxxxxKxPxfxR- - - -= =- -= =!4)()()4(xfxK = =與與 Lagrange 分析完全類分析完全類似似3 Hermite Interpola

41、tion一般地，已知一般地，已知 x0 , , xn 處有處有 y0 , , yn 和和 y0 , , yn ，求，求 H2n+1(x) 滿足滿足 H2n+1(xi) = yi ， H2n+1(xi) = yi。解：解：設(shè)設(shè)=ni)()()(=0iix x yixH2n+1n=0iyi其中其中 i(xj) = ij , i(xj) = 0, (xj) = 0, (xj) = ij i i(x)有根有根 x0 , , xi , , xn且都是且都是2重根重根 - - -= =ijjijixxxxxl)()()(由余下條件由余下條件 i (xi) = 1 和和 i(xi) = 0 可解可解Ai 和

42、和 Bi 2( )1 2 ( )()( )iiiiixl xxxlx=- (x) i有根有根 x0 , , xn, 除了除了xi 外都是外都是2重根重根 i)()(iili2(x)xxCx- -= =又又 i (xi) = 1 Ci = 1 i)(x)(ili2(x)xx- -= =設(shè)設(shè),.210baCfbxxxann = = = =則則20)22()()!22()()( - - = = = niixnnxxnfxR 這樣的這樣的Hermite 插值唯插值唯一一 i)()()(2xlBxAxiii = = i類似的，類似的，Th. Th. 設(shè)設(shè)f(x)C 2n+2a,b,則則 a,b, s.t

43、. 滿足下面插值條件滿足下面插值條件3 Hermite InterpolationQuiz: 給定給定 xi = i +1, i = 0, 1, 2, 3, 4, 5. 下面哪個(gè)是下面哪個(gè)是 i(x)的圖像？的圖像？x0-10.5123456yxy0-10.5123456斜率斜率=1 求求Hermite多項(xiàng)式的基本步驟：多項(xiàng)式的基本步驟： 寫出相應(yīng)于條件的寫出相應(yīng)于條件的 i(x)、 i(x) 的組合形式；的組合形式； 對(duì)每一個(gè)對(duì)每一個(gè) i(x)、 i(x)找出盡可能多的條件給出的根；找出盡可能多的條件給出的根； 根據(jù)多項(xiàng)式的總階數(shù)和根的個(gè)數(shù)寫出表達(dá)式；根據(jù)多項(xiàng)式的總階數(shù)和根的個(gè)數(shù)寫出表達(dá)式；

44、 根據(jù)尚未利用的條件解出表達(dá)式中的待定系數(shù)；根據(jù)尚未利用的條件解出表達(dá)式中的待定系數(shù)； 最后完整寫出最后完整寫出H2n+1(x)。HW: p.59#17-19注：待定系數(shù)法仍適用，但插值節(jié)點(diǎn)多注：待定系數(shù)法仍適用，但插值節(jié)點(diǎn)多時(shí)比較麻煩。時(shí)比較麻煩。4 分段低次插值分段低次插值 /* piecewise polynomial approximation */Remember what I have said? Increasing the degree of interpolating polynomial will NOT guarantee a good result, since hig

45、h-degree polynomials are oscillating.例：例：在在- -5, 5上考察上考察 的的Ln(x)。取。取211)(xxf=),., 0(105niinxi= = - -= = -5 -4 -3 -2 -1 0 1 2 3 4 5 -0.5 0 0.5 1 1.5 2 2.5 n 越大，越大，端點(diǎn)附近抖動(dòng)端點(diǎn)附近抖動(dòng)越大，稱為越大，稱為Runge 現(xiàn)象現(xiàn)象Ln(x) f (x) 分段分段低次低次插值插值4 Piecewise Polynomial Approximation 分段線性插值分段線性插值 /* piecewise linear interpolatio

46、n */在每個(gè)區(qū)間在每個(gè)區(qū)間 上，用上，用1階多項(xiàng)式階多項(xiàng)式 (直線直線) 逼近逼近 f (x):,1 iixx11111)()( - - - - - -= = iiiiiiiiyxxxxyxxxxxPxf,for 1 iixxx記記 ，易證：當(dāng)，易證：當(dāng) 時(shí)，時(shí)，|max1iixxh- -= = 0h)()(1xfxPh一致一致失去了原函數(shù)的光滑性。失去了原函數(shù)的光滑性。 分段分段Hermite插值插值 /* Hermite piecewise polynomials */給定給定nnnyyyyxx ,.,;,.,;,.,000在在 上利用兩點(diǎn)的上利用兩點(diǎn)的 y 及及 y 構(gòu)造構(gòu)造3次次He

47、rmite函數(shù)函數(shù),1 iixx導(dǎo)數(shù)一般不易得到。導(dǎo)數(shù)一般不易得到。How can we make a smooth interpolation without asking too much from f ?Headache 5 三次樣條三次樣條 /* Cubic Spline */定義定義設(shè)設(shè) 。三次樣條函數(shù)三次樣條函數(shù) , 且在每個(gè)且在每個(gè) 上為上為三次多項(xiàng)式三次多項(xiàng)式 /* cubic polynomial */。若它同。若它同時(shí)還滿足時(shí)還滿足 ，則稱為，則稱為 f 的的三次樣條插值函三次樣條插值函數(shù)數(shù) /* cubic spline interpolant */.bxxxan= =

48、= =.10,)(2baCxS ,1 iixx),., 0(),()(nixfxSii= = =注：注：三次樣條與分段三次樣條與分段 Hermite 插值的根本區(qū)別在于插值的根本區(qū)別在于S(x)自自身光滑身光滑，不需要知道，不需要知道 f 的導(dǎo)數(shù)值（除了在的導(dǎo)數(shù)值（除了在2個(gè)端點(diǎn)可能需個(gè)端點(diǎn)可能需要）；而要）；而Hermite插值依賴于插值依賴于f 在所有插值點(diǎn)的導(dǎo)數(shù)值。在所有插值點(diǎn)的導(dǎo)數(shù)值。f(x)H(x)S(x)5 Cubic Spline 構(gòu)造三次樣條插值函數(shù)的構(gòu)造三次樣條插值函數(shù)的三彎矩法三彎矩法 /* method of bending moment */在在 上，記上，記,1jjx

49、x- -,1- - -= =jjjxxh,for )()(1jjjxxxxSxS- - = =對(duì)每個(gè)對(duì)每個(gè)j, 此為此為3次多項(xiàng)式次多項(xiàng)式則則 Sj”(x) 為為 次多項(xiàng)式，需次多項(xiàng)式，需 個(gè)點(diǎn)的值確定之。個(gè)點(diǎn)的值確定之。12設(shè)設(shè) Sj”(xj- -1) = Mj- -1， Sj”(xj) = Mj 對(duì)應(yīng)力學(xué)中的對(duì)應(yīng)力學(xué)中的梁彎矩梁彎矩，故名，故名對(duì)于對(duì)于x xj- -1, xj 可可得到得到Sj”(x) =jjjjjjhxxMhxxM11- - - - - -積分積分2次，可得次，可得 Sj(x) 和和 Sj(x) ：jjjjjjjAhxxMhxxM - - - - - - - -2)(2)

50、(21121Sj(x) =jjjjjjjjBxAhxxMhxxM - - - - - -6)(6)(3131Sj(x) =利用已知利用已知Sj(xj- -1) = yj- -1 Sj(xj) = yj 可解可解5 Cubic SplinejjjjjjjhMMhyyA611- - - - - -= =jjjjjjjjjjjjhxxhMyhxxhMyBxA12211)6()6(- - - - - - - - -= = 下面解決下面解決 Mj ： 利用利用S 在在 xj 的的連續(xù)性連續(xù)性xj- -1, xj : Sj(x) =jjjjjjjjjjjhMMxxfhxxMhxxM6,2)(2)(1121

51、21- - - - - - - - - - - -1111211216,2)(2)( - - - - - - - -jjjjjjjjjjjhMMxxfhxxMhxxMxj , xj+1: Sj+1(x) =利用利用Sj(xj) = Sj+1(xj)，合并關(guān)于，合并關(guān)于Mj- -1、 Mj、 Mj+1的同類項(xiàng)，并的同類項(xiàng)，并記記 , , , 整理整理后得到：后得到：11jjjjhhh=l l1jj-=l lm m),(6111jjjjjjjxxfxxfhhg-=211gMMMjjjjjj=-l lm m j 1n- -1即：有即：有 個(gè)未知數(shù)，個(gè)未知數(shù)， 個(gè)方程。個(gè)方程。n- -1n+1 = =

52、 - - - -110111122nnnnggMMl lm ml lm m還需還需2個(gè)個(gè)邊界條件邊界條件 /* boundary conditions */5 Cubic Spline 第第1類邊條件類邊條件 /* clamped boundary */： S(a) = y0， S(b) = yna , x1 : S1(x) =1011012112106,2)(2)(hMMxxfhaxMhxxM- - - - - - - -010110),(62gy0 xxfhMM= =- -= = nnnnnngxxfynhMM= =- -= = - - -),(6211類似地利用類似地利用 xn- -1,

53、 b 上的上的 Sn(x) 第第2類邊條件：類邊條件： S”(a) = y0” = M0， S”(b) = yn” = Mn這時(shí)：這時(shí)：nnnygyg = = = = = =2,0;2,0000m ml l特別地，特別地，M0 = Mn = 0 稱為稱為自由邊界自由邊界 /* free boundary */，對(duì)應(yīng)的對(duì)應(yīng)的樣條函數(shù)稱為樣條函數(shù)稱為自然樣條自然樣條 /* Natural Spline */。 第第3類邊條件類邊條件 /* periodic boundary */ ： 當(dāng)當(dāng) f 為為周期函數(shù)周期函數(shù)時(shí)，時(shí)， yn = y0 ， S(a+) = S(b- -) M0 = Mn = =

54、 - - -nnnnnnggMM111122112222m ml ll lm ml lm mm ml l5 Cubic Spline注：注：另有另有三轉(zhuǎn)角法三轉(zhuǎn)角法(p.49-53)得到樣條函數(shù)，即設(shè)得到樣條函數(shù)，即設(shè) Sj(xj) = mj，則易知，則易知xj- -1, xj 上的上的Sj(x) 就是就是Hermite函數(shù)函數(shù)。再。再利用利用S”的連續(xù)性，可導(dǎo)出關(guān)于的連續(xù)性，可導(dǎo)出關(guān)于mj 的方程組，加上邊的方程組，加上邊界條件即可解。界條件即可解。Cubic Spline 由由boundary conditions 唯一唯一確定。確定。收斂性：收斂性：若若 ，且，且 ，則，則,baCf Chhiiminmax一致一致S(x) f(x)0maxihas即即:提高精度只須提高精度只須增加節(jié)點(diǎn)增加節(jié)點(diǎn), 而無(wú)須提高樣條階數(shù)。而無(wú)須提高樣條階數(shù)。穩(wěn)定性：穩(wěn)定性：只要邊條件保證只要邊條件保證 | m m0 |, | l l0 |, | m mn |, | l l n | 2，則方程組系數(shù)陣為則方程組系數(shù)陣為SDD陣陣，保證數(shù)值穩(wěn)定。保證數(shù)值穩(wěn)定。HW: p.59-60,#19-255 Cubic Spline Sketch of the Algorithm: Cubic Spline 計(jì)算計(jì)算 m mj , l l j , gj ; 計(jì)算計(jì)算 Mj (