1、第5章 數(shù)值積分與數(shù)值微分方法關(guān)于定積分計(jì)算，已經(jīng)有較多方法，如公式法、分步積分法等，但實(shí)際問題中，經(jīng)常出現(xiàn)不能用通常這些積分方法計(jì)算的定積分問題。怎樣把這些通常方法失效的定積分在一定精度下快速計(jì)算出來，特別是通過計(jì)算機(jī)編程計(jì)算出來就是本章研究的內(nèi)容。此外，怎樣根據(jù)函數(shù)在若干個(gè)點(diǎn)處的函數(shù)值去求該函數(shù)的導(dǎo)數(shù)近似值也是本章介紹的內(nèi)容。本章涉及的方法有Newton-Cotes求積公式、Gauss求積公式、復(fù)化求積公式、Romberg求積公式和數(shù)值微分。5.1 引 例人造地球衛(wèi)星軌道可視為平面上的橢圓。我國的第一顆人造地球衛(wèi)星近地點(diǎn)距離地球表面439km，遠(yuǎn)地點(diǎn)距地球表面2384km，地球半徑為637

2、1km，求該衛(wèi)星的軌道長度。本問題可用橢圓參數(shù)方程來描述人造地球衛(wèi)星的軌道，式中a, b分別為橢圓的長短軸，該軌道的長度L就是如下參數(shù)方程弧長積分但這個(gè)積分是橢圓積分，不能用解析方法計(jì)算。5.2問題的描述與基本概念要想用計(jì)算機(jī)來計(jì)算，應(yīng)對其做離散化處理。注意到定積分是如下和式的極限要離散化，做1） 去掉極限號(hào)2） 將取為具體的值3） 為減少離散化帶來的誤差，將用待定系數(shù)代替于是就得到定義5.1 若存在實(shí)數(shù)且任取都有 （5.1）則稱式(5.1)為一個(gè)數(shù)值求積公式。式中稱為求積系數(shù)，稱為求積節(jié)點(diǎn)，而稱 （5.2）為求積余項(xiàng)或求積公式（5.1）的截?cái)嗾`差。從定義可以看到，數(shù)值求積公式依賴于求積節(jié)點(diǎn)個(gè)

3、數(shù)n、求積節(jié)點(diǎn)和求積系數(shù)，這三個(gè)量有一個(gè)發(fā)生變化，則產(chǎn)生不同的求積公式。定義5.2 若求積公式對所有不超過m次的多項(xiàng)式有求積余項(xiàng)，而對某一個(gè)m+1次多項(xiàng)式有，則稱該求積公式的代數(shù)精度為m。一般，一個(gè)求積公式的代數(shù)精度越大，則該求積公式越好。確定代數(shù)精度的方法依次取代入公式并驗(yàn)證是否成立。若第一個(gè)使不成立的k值為m，則對應(yīng)的代數(shù)精度為m-1。例 5.1確定求積公式的代數(shù)精度。解 取代入求積公式有易驗(yàn)證，但，故本題求積公式代數(shù)精度為3。例5.2確定下面求積公式的參數(shù)A，B，C，使它具有盡可能高的代數(shù)精度，并指出相應(yīng)的代數(shù)精度。解本題要先求出具體的求積公式，然后再判斷所求公式的代數(shù)精度。公式有3個(gè)待

4、定參數(shù)，h不是求積公式的參數(shù)，故利用3個(gè)條件得到的3個(gè)等式關(guān)系就可以解決求出具體求積公式的問題。依次取代入求積公式并取等號(hào)，有解之得故所求的求積公式為為確定其代數(shù)精度，再取代入求出的公式繼續(xù)計(jì)算，有，故所求的求積公式具有二階代數(shù)精度。5.3 插值型求積公式借助多項(xiàng)式插值函數(shù)來構(gòu)造的求積公式稱為插值型求積公式。一般選用不同的插值公式就可以得到不同的插值型求積公式。基本思想利用被積函數(shù)的插值函數(shù)代替做定積分的近似計(jì)算來構(gòu)造求積公式。1.構(gòu)造原理考慮在n個(gè)節(jié)點(diǎn)上的n-1次Lagrange插值多項(xiàng)式與的余項(xiàng)，有這里。兩邊取積分,有記 （5.3）則有 （5.4）若舍去，得求積公式求積系數(shù)的求積公式就是插

5、值型求積公式。插值型求積公式的求積余項(xiàng)當(dāng)為次數(shù)小于n次的多項(xiàng)式時(shí)，有，對應(yīng)的。因此插值型求積公式的代數(shù)精度至少為n-1。若取，代入式(5.4)，可得插值型求積公式的求積系數(shù)之和為下面具體介紹常用的幾個(gè)插值型求積公式。2.Newton-Cotes求積公式1)n點(diǎn)的Newton-Cotes公式的構(gòu)造將求積節(jié)點(diǎn)取為a,b上的等距節(jié)點(diǎn)做積分變量變換:則當(dāng)時(shí),有,于是有插值型求積公式的求積系數(shù)為記，則有常稱為Cotes系數(shù)，易驗(yàn)證通常稱 （5.6）為n點(diǎn)的Newton-Cotes公式。由于求積節(jié)點(diǎn)是等距的，因此也稱式（5.6）為等距節(jié)點(diǎn)求積公式。利用可以得出下面常用的Newton-Cotes公式A) 2

6、 點(diǎn)的Newton-Cotes公式（5.7）這正是我們熟悉的梯形公式。B) 3點(diǎn)的Newton-Cotes公式為 （5.8）稱它為Simpson公式或拋物線公式。表5.1 部分Cotes系數(shù)n23456789例5.3 試分別用梯形公式和Simpson公式計(jì)算解 用梯形公式計(jì)算，有用Simpson公式計(jì)算，有2）n點(diǎn)Newton-Cotes公式的代數(shù)精度定理5.1 當(dāng)求積節(jié)點(diǎn)個(gè)數(shù)n為奇數(shù)時(shí)，對應(yīng)的Newton-Cotes求積公式的代數(shù)精度至少為n。證明由于是插值型求積公式，故有對有記 ，易知，故是奇函數(shù)，得，得證。3） 梯形公式與Simpson公式的余項(xiàng)引理 5.1(積分中值定理)設(shè)，且在上不變

7、號(hào)，則有梯形公式余項(xiàng)為在a,b不變號(hào)，有2階連續(xù)導(dǎo)數(shù)，由引理 5.1，有梯形公式余項(xiàng) （5.9）拋物線公式的余項(xiàng) （5.10）4）Newton-Cotes公式的數(shù)值穩(wěn)定性設(shè)計(jì)算函數(shù)時(shí)產(chǎn)生舍入誤差為實(shí)際在計(jì)算機(jī)中參加計(jì)算的是的近似值故Newton-Cotes公式在計(jì)算機(jī)中產(chǎn)生的誤差為若記，則有由Cotes表5.1,當(dāng)時(shí),，從而有說明此時(shí)計(jì)算舍入誤差可以控制，從而Newton-Cotes公式是數(shù)值穩(wěn)定的。但當(dāng)n>8時(shí)，有正有負(fù)， 隨n增大而增大，從而導(dǎo)致舍入誤差增加。故n>8時(shí)，Newton-Cotes公式是數(shù)值不穩(wěn)定的。因而一般不用n>8的Newton-Cotes公式來做定積分

8、計(jì)算。3.Gauss求積公式1）Gauss求積公式的構(gòu)造與概念n點(diǎn)的插值型求積公式的代數(shù)精度至少是n-1 ,那么是否還能提高其代數(shù)精度呢？若能，其代數(shù)精度最大能是多少？為回答這個(gè)問題，觀察一下插值求積公式的構(gòu)造方法，發(fā)現(xiàn)其至少具有n-1次代數(shù)精度的結(jié)論是在限定求積節(jié)點(diǎn)的前提下得出的，若讓求積節(jié)點(diǎn)也可以自由取值，則就給提高代數(shù)精度創(chuàng)造了條件。為使問題討論更具一般性，這里考慮帶權(quán)的定積分求積公式 （5.11）式中是已知的非負(fù)函數(shù)，為區(qū)間a,b上的權(quán)函數(shù)，a,b可以取為。顯然在式（5.11）中，若，就是前面討論的求積公式。定理5.2 求積公式(5.11)的代數(shù)精度最大為。證明 設(shè)是求積公式(5.11

9、)的任意一組求積節(jié)點(diǎn)，用此節(jié)點(diǎn)構(gòu)造插值型求積公式(5.11)，并令取2n次多項(xiàng)式代入公式(5.11), 考慮它的求積余項(xiàng)，有因?yàn)槭谴味囗?xiàng)式，由代數(shù)精度定義得的代數(shù)精度不大于。為證求積公式的代數(shù)精度可以為，設(shè)是任意一個(gè)次多項(xiàng)式，用去除，由多項(xiàng)式除法有式中都是次數(shù)不大于的多項(xiàng)式，于是有(5.12)因?yàn)樯厦娴那蠓e公式(5.11)是具有n個(gè)節(jié)點(diǎn)的插值型求積公式，故其代數(shù)精度不小于，從而有要讓式（5.11）成為等式，由式（5.12）應(yīng)有 （5.13）式（5.13）要求對固定的n次多項(xiàng)式和任意次數(shù)不超過的多項(xiàng)式都成立，這樣可以用的這種任意性，選擇求積節(jié)點(diǎn)。由正交多項(xiàng)式理論可知：使式（5.13）成立的點(diǎn) 是

10、存在唯一的，且都在內(nèi)，它就是在關(guān)于權(quán)的次正交多項(xiàng)式的零點(diǎn)。于是選取這些點(diǎn)作為求積公式的求積節(jié)點(diǎn)，并構(gòu)造對應(yīng)的插值型求積公式，就得到具有次迭代精度的求積公式了。關(guān)于插值型求積公式的結(jié)論定理5.3點(diǎn)插值型求積公式的代數(shù)精度至少是，至多為。定義5.3 若點(diǎn)的求積公式具有次代數(shù)精度，則稱該求積公式為Gauss型求積公式，對應(yīng)的求積節(jié)點(diǎn)和求積系數(shù)分別稱為Gauss點(diǎn)和Gauss系數(shù)。例5.4確定參數(shù)，使求積公式成為Gauss求積公式。解 注意到被積函數(shù)中有一因式與求積公式右端無明顯的關(guān)系，可將其視為權(quán)函數(shù)。為確定四個(gè)參數(shù)，依次取代入公式并將近似號(hào)取為等號(hào)，得聯(lián)立方程組解出由于本題求積節(jié)點(diǎn)個(gè)數(shù)為，其最高代

11、數(shù)精度，故所求的求積公式代數(shù)精度至少是3，故它是Gauss公式。2) Gauss型求積公式的求積余項(xiàng)設(shè),取在Gauss點(diǎn)上的次Hermite插值多項(xiàng)式,由Hermite插值余項(xiàng)公式，有兩邊積分，有由Gauss型求積公式的代數(shù)精度為及積分中值定理有，Gauss求積余項(xiàng)為3) Gauss型求積公式的數(shù)值穩(wěn)定性對任意i, 取是關(guān)于Gauss點(diǎn)的次Lagrange插值函數(shù)，有,由Gauss公式的代數(shù)精度為，而是次多項(xiàng)式，有由i的任意性可知，Gauss求積系數(shù)。若在Gauss公式中取，可得類似Newton-Cotes穩(wěn)定性處理方法，有在計(jì)算機(jī)計(jì)算時(shí)的舍入誤差這說明Gauss型求積公式是穩(wěn)定的。在Gaus

12、s公式中，不同的權(quán)函數(shù)和不同積分區(qū)間，對應(yīng)不同形式的Gauss公式，但最基本和常用的Gauss公式有4種。4） 常用的Gauss型求積公式（1） Gauss-Legendre求積公式權(quán)函數(shù)，積分區(qū)間為, Gauss點(diǎn)為次Legendre正交多項(xiàng)式的零點(diǎn)，Gauss-Legendre求積公式為Gauss-Legendre求積余項(xiàng)為Gauss-Legendre求積公式的Gauss點(diǎn)與系數(shù)表 5.2nx kA knx kA k10.000 000 02.000 000 040.861 136 30.347 854 820.577 350 31.000 000 00.339 981 00.652 14

13、5 230.000 000 00.888 888 950.000 000 00.568 888 90.774 596 70.555 555 60.906 179 80.236 926 90.538 469 30.478 628 7Gauss-Legendre求積公式可以計(jì)算任何有限積分區(qū)間的定積分，計(jì)算之前先作變量代換將積分區(qū)間變到，有然后再對用Gauss-Legendre求積公式。（2） Gauss-Chebyshev求積公式權(quán)函數(shù),求積區(qū)間為，Gauss點(diǎn)為次Chebyshev正交多項(xiàng)式的零點(diǎn)，Gauss-Chebyshev求積公式為式中Gauss點(diǎn)與系數(shù)Gauss-Chebyshev求積

14、余項(xiàng)（3） Gauss-Laguerre求積公式權(quán)函數(shù),積分區(qū)間為，Gauss點(diǎn)為次Laguerre正交多項(xiàng)式的零點(diǎn)，Gauss-Laguerre求積公式為Gauss-Laguerre求積余項(xiàng)為Gauss-Laguerre求積公式的Gauss點(diǎn)和系數(shù)也有事先計(jì)算好的表。（4）Gauss-Hermite求積公式權(quán)函數(shù),積分區(qū)間, Gauss點(diǎn)為n次Hermite正交多項(xiàng)式的零點(diǎn)，Gauss-Hermite求積公式為Gauss-Hermite求積余項(xiàng)為Gauss-Hermite求積公式的Gauss點(diǎn)及系數(shù)有表。例5.5用兩點(diǎn)Gauss公式求定積分的近似值。解 本題為有限區(qū)間的定積分，可用兩點(diǎn)Gau

15、ss-Legendre 求積公式計(jì)算。做積分換元，將其化為-1,1上的定積分，即令有本題準(zhǔn)確值為，可見精度很高。5.4復(fù)化求積公式Newton-Cotes公式在n>8時(shí)數(shù)值不穩(wěn)定，因此不能用增加求積節(jié)點(diǎn)的方法來提高計(jì)算精度。實(shí)用中常用復(fù)化求積公式來求積區(qū)間a,b上的定積分，以獲得滿足給定計(jì)算精度要的定積分值。常用的復(fù)化求積公式有復(fù)化梯形公式和復(fù)化Simpson公式。 基本思想將求積區(qū)間a,b分成若干個(gè)小區(qū)間，然后在每個(gè)小區(qū)間上采用數(shù)值穩(wěn)定的Newton-Cotes公式求小區(qū)間上的定積分，最后把所有小區(qū)間上的計(jì)算結(jié)果相加起來作為原定積分的近似值。1. 復(fù)化梯形公式1）復(fù)化梯形公式的構(gòu)造原理

16、取等距節(jié)點(diǎn)將積分區(qū)間a,bn等分,在每個(gè)小區(qū)間上用梯形公式做近似計(jì)算，就有得求積公式-復(fù)化梯形公式： (5.14)2）復(fù)化梯形公式的余項(xiàng)記(5.15)故復(fù)化梯形公式的求積余項(xiàng)由此可知,復(fù)化梯形公式的代數(shù)精度是1。若，對給定計(jì)算精度，令得出說明利用復(fù)化求積公式能得到計(jì)算誤差小于的定積分近似值。2. 復(fù)化Simpson公式1） 復(fù)化Simpson公式的構(gòu)造原理取等距節(jié)點(diǎn)將積分區(qū)間a,b n等分,在每個(gè)小區(qū)間上用Simpson公式做近似計(jì)算，再累加起來就有式中，得復(fù)化Simpson公式 (5.16)2）復(fù)化Simpson公式的余項(xiàng)記有復(fù)化Simpson公式的求積余項(xiàng)從復(fù)化simpson公式的求積余項(xiàng)

17、可以看出復(fù)化simpson公式的代數(shù)精度是3，它在代數(shù)精度和計(jì)算精度上都比復(fù)化梯形公式好。復(fù)化Simpson公式也稱為復(fù)化拋物線公式。例5.6分別用復(fù)化梯形公式和復(fù)化Simpson公式計(jì)算，要求誤差不超過。解 為較快得結(jié)果，積分區(qū)間分割數(shù)按進(jìn)行。數(shù)值計(jì)算結(jié)果列表，其中代表求積余項(xiàng)。N復(fù)化梯形公式復(fù)化Simpson公式2-17.389 2595.32-11.592 840-0.47822-13.336 0231.27-11.984 94423-12.382 1620.312-12.064 20924-12.148 004-12.069 95125-12.089 74226-12.075 1942

18、7-12.071 55828-12.070 649本題積分的準(zhǔn)確值為，可見復(fù)化梯形公式和復(fù)化Simpson公式能求出精度較高的解。例5.7考慮用復(fù)化Simpson公式計(jì)算要使誤差小于，那么求積區(qū)間0,1應(yīng)分成多少個(gè)子區(qū)間。以此計(jì)算積分近似值。解 復(fù)化Simpson公式的求積余項(xiàng)為式中。為估計(jì)誤差，要計(jì)算。注意到，故由此得從而有解出，故要求出滿足計(jì)算精度要求的定積分值，只要將0,1分成4個(gè)子區(qū)間即可，此時(shí)可算出5.5 Romberg 求積方法Romberg 求積方法是對復(fù)化梯形公式用加速技術(shù)得到的一種求積方法，它也稱為逐次分半加速收斂法。基本思想將Richardson 外推算法應(yīng)用于復(fù)化梯形公式

19、中，用產(chǎn)生的加速數(shù)列來求定積分值。1.構(gòu)造原理定理 5.4 設(shè)函數(shù)逼近數(shù)的余項(xiàng)為式中都是與h無關(guān)的常數(shù)，且時(shí)，,則由 ，q為任意常數(shù)定義的函數(shù)也逼近，且有式中都是與h無關(guān)的非零常數(shù)。證明 用qh 代替余項(xiàng)式（5.18）的變元h，有(5.19)用乘式（5.18）并與式(5.19)相加，整理后，有(5.20)因?yàn)?<q<1，故，用同除式(5.20)，有令即得所證。由微積分收斂階概念可知的收斂階比的收斂階大，故比逼近的速度快。通常稱如上的方法為Richardson 外推法。顯然這種外推可以不斷做下去以獲得逼近更快的函數(shù)，一般有這樣，每用一次Richardson外推，就使逼近階提高一個(gè)等級(jí)

20、，從而達(dá)到加速的目的。分析記稱為復(fù)化梯形公式的T形值，可以證明：此式說明逼近的階為。利用Richardson 外推法對做加速。因?yàn)?有逼近的階變?yōu)椤Ｈ羧。瑒t有同理對再做一次Richardson 加速，有逼近I的階為。一般，有經(jīng)Richardson加速的求定積分的序列為 (5.21)逼近I的階變?yōu)椤Ｗ⒁獾剑苯佑?jì)算可知就是復(fù)化Simpson公式的，即有 (5.22)類似地可得 (5.23) (5.24)這樣可以用復(fù)化梯形公式計(jì)算出序列。再逐次用公式（5.22），（5.23），（5.24）對其進(jìn)行加工得到收斂更快的序列，等等。通常將序列，和依次稱為梯形序列、Simpson序列、Cotes序列和Ro

21、mberg序列。若繼續(xù)外推，舍入誤差增大將使新序列體現(xiàn)不出更明顯的加速效果，因此一般只外推到Romberg序列。用Romberg序列求定積分的算法稱為Romberg求積方法。3Romberg求積方法的計(jì)算過程k01234因?yàn)榕c有相重的點(diǎn)，計(jì)算時(shí)還可以利用下面關(guān)系式來減小計(jì)算量。在給定計(jì)算誤差界后，可用來作為終止Romberg計(jì)算的條件。例5.8 用Romberg算法計(jì)算計(jì)算結(jié)果精確到。解 由Romberg求積方法，得計(jì)算結(jié)果如下表k00.920 73550.946 14590.946 08300.946 083110.939 79330.946 08690.946 08310.946 0830

22、20.944 51350.946 08340.946 083030.945 69090.946 083040.9845 9850由 ，故有。5.6 數(shù)值微分根據(jù)函數(shù)在若干個(gè)點(diǎn)處的函數(shù)值去求該函數(shù)的導(dǎo)數(shù)近似值稱為數(shù)值微分，所求導(dǎo)數(shù)的近似值常稱為數(shù)值導(dǎo)數(shù)。設(shè)是已知的個(gè)點(diǎn)處的函數(shù)值，下面分別討論兩種求數(shù)值微分的方法。1.利用次多項(xiàng)式插值函數(shù)求數(shù)值導(dǎo)數(shù)用插值函數(shù)代替被插值函數(shù)可以用來計(jì)算被插函數(shù)的近似值；用插值函數(shù)代替定積分的被積函數(shù)可以用來計(jì)算定積分的近似值；用插值函數(shù)的導(dǎo)數(shù)能否用來求被插函數(shù)導(dǎo)數(shù)值？定義5.4設(shè)是的次插值多項(xiàng)式，稱用插值函數(shù)的導(dǎo)數(shù)來計(jì)算被插函數(shù)導(dǎo)數(shù)的公式為的插值型求導(dǎo)公式。由插值理

23、論，有特別有一階數(shù)值導(dǎo)數(shù)的余項(xiàng)關(guān)系在非節(jié)點(diǎn)處的余項(xiàng)由于涉及對未知中值函數(shù)的求導(dǎo)，使余項(xiàng)不易描述和控制，但在節(jié)點(diǎn)處有此余項(xiàng)易得出誤差的界限，因此插值型求導(dǎo)公式多用于求在節(jié)點(diǎn)處的數(shù)值導(dǎo)數(shù)。下面給出幾個(gè)常用的數(shù)值微分公式。1) 兩點(diǎn)公式給定點(diǎn)，有線性插值多項(xiàng)式記，有不帶余項(xiàng)的數(shù)值微分公式在節(jié)點(diǎn)處的余項(xiàng)帶余項(xiàng)的兩點(diǎn)數(shù)值微分公式前差公式 后差公式 2) 三點(diǎn)公式只給出等距節(jié)點(diǎn)的公式形式，對非等距節(jié)點(diǎn)有類似的公式。給定點(diǎn)集，有二次插值多項(xiàng)式，類似前面的做法，可得及帶余項(xiàng)的三點(diǎn)數(shù)值微分公式前差公式中心公式后差公式中間的節(jié)點(diǎn)使用中心公式，左邊界節(jié)點(diǎn)用前插公式，右邊界點(diǎn)用后插公式。在求數(shù)值微分公式余項(xiàng)表達(dá)式時(shí)，有時(shí)從插值函數(shù)的余項(xiàng)表達(dá)式不易得出數(shù)值微分公式的余項(xiàng)，但借助Taylor中值定