龍貝格積分實(shí)驗(yàn)報(bào)告

1、二、Romberg積分法1.變步長Romberg積分法的原理復(fù)化求積方法對于提高精度是行之有效的方法，但復(fù)化公式的一個(gè)主要缺點(diǎn)在于要事先估計(jì)出部長。若步長過大，則精度難于保證；若步長過小，則計(jì)算量又不會(huì)太大。而用復(fù)化公式的截?cái)嗾`差來估計(jì)步長，其結(jié)果是步長往往過小，而且和在區(qū)間上的上界的估計(jì)是較為困難的。在實(shí)際計(jì)算中通常采用變步長的方法，即把步長逐次分半（也就是把步長二等分），直到達(dá)到某種精度為止，這種方法就是Romberg積分法的思想。在步長的逐步分半過程中，要解決兩個(gè)問題：1. 在計(jì)算出后，如何計(jì)算,即導(dǎo)出和之間的遞推公式；2. 在計(jì)算出后，如何估計(jì)其誤差，即算法的終止的準(zhǔn)則是什么。首先推導(dǎo)

2、梯形值的遞推公式，在計(jì)算時(shí)，需要計(jì)算個(gè)點(diǎn)處的函數(shù)值在計(jì)算出后，在計(jì)算時(shí)，需將每個(gè)子區(qū)間再做二等分，共新增個(gè)節(jié)點(diǎn)。為了避免重復(fù)計(jì)算，計(jì)算時(shí)，將已計(jì)算的個(gè)點(diǎn)的數(shù)值保留下來，只計(jì)算新增個(gè)節(jié)點(diǎn)處的值。為此，把表示成兩部分之和，即由此得到梯形值遞推公式因此由復(fù)化梯形公式的截?cái)嗾`差有若變化不大時(shí)，即，則有式（2）表明，用作為定積分的近似值，其誤差大致為，因此其終止條件為其中是預(yù)先給定的精度。2.Romberg積分公式將上述方法不斷推廣下去，可以得到一個(gè)求積分的序列，而且這個(gè)序列很快收斂到所求的定積分。記，將區(qū)間等分的梯形值。，將區(qū)間等分的Simpson，將區(qū)間等分的Cotes。，將區(qū)間等分的Romberg

3、。由其可構(gòu)造一個(gè)序列，次序列稱為Romberg序列，并滿足如下遞推關(guān)系：以上遞推公式就是Romberg積分遞推公式。3.Romberg積分程序1. 置，精度要求，；2. 計(jì)算；3. 置，并計(jì)算；4. 置5. 計(jì)算；6. 若 ，則轉(zhuǎn)（7）；否則置，轉(zhuǎn)（5）；7. 若，則停止計(jì)算（輸出），否則轉(zhuǎn)（3）。4.Romberg積分法的應(yīng)用function T,n = romb(f,a,b,eps)double R;if nargin<4,eps=1e-8;endh=b-a;R(1,1)=(h/2)*(feval(f,a)+feval(f,b);n=1;J=0;err=1;while (err>

4、;eps) J=J+1;h=h/2;S=0; for i=1:n x=a+h*(2*i-1); S=S+feval(f,x); end R(J+1,1)=R(J,1)/2+h*S; for k=1:J R(J+1,k+1)=(4k*R(J+1,k)-R(J,k)/(4k-1); end err=abs(R(J+1,J+1)-R(J+1,J); n=2*n;endR;T=R(J+1,J+1) End其中輸入項(xiàng)：f為被積函數(shù)，ab為積分區(qū)間的端點(diǎn)值，ep為積分精度；輸出項(xiàng)：T是逐次積分表值，n是迭代次數(shù)，R是最后積分值。4.1程序調(diào)用 可以將被積分函數(shù)編成函數(shù)文件，也可以直接使用內(nèi)聯(lián)函數(shù)來表示被積

5、分函數(shù)，示例如下:>>f=inline('1/(1+x.2)','x');>> T,n,R=romb(f,2,9,1e-9)運(yùn)行后得出其迭代次數(shù)，最終積分結(jié)果以及龍貝格積分矩陣如表2-1所示，迭代次數(shù)N=64，最終的積分值R=0.3530.0.7427 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 0.4833 0.3969 0.0000 0.0000 0.0000 0.0000 0.0000 0.3905 0.3596 0.3571 0.0000 0.0000 0.0000 0.0000 0.362

6、8 0.3536 0.3532 0.3532 0.0000 0.0000 0.0000 0.3555 0.3530 0.3530 0.3530 0.3530 0.0000 0.0000 0.3536 0.3530 0.3530 0.3530 0.3530 0.3530 0.0000 0.3531 0.3530 0.3530 0.3530 0.3530 0.3530 0.3530 表2-1 龍貝格積分矩陣3.課本例題求解1 當(dāng)?shù)萫p=1e-9的條件下，迭代次數(shù)N=32，迭代結(jié)果R=0.6931表2-2 式1對應(yīng)的龍貝格積分矩陣0.7500 0.0000 0.0000 0.0000 0.00

7、00 0.0000 0.7083 0.6944 0.0000 0.0000 0.0000 0.0000 0.6970 0.6933 0.6932 0.0000 0.0000 0.0000 0.6941 0.6932 0.6931 0.6931 0.0000 0.0000 0.6934 0.6931 0.6931 0.6931 0.6931 0.0000 0.6932 0.6931 0.6931 0.6931 0.6931 0.6931 2 當(dāng)?shù)萫p=1e-9的條件下，迭代次數(shù)N=32，迭代結(jié)果R=0.2722.表2-3 式2對應(yīng)的龍貝格積分矩陣0.1733 0.0000 0.0000 0

8、.0000 0.0000 0.0000 0.2488 0.2740 0.0000 0.0000 0.0000 0.0000 0.2665 0.2723 0.2722 0.0000 0.0000 0.0000 0.2708 0.2722 0.2722 0.2722 0.0000 0.0000 0.2718 0.2722 0.2722 0.2722 0.2722 0.0000 0.2721 0.2722 0.2722 0.2722 0.2722 0.2722 3對于積分，由于積分下限0為其奇點(diǎn)，理論上無法進(jìn)行數(shù)值積分，本題中近似取下限為1*10-9來進(jìn)行計(jì)算。當(dāng)?shù)萫p=1e-9的條件下，迭代

9、次數(shù)N=16，迭代結(jié)果R=0.2722.表2-4 式3對應(yīng)的龍貝格積分矩陣0.8466 0.0000 0.0000 0.0000 0.0000 0.8288 0.8228 0.0000 0.0000 0.0000 0.8241 0.8225 0.8225 0.0000 0.0000 0.8229 0.8225 0.8225 0.8225 0.0000 0.8226 0.8225 0.8225 0.8225 0.8225 4.對于積分，同樣積分下限0為積分函數(shù)的奇點(diǎn)，理論上無法進(jìn)行數(shù)值積分運(yùn)算，本題中仍取積分下限近似為1*10-9進(jìn)行計(jì)算。當(dāng)?shù)萫p=1e-9的條件下，迭代次數(shù)N=16，迭代結(jié)果R=1.3708.表2-5 式4對應(yīng)的龍貝格積分矩陣1.2854 0.0000 0.0000 0.0000 0.