1、數值分析(研究生)第三章數值積分與數值微分二第三章第三章 數值積分與數值微分（二）數值積分與數值微分（二）第五節第五節 Romberg求積算法求積算法第六節第六節 Gauss求積公式求積公式第七節第七節 數值微分數值微分數值分析(研究生)第三章數值積分與數值微分二一、梯形公式的遞推公式及事后估計法一、梯形公式的遞推公式及事后估計法 上一節介紹的復化求積方法對提高精度是行之有效的，但上一節介紹的復化求積方法對提高精度是行之有效的，但在使用求積公式之前必須給出合適的步長，步長取得太大精度在使用求積公式之前必須給出合適的步長，步長取得太大精度難以保證，步長太小則會導致計算量的增加，而事先給出一個難以

2、保證，步長太小則會導致計算量的增加，而事先給出一個恰當的步長又往往是困難的恰當的步長又往往是困難的 實際計算中常常采用變步長的計算方案，即在步長逐次分實際計算中常常采用變步長的計算方案，即在步長逐次分半半(即步長二分即步長二分)的過程中，反復利用復化求積公式進行計算，的過程中，反復利用復化求積公式進行計算，直至所求得的積分值滿足精度要求為止直至所求得的積分值滿足精度要求為止 設將求積區間設將求積區間a，b分成分成n等分，則一共有等分，則一共有n+1個分點，按個分點，按梯形公式計算積分值梯形公式計算積分值Tn，需要提供，需要提供n+1個函數值如果將求積個函數值如果將求積區間再二分一次，則分點增至

3、區間再二分一次，則分點增至2n+1個，我們將二分前后兩個積個，我們將二分前后兩個積分值聯系起來加以考察分值聯系起來加以考察5 龍貝格求積公式龍貝格求積公式數值分析(研究生)第三章數值積分與數值微分二 注意到每個子區間注意到每個子區間xk，xk+1經過二分只增加了一個分點經過二分只增加了一個分點xk+1/2( xk+xk+1)/2，用復化梯形公式求得該子區間上的積分，用復化梯形公式求得該子區間上的積分值為值為 101021102110122)12(221)(221)(2)()(4nknnkknnkknkkknhkafhTxfhTxfhxfxfhT數值分析(研究生)第三章數值積分與數值微分二.si

4、n110dxxxI 計算積分值計算積分值例例.9207355. 0)1()0(211 , 0sin)(1 ffTxxxf上上應應用用梯梯形形公公式式，有有先先在在對對解解.9397933. 0)21(2121,9588510. 0)21(12 fTTf代代入入遞遞推推公公式式，有有將將區區間間二二等等分分，求求出出.9445135. 0)43()41(4121,9088516. 0)43(,9896158. 0)41(24 ffTTff代代入入遞遞推推公公式式，有有求求出出新新分分點點上上的的函函數數值值進進一一步步二二分分求求積積區區間間，數值分析(研究生)第三章數值積分與數值微分二).2k

5、nk 區間等分數區間等分數代表二分次數，代表二分次數，算結果見下表（表中算結果見下表（表中這樣不斷二分下去，計這樣不斷二分下去，計9460831. 09460830. 09460827. 09460815. 09460769. 01098769460596. 09459850. 09456909. 09445135. 09397933. 054321nnTkTk.1025.107計計算算量量的的問問題題何何提提高高收收斂斂速速度度、節節省省因因此此，接接下下來來要要研研究究如如斂斂速速度度太太慢慢量量很很大大，收收個個分分點點的的函函數數值值，計計算算即即要要提提供供次次需需要要二二分分區區間

6、間位位有有效效數數字字的的精精度度要要達達到到式式計計算算積積分分它它表表明明，用用復復化化梯梯形形公公I數值分析(研究生)第三章數值積分與數值微分二二、龍貝格算法二、龍貝格算法).,()(212);,()(12)(222bafhabTIbafhabTIfRnnn 有有：根據復化梯形公式的余項表達式根據復化梯形公式的余項表達式. )(31.41)()(222nnnnnTTTITITIff 整整理理后后可可得得：，則則有有假假定定 可見，利用兩種步長計算的結果能估計截斷誤差可見，利用兩種步長計算的結果能估計截斷誤差.這種利用計這種利用計算結果估計誤差的方法稱為算結果估計誤差的方法稱為事后估計法事

7、后估計法. 若將該截斷誤差加到計算結果中，就可得出若將該截斷誤差加到計算結果中，就可得出“改進梯形改進梯形求積公式求積公式”：nnnnnTTTTTT3134)(31222 數值分析(研究生)第三章數值積分與數值微分二改進梯形求積公式的右邊實際是改進梯形求積公式的右邊實際是nnknkkknkknkknkknnnkknnnSbfxfxfafhxfhbfxfafhxfhTTxfhTTT 101121102111102110212)()(2)(4)(6)(2)()(2)(231)(231)(221431)4(31這就是說用梯形法二分前后的兩個積分值這就是說用梯形法二分前后的兩個積分值Tn與與T2n的線

8、性組的線性組合的結果得到復化辛普森法求積公式合的結果得到復化辛普森法求積公式)1 . 5(141144313422nnnnnTTTTS 數值分析(研究生)第三章數值積分與數值微分二 類似的情況，用辛普森法二分前后的兩個積分值類似的情況，用辛普森法二分前后的兩個積分值Sn與與S2n的的線性組合的結果可得到復化柯特斯求積公式線性組合的結果可得到復化柯特斯求積公式)2 . 5(151151614114422222nnnnnSSSSC 重復同樣的手續，用柯特斯法二分前后的兩個積分值重復同樣的手續，用柯特斯法二分前后的兩個積分值Cn與與C2n的線性組合的結果可得到龍貝格的線性組合的結果可得到龍貝格(Ro

9、mberg)求積公式求積公式)3 . 5(631636414114423233nnnnnCCCCR 我們在變步長的過程中運用加速公式我們在變步長的過程中運用加速公式(5.1)、(5.2)、(5.3)，就能將粗糙的梯形值就能將粗糙的梯形值Tn逐步加工成精度較高的辛普森值逐步加工成精度較高的辛普森值Sn 、柯、柯特斯值特斯值Cn和龍貝格值和龍貝格值Rn .龍貝格求積算法可用下表來表示：龍貝格求積算法可用下表來表示：數值分析(研究生)第三章數值積分與數值微分二數值分析(研究生)第三章數值積分與數值微分二 例例2 用龍貝格方法計算橢圓用龍貝格方法計算橢圓 x2/4 + y2 l 的周長，使結果具有的周

10、長，使結果具有五位有效數字五位有效數字 分析分析 為便于計算，先將橢圓方程采用參數形式表示，再根據弧為便于計算，先將橢圓方程采用參數形式表示，再根據弧長公式將橢圓周長用積分形式表示由于計算結果要求具有五位有長公式將橢圓周長用積分形式表示由于計算結果要求具有五位有效數字，因此需要估計所求積分值有幾位整數，從而確定所求積分效數字，因此需要估計所求積分值有幾位整數，從而確定所求積分值的絕對誤差限最后再應用龍貝格方法計算積分值的絕對誤差限最后再應用龍貝格方法計算積分 解解 令令 x 2cosq q，y sinq q 則橢圓的周長為則橢圓的周長為Iyxl4d sin314d42022022 q qq q

11、q qq qq q.10125. 01081)(1021)(4422d sin3124451202 fRIfRIlI斷斷誤誤差差為為的的截截，故故計計算算的的截截斷斷誤誤差差為為五五位位有有效效數數字字，則則需需有有有有一一位位整整數數，要要求求結結果果，因因此此由由于于 q qq q 數值分析(研究生)第三章數值積分與數值微分二kk2kT212 kS22 kC32 kR4322 kkRR數值分析(研究生)第三章數值積分與數值微分二三、理查森三、理查森(Richardson)外推加速法外推加速法 上面討論說明由梯形公式出發，將區間上面討論說明由梯形公式出發，將區間a，b逐次二分可逐次二分可提高

12、求積公式的精度，上述加速過程還可繼續下去，其理論依提高求積公式的精度，上述加速過程還可繼續下去，其理論依據是梯形公式的余項展開，即據是梯形公式的余項展開，即 llllhhhIhThhhIhT242212422121642,)( .)(122nabhbafhabTIn ， 22hTTn若記若記Tn = T(h)，當區間，當區間a，b分為分為2n等分時，有等分時，有 ，則，則可見可見I = T(h)的誤差為的誤差為O(h2)階階.若記若記 ，則，則 3)(24)(1hThThT 數值分析(研究生)第三章數值積分與數值微分二 64162)(6241162411hhIhThhIhT 顯然顯然T1(h)

13、與積分值與積分值 I 近似的階為近似的階為O(h4) .這樣構造的這樣構造的， 2)(11hThT就是辛普森公式序列就是辛普森公式序列Sn，S2n， . .若令若令 ，)(15121516)(112hThThT 則又可進一步從余項中消去則又可進一步從余項中消去 h4 項，這樣構造出的項，這樣構造出的 ，其實，其實就是柯特斯公式序列，它與積分值就是柯特斯公式序列，它與積分值 I 的逼近階為的逼近階為O(h6) . 如此繼如此繼續下去，每加速一次，誤差的量級便提高續下去，每加速一次，誤差的量級便提高2階，一般地，若記階，一般地，若記T0(h) = T(h)，經過，經過m (m = 1，2，)次加速

14、后，則有次加速后，則有 )(2hT)(1412144)(11hThThTmmmmmm 數值分析(研究生)第三章數值積分與數值微分二 )21(141144)(1)1(1)()(0)()(0，次次加加速速值值，可可得得的的表表示示序序列列以以次次后后求求得得的的梯梯形形值值，且且表表示示二二分分設設以以 kTTTmTTkTkmmkmmmkmkkmk.)3 . 5()2 . 5()1 . 5(321數數表表出出三三角角形形數數表表根根據據公公式式可可以以逐逐行行構構造造、公公式式，即即可可得得到到加加速速、算算法法，若若取取上上式式也也稱稱為為龍龍貝貝格格求求積積Tm 可以證明，如果可以證明，如果

15、f (x) 充分光滑，那么充分光滑，那么T 數表每一列的元數表每一列的元素及對角線元素均收斂到所求的積分值素及對角線元素均收斂到所求的積分值 I ，即，即ITmITkmmkmk )()(lim)(lim，固固定定數值分析(研究生)第三章數值積分與數值微分二.,)()(34200的的外外推推公公式式計計算算立立近近似似試試用用理理查查森森外外推推原原理理建建無無關關，為為常常數數與與步步長長、其其中中的的公公式式：是是計計算算函函數數值值設設例例JhbabhahhuJJhu 數值分析(研究生)第三章數值積分與數值微分二 機械求積公式機械求積公式 含有含有2n+2個待定參數個待定參數xk、Ak(k

16、0，1，n)當當 xk 為等距節點時得到的插值求積公式的代數精為等距節點時得到的插值求積公式的代數精度至少為度至少為n次，如果適當選取次，如果適當選取 xk (k0，1，n)，有可能使求積公式，有可能使求積公式具有具有 2n+1 次代數精度，這類求積公式稱為次代數精度，這類求積公式稱為高斯高斯(Gauss)求積公式求積公式. 為使問題更具一般性，我們研究帶權積分為使問題更具一般性，我們研究帶權積分 ，這，這里里r r (x)為權函數，類似機械求積公式，它的求積公式為為權函數，類似機械求積公式，它的求積公式為 Ak(k0，1，n)為不依賴于為不依賴于f (x)的求積系數，的求積系數，xk (k0

17、，1，n)為為求積節點，可適當選取求積節點，可適當選取 xk 及及 Ak (k0，1，n)使求積公式使求積公式(6.1)具有具有2n+1次代數精度次代數精度 bankkkxfAxxf0)(d)( baxxxfId)()(r r)1 . 6()(d)()(0 nkkkbaxfAxxxfr r6 6 高斯求積公式高斯求積公式數值分析(研究生)第三章數值積分與數值微分二一、高斯點一、高斯點 定義定義4 如果求積公式如果求積公式(6.1)具有具有2n+1次代數精度，則稱其節點次代數精度，則稱其節點 xk (k0，1，n)為高斯點，相應公式為高斯點，相應公式(6.1)稱為稱為高斯求積公式高斯求積公式 .

18、 根據定義要使根據定義要使(6.1)具有具有2n+1次代數精度，只要取次代數精度，只要取f(x)xm，對對m0，1，2n+1，(6.1)精確成立，則得精確成立，則得 當給定權函數當給定權函數r r (x)，求出右端積分，則可由，求出右端積分，則可由(6.2)解得解得 Ak 及及 xk (k0，1，n) bamnkmkknmxxxxA)2.6(.12,1,0d)(0r r 求解非線性方程組求解非線性方程組(6.2)較復雜，通常較復雜，通常n2就很難求解故一就很難求解故一般不通過解方程般不通過解方程(6.2)求求 xk 及及 Ak (k0，1，n)，而從分析高，而從分析高斯點的特性來構造高斯求積公

19、式斯點的特性來構造高斯求積公式數值分析(研究生)第三章數值積分與數值微分二 定理定理5 插值型求積公式插值型求積公式(6.1)的節點的節點 ax0 xlxnb是是高斯點的充分必要條件是以這些節點為零點的多項式高斯點的充分必要條件是以這些節點為零點的多項式 與任何次數不超過與任何次數不超過n的多項式的多項式P(x)帶權帶權r r (x) 正交，即正交，即)()()(101nnxxxxxxx )3 . 6(0)(d)()()(1 banxxxxxPr r 定理表明在定理表明在a，b上帶權上帶權r r (x)的的n+1次正交多項式的零點就是次正交多項式的零點就是求積公式求積公式(6.1)的高斯點，有

20、了求積節點的高斯點，有了求積節點 xk (k0，l，n)，再，再利用利用(6.2)對對m0，l，n 成立，則得到一組關于求積系數成立，則得到一組關于求積系數A0，A1，An 的線性方程解此方程則得的線性方程解此方程則得Ak(k0，1，n) . 也也可直接由可直接由x0，x1，xn 的插值多項式求出求積系數的插值多項式求出求積系數Ak(k = 0, 1, ,n) .數值分析(研究生)第三章數值積分與數值微分二二、高斯求積公式的余項二、高斯求積公式的余項 利用利用 f (x)在節點在節點xk(k = 0，1，n)的的埃爾米特插值埃爾米特插值 H2n+1 (x)，即，即于是于是 ，兩端乘，兩端乘r

21、r (x)，并由，并由a到到b積分，積分，則得則得其中右端第一項積分對其中右端第一項積分對2n+1次多項式精確成立，故次多項式精確成立，故 由于由于 0，故由積分中值定理得，故由積分中值定理得(6.1)的余項為的余項為.,1,0, )()(, )()(1212nkxfxHxfxHkknkkn )()!22()()()(21)22(12xnfxHxfnnn bannnkkknxxxnfxfAIfR.d)()()!22()()(21)22(0r r )()(21xxnr r )4.6(.d)()()!22()(21)22( bannnxxxnffRr r bannbafRxxxHxxxfId)()

22、(d)()(12r rr r數值分析(研究生)第三章數值積分與數值微分二 與積分相反，數值微分非常困難與積分相反，數值微分非常困難. 積分描述了一個函數的積分描述了一個函數的整體或宏觀性質，而微分則描述一個函數在一點處的斜率，這整體或宏觀性質，而微分則描述一個函數在一點處的斜率，這是函數的微觀性質是函數的微觀性質. 因此積分對函數的形狀在小范圍內的改變因此積分對函數的形狀在小范圍內的改變不敏感。而微分卻很敏感不敏感。而微分卻很敏感. 一個函數小的變化，容易產生相鄰一個函數小的變化，容易產生相鄰點的斜率的大的改變點的斜率的大的改變. 由于微分這個固有的困難，所以應盡可能避免數值微分，由于微分這個

23、固有的困難，所以應盡可能避免數值微分，特別是對實驗獲得的數據進行微分特別是對實驗獲得的數據進行微分. 在這種情況下，最好用最在這種情況下，最好用最小二乘曲線擬合這種數據，然后對所得到的多項式進行微分小二乘曲線擬合這種數據，然后對所得到的多項式進行微分. 或用另一種方法，對該數據進行三次樣條擬合，然后尋找該樣或用另一種方法，對該數據進行三次樣條擬合，然后尋找該樣條函數的微分條函數的微分. 7 7 數值微分數值微分數值分析(研究生)第三章數值積分與數值微分二一、中點方法與誤差分析一、中點方法與誤差分析 數值微分就是用函數值的線性組合近似函數在某點的導數數值微分就是用函數值的線性組合近似函數在某點的

24、導數 值按導數定義可以簡單地用差商近似導數，這樣立即得到幾值按導數定義可以簡單地用差商近似導數，這樣立即得到幾種數值微分公式種數值微分公式hhafhafafhhafafafhafhafaf2)()()()()()()()()( 其中其中h為一增量稱為步長后一種數值微分方法稱為中點為一增量稱為步長后一種數值微分方法稱為中點方法、它是前兩種方法的算術平均但它的誤差階卻由方法、它是前兩種方法的算術平均但它的誤差階卻由O(h)提提高到高到O(h2) 上面所給出的三個公式是很實用的尤其是中點公上面所給出的三個公式是很實用的尤其是中點公式更為常用式更為常用數值分析(研究生)第三章數值積分與數值微分二為要利

25、用中點公式為要利用中點公式hhafhafhG2)()()( )(af 計算導數計算導數的近似值，首先須選取合適的步長為此需要進行誤差分析的近似值，首先須選取合適的步長為此需要進行誤差分析分別將分別將在在 x=a 處做泰勒展開有處做泰勒展開有)(haf )(! 5)(! 4)(! 3)(! 2)()()()5(5)4(432afhafhafhafhafhafhaf代入代入G(h)得得 )(! 5)(! 3)()()5(42afhafhafhG 由此得知，從截斷誤差的角度看，步長越小，計算結果越由此得知，從截斷誤差的角度看，步長越小，計算結果越準確準確.且且MhhGaf6)()(2 其中其中.)(

26、maxxfMhax 數值分析(研究生)第三章數值積分與數值微分二 再考察舍入誤差按中點公式計算，當再考察舍入誤差按中點公式計算，當h很小時，因很小時，因 f (a+h)與與 f (a-h)很接近，直接相減會造成有效數字的嚴重損失很接近，直接相減會造成有效數字的嚴重損失(參看第參看第1章第章第4節節)因此，從舍入誤差的角度來看，步長不宜太小因此，從舍入誤差的角度來看，步長不宜太小則計算則計算 當當 f (a+h)及及 f (a-h)分別有舍入誤差分別有舍入誤差 1及及 2時，若令時，若令 21,max )(af 的舍入誤差上界為的舍入誤差上界為hhaGafaf 2)()()(21步長步長h不宜太

27、大，也不宜太小其最優步長應為不宜太大，也不宜太小其最優步長應為 它表明它表明h越小，舍入誤差越小，舍入誤差)(af 越大，故它是病態的用中點越大，故它是病態的用中點公式計算公式計算)(af 的誤差上界為的誤差上界為,6)(2hMhhE 要使誤差要使誤差E(h)最小，最小，3opt/3Mh 數值分析(研究生)第三章數值積分與數值微分二二、插值型的求導公式二、插值型的求導公式 x x0 x1 x2 xn y y0 y1 y2 yn 對于列表函數對于列表函數 y = f (x)： 運用插值原理，可以建立插值多項式運用插值原理，可以建立插值多項式 y = Pn(x)作為它的近作為它的近似由于多項式的求

28、導比較容易，我們取似由于多項式的求導比較容易，我們取統稱插值型的求導公式統稱插值型的求導公式)1. 7()()(xPxfn )(xPn )(xf 的近似值，這樣建立的數值公式的近似值，這樣建立的數值公式的值作為的值作為 依據插值余項定理，求導公式依據插值余項定理，求導公式(7.1)的余項為的余項為)(dd)!1()()()!1()()()()1(11)1( nnnnnfxnxxnfxPxf式中式中. )()(01 niinxxx 數值分析(研究生)第三章數值積分與數值微分二 如果我們限定求某個節點如果我們限定求某個節點 xk 上的導數值，那么上面的第二上的導數值，那么上面的第二項因式變為零，這

29、時有余項公式項因式變為零，這時有余項公式)2 . 7()()!1()()()(1)1(knnknkxnfxPxf 下面我們僅僅考察節點處的導數值為簡化討論，假定所給下面我們僅僅考察節點處的導數值為簡化討論，假定所給的節點是等距的的節點是等距的1兩點公式兩點公式 設已給出兩個節點設已給出兩個節點 x0， x1 上的函數值上的函數值 f (x0)，f (x1)，做線性，做線性插值得公式插值得公式)()()(101001011xfxxxxxfxxxxxP 對上式兩端求導，記對上式兩端求導，記x1 x0 = h，有，有)()(1)(101xfxfhxP 數值分析(研究生)第三章數值積分與數值微分二)(

30、)(1)()()(1)(01110101xfxfhxPxfxfhxP 于是有下列求導公式：于是有下列求導公式：)(2)()(1)()(2)()(1)(011010 fhxfxfhxffhxfxfhxf 而利用余項公式而利用余項公式(7.2)知，帶余項的兩點公式是：知，帶余項的兩點公式是：數值分析(研究生)第三章數值積分與數值微分二 2三點公式三點公式 設已給出三個節點設已給出三個節點x0，xl = x0+h，x2 = x0+2h上的函數值，上的函數值，做二次插值做二次插值)()()()()()()()()()(2120210121012002010212xfxxxxxxxxxfxxxxxxxx

31、xfxxxxxxxxxP 令令 x = x0+th，則，則 )3 . 7()()12()()44()()32(21)()()1(21)()2()()2)(1(21)(2100221002xftxftxfththxPxfttxfttxfttthxP 數值分析(研究生)第三章數值積分與數值微分二這里撇號這里撇號()表示對變量表示對變量x求導數上式分別取求導數上式分別取t = 0，1，2，得，得到三種三點公式：到三種三點公式： .)(3)(4)(21)(;)()(21)(;)()(4)(321)(21022201221002xfxfxfhxPxfxfhxPxfxfxfhxP 而帶余項的三點求導公式如

32、下：而帶余項的三點求導公式如下： . )(3)(3)(4)(21)()4 . 7();(6)()(21)(;)(3)()(4)(321)(22102220122100 fhxfxfxfhxffhxfxfhxffhxfxfxfhxf 數值分析(研究生)第三章數值積分與數值微分二公式公式(6.6)是我們所熟悉的中點公式在三點公式中，它由于少是我們所熟悉的中點公式在三點公式中，它由于少用了一個函數值用了一個函數值 f(x1) 而引人注目而引人注目 用插值多項式用插值多項式Pn(x)作為作為 f (x)的近似函數，還可以建立高階的近似函數，還可以建立高階數值微分公式：數值微分公式：,2,1, )()()()( kxPxfkmk例如，將式例如，將式(6.5)再對再對t求導一次，有求導一次，有 )()(2)(1)(210202xfxfxfhthxP 于是有于是有 )()(2)(1)(111212hxfxfhxfhxP 而帶余項的二階三點公式如下：而帶余項的二階三點公式如下： )5 . 7()(12)()(2)(1)()4(211121 fhhxfxfhxfhxf 數值分析(研究生)第三章數值積分與數值微分二例例4 4數值分析(研究生)第三章數值積分與數值微分二數值分析(研究