2023/2/2數值積分計算吳鵬（rocwoods）rocwoods@126.comMATLAB從零到進階2023/2/2主要內容矩形區域積分以及離散數據積分

含參數的積分一般區域二重、三重積分一般區域n重積分

蒙特卡洛法計算n重積分

2023/2/2第一節矩形區域積分以及離散數據積分2023/2/2一、矩形區域積分積分上下限是常數，分單重、二重、三重三種情況1.單重quad(自適應Simpson積分)quadl(自適應Gauss-Lobatto積分，最常用)quadgk(自適應Gauss-Kronrod積分，尤其適合震蕩積分、含奇點的積分；從R2007b開始支持)

具體用法參見書中【例10.1-1】2.二重情形一般二重積分可以用下面兩個函數解決：dblquadquad2d(R2009a開始支持，不僅可以求矩形區域二重積分，還可以求一般區域二重積分)具體用法參見書中【例10.1-2】2023/2/2一、矩形區域積分3.三重情形MATLAB中可以利用triplequad函數求解長方體區域的三重積分。具體用法參見書中【例10.1-3】4.向量化積分所謂向量化積分，是指被積函數含有參數，需要對參數的一系列值求出相應的積分。向量話積分可以用quadv函數。具體用法參見書中【例10.1-4】2023/2/2二、離散數據積分離散數據積分的函數MATLAB只有trapz函數，是針對一重情形的，如果是二重情形，需要自己編寫，其中會用到trapz函數。用法示例參見書中【例10.1-5】以及【例10.1-6】2023/2/2第二節含參數的積分2023/2/2

在解決實際問題的時候我們常要求解一些帶參數的積分問題，這涉及到參數的傳遞，總的說來有以下幾種方法可以求解：用inline+num2str方法(MATLAB7.0以前版本使用，7.0以后不推薦使用)利用匿名函數實現利用嵌套函數實現利用積分函數本身傳遞參數 這幾種方法的具體實現辦法可以參見【例10.2-1】。2023/2/2第三節一般區域二重、三重積分2023/2/2一、概要

一般區域上的二重積分，低版本的MATLAB向來支持不好，7.X版本之前不能通過簡單的形式直接求取一般區域上的的二重積分，往往要借助廣泛流傳的NIT工具箱（數值積分工具箱）來實現。后來到了7.X版本，MATLAB引入了匿名函數結構，dblquad的被積函數可以是匿名函數的形式，利用匿名函數結構，通過適當改寫被積函數。dblquad可以求解一般區域上的二重積分。但這種改寫就是將積分區域表示成邏輯表達式形式并與被積函數相乘，本質上還是計算的矩形區域的積分。 仔細分析dblquad的幫助文檔中求一般區域二重積分的例子就會發現，這種辦法是把原被積函數外推到了一個矩形區域上來實現的，該矩形區域包含了被積區域，在被積區域上，外推函數取值和原函數一樣，而在矩形區域內被積區域外的那部分，外推函數取值都為0.這樣外推函數在矩形區域上積分的結果就等于原被積函數在被積區域的結果。這種數學上簡單的變換導致的結果就是，按這種辦法，計算量大增，特別是被積區域非常不規則并布滿整個矩形區域的時候。2023/2/2一、概要

從MATLABR2009a版本起，MATLAB終于有了專門求解一般區域二重積分的函數——quad2d，該函數采用的自適應積分算法基于L.F.Shampine的文章："VectorizedAdaptiveQuadratureinMATLAB,"JournalofComputationalandAppliedMathematics,211,2008。該函數求解一般區域二重積分的效率要遠高于上面提到的dblquad的方法。 雖然quad2d可以求解一般區域二重積分，但還是不能直接求解一般區域三重積分，而NIT工具箱也沒有一般區域三重積分的計算函數。這些需要靈活運用MATLAB相關函數來構建。2023/2/2二、一般區域二重積分的計算dblquad求解一般區域二重三重積分的思路方法，就是將被積函數“延拓”到矩形或者長方體區域，但是這種方法不可避免引入很多乘0運算浪費時間。而quad2d雖然效率很高，但是R2009a之前沒有這個函數。因此，要想在2009a之前的MATLAB7版里里計算一般區域二重三重積分，需要開辟新的方法。新的方法是調用已有的MATLAB函數求解，在求一般區域二重積分時，效率和quad2d相比有一些差距，但是相對于“延拓”函數的做法，效率大大提高了

。【例10.3-2】和【例10.3-3】給出了這一方法的應用實例，可以看出該方法在R2009a以前的版本中不失為一種方法，效率要明顯高于dblquad幫助文檔里推薦的做法。更重要的是，這給我們求解一般區域三重積分提供了一種途徑。2023/2/2三、一般區域三重積分的計算模板1：quadl(@(x)arrayfun(@(xx)quad2d(被積函數f(xx,y,z)關于y,z變量的函數句柄,y積分下限y1(xx),y積分上限y2(xx),z積分下限z1(xx,y),z積分上限z2(xx,y)),x),x積分下限值,x積分上限值)模板2：quad2d(@(x,y)arrayfun(@(xx,yy)quadl(被積函數f(xx,yy,z)關于z變量的函數句柄,z積分下限z1(xx,yy),z積分上限z2(xx,yy)),x,y),x積分下限值,x積分上限值,y積分下限y1(x),y積分上限y2(x))模板3：quadl(@(x)arrayfun(@(xx)quadl(@(y)arrayfun(@(yy)quadl(被積函數f(xx,yy,z)關于z變量的函數句柄,z積分下限z1(xx,yy),z積分上限z2(xx,yy)),y),y積分下限y1(xx),y積分上限y2(xx)),x),x積分下限值,x積分上限值)模板使用說明：x,y,z是積分變量，模板中除了用語言描述的參量用相應表達式替換掉外，其余結構保持不變。具體使用示例可以參考【例10.3-4】

2023/2/2第四節一般區域n重積分的計算2023/2/2介紹求解一般區域n重積分的函數——nIntegratef=nIntegrate(fun,low,up)

f為函數的返回值是n重積分積分結果。fun是被積函數字符串形式，不同的變量依次以x1,x2,...xn表示，（需要注意的是，必須以x1,x2,...xn這種形式表示，其余像y1,y2,...yn或是其他表示方法都不行）。low和up都是長度為n的cell數組，low存儲從外到內各重積分的積分下限函數，up存儲從外到內各重積分的積分上限函數（都是字符串形式）。low和up內的函數表示都要遵循一些原則，這些原則在程序注釋里進行了說明.nIntegrate函數源代碼以及使用示例參考書中相關內容以及【例10.4-1】至【例10.4-4】2023/2/2第五節蒙特卡洛法計算n重積分2023/2/2 一些積分重數比較高的積分，按照傳統方法進行計算的話，運算量會非常大，求解時間隨著重數的增加而迅速增加。在10.4節雖然給出了一般區域n重積分的通用程序，但是在實際計算中，n不會很大，一般n超過5后，求解時間就會長到無法接受。對于更高維的積分，一個非常有效的求解方法就是蒙特卡洛積分法。 蒙特卡洛方法又稱隨機抽樣法或統計試驗方法，它用于求積分時，與積分重數無關，這點非常重要。雖然四維以下的積分用蒙特卡洛法效率可能不如傳統的一些數值積分方法，但是維數高的時候，蒙特卡洛法比傳統方法要有效的多，而且實現起來也非常容易。可以說，計算高維積分是蒙特卡洛方法最成功和典型的應用。 實際應用中，有多種蒙特卡洛方法可以計算n重積分，比較常用的是基本蒙特卡洛法和等分布序列的蒙特卡洛法。2023/2/2