二重積分的數(shù)值方法

1、精選優(yōu)質(zhì)文檔-傾情為你奉上數(shù)值分析課程設(shè)計(jì)報(bào) 告 專業(yè)： 學(xué)號(hào)： 學(xué)生姓名： 指導(dǎo)教師： 一、題目 數(shù)值積分中二重積分探究。二、理論數(shù)值積分就是用數(shù)值方法近似計(jì)算定積分。其原理很簡(jiǎn)單，就是將積分核用插值多項(xiàng)式替代，用多項(xiàng)式的結(jié)果近似定積分的值。一般常用的方法是，將積分區(qū)間等分為個(gè)子區(qū)間，即取步長，子區(qū)間端點(diǎn)為（k=0,1,，n），在每個(gè)子區(qū)間上套用插值積分公式，再將個(gè)區(qū)間的結(jié)果累加起來。比較常用的有梯形公式，是在每個(gè)子區(qū)間上用1階多項(xiàng)式（即直線段）近似并積分的結(jié)果： 另外一種在實(shí)際應(yīng)用中很受歡迎的方法是，在每個(gè)子區(qū)間上用2階多項(xiàng)式（即拋物線）近似并積分，得到著名的辛普森（Simpson）公式：

2、 其中。 三、方法、算法與程序設(shè)計(jì).辛普森公式求二重積分 考慮二重積分，它是曲面與平面區(qū)域R圍成的體積，對(duì)于矩形區(qū)域，可將它寫成累次積分 。若用復(fù)合辛普森公式，可分別將,分成N,M等份，步長,，先對(duì)積分，應(yīng)用復(fù)合辛普森公式，令，則從而得 。對(duì)每個(gè)積分再分別用復(fù)合辛普森公式即可求得積分值。MATLAB程序見附錄1，MATLAB中自帶自適應(yīng)辛普森公式dblquad()，對(duì)于變量區(qū)域同樣適用。對(duì)于變量區(qū)域，寫成累次積分的形式：進(jìn)行數(shù)值計(jì)算的表達(dá)式為： 上面的表達(dá)式中、表示權(quán)重，取決于一維積分方法。我們常用復(fù)合辛普森公式，先對(duì)內(nèi)積分進(jìn)行計(jì)算，在計(jì)算外積分，與矩形區(qū)域情況基本一致。高斯求積公式求二重積分

3、在高斯求積公式中，若取權(quán)函數(shù)，區(qū)間為，則得公式 。勒讓德多項(xiàng)式是區(qū)間上的正交多項(xiàng)式，因此，勒讓德多多項(xiàng)式的零點(diǎn)就是求積公式的高斯點(diǎn)。 若取的零點(diǎn)做節(jié)點(diǎn)構(gòu)造求積公式 ；若取的零點(diǎn)構(gòu)造求積公式 ；當(dāng)時(shí)，求積公式為 同樣先用高斯求積公式求內(nèi)積分，再求外積分，可得二重積分值。 四、算例、應(yīng)用實(shí)例算例：計(jì)算二重積分。（1）若區(qū)域，試分別用復(fù)合辛普森公式（取n=4）及高斯求積公式（取n=4）求積分。（2）若區(qū)域用復(fù)合辛普森公式（取n=4）求積分。解：（1） = 對(duì)各個(gè)積分應(yīng)用復(fù)合辛普森公式。 也可應(yīng)用MATLAB中的函數(shù)進(jìn)行計(jì)算，程序見附錄2。 先將區(qū)域變換為區(qū)域，其中 ，等價(jià)于，有 。對(duì)于取時(shí)的高斯求積

4、公式節(jié)點(diǎn)及系數(shù)，即， 用的高斯積分公式計(jì)算積分I， （2） ，等分為4等份，對(duì)應(yīng)值為的值，用節(jié)點(diǎn)應(yīng)用辛普森公式對(duì)內(nèi)積分求積，再用復(fù)合辛普森公式對(duì)外積分求積，也可用MATLAB中的函數(shù)實(shí)現(xiàn)，結(jié)果和程序如下（附錄3）。五、參考文獻(xiàn)【1】 數(shù)值分析 李慶揚(yáng)，王朝能，易大義 清華大學(xué)出版社【2】 數(shù)值分析課程設(shè)計(jì) 陳越，童若鋒 浙江大學(xué)出版社【3】 MATLAB教程 張志涌 北京航空航天大學(xué)出版社六、附錄附錄1：function q=DblSimpson(f,a,A,b,B,m,n)if(m=1 && n=1)  %辛普森公式q=(B-b)*(A-a)/9)*(su

5、bs(sym(f),findsym(sym(f),a,b)+.subs(sym(f),findsym(sym(f),a,B)+.subs(sym(f),findsym(sym(f),A,b)+.subs(sym(f),findsym(sym(f),A,B)+.4*subs(sym(f),findsym(sym(f),(A-a)/2,b)+.4*subs(sym(f),findsym(sym(f),(A-a)/2,B)+.4*subs(sym(f),findsym(sym(f),a,(B-b)/2)+. 4*subs(sym(f),findsym(sym(f),A,(B-b)/2)+.16*su

6、bs(sym(f),findsym(sym(f),(A-a)/2,(B-b)/2);else   %復(fù)合辛普森公式q=0;for i=0:n-1for j=0:m-1x=a+2*i*(A-a)/2/n;y=b+2*j*(B-b)/2/m;            x1=a+(2*i+1)*(A-a)/2/n;y1=b+(2*j+1)*(B-b)/2/m;x2=a+2*(i+1)*(A-a)/2/n;y2=b+2*(j+1)*(B-b)/2/m;q=q+subs(sym

7、(f),findsym(sym(f),x,y)+.subs(sym(f),findsym(sym(f),x,y2)+.subs(sym(f),findsym(sym(f),x2,y)+.subs(sym(f),findsym(sym(f),x2,y2)+.4*subs(sym(f),findsym(sym(f),x,y1)+.4*subs(sym(f),findsym(sym(f),x2,y1)+.4*subs(sym(f),findsym(sym(f),x1,y)+.4*subs(sym(f),findsym(sym(f),x1,y2)+.16*subs(sym(f),findsym(sym(f),x1,y1);endendendq=(B-b)*(A-a)/36/m/n)*q;附錄2：fun=inline(exp(x*y);d=dblquad(fun,0,1,0,1);附錄3：f = inline('exp(-x*y),'x','y');xlower = inlin