這裡討論的計算方法指的是利用現有的matlab函式來求解,而不是根據具體的數值計算方法來編寫相應程式。目前最新版的2009a有關於一般區域二重積分的計算函式quad2d,但沒有一般區域三重積分的計算函式,而nit工具箱似乎也沒有一般區域三重積分的計算函式。
本貼的目的是介紹一種在7.x版本matlab(不一定是2009a)裡求解一般區域二重三重積分的思路方法。需要說明的是,在matlab的dblquad幫助文件裡已經討論了一種求解一般區域二重三重積分的思路方法,就是將被積函式「延拓」到矩形或者長方體區域,但是這種方法不可避免引入很多乘0運算浪費時間。因此,新的思路將避免這些。由於是呼叫已有的matlab函式求解,在求一般區域二重積分時,效率和2009a的quad2d相比有一些差距,但是相對於"延拓"函式的做法,效率大大提高了。下面結合一些簡單例子說明下計算方法。
譬如二元函式f(x,y) = x*y,y從sin(x)積分到cos(x),x從1積分到2,這個積分可以很容易用符號積分算出結果
syms x
yint(int(x*y,y,sin(x),cos(x)),1,2) ]
結果是
-1/2*cos(1)*sin(1)-1/4*cos(1)^2+cos(2)*sin(2)+1/4*cos(2)^2 =
-0.635412702399943
複製**
如果你用的是2009a,你可以用
quad2d(@(x,y)
x.*y,1,2,@(x)sin(x),@(x)cos(x),'abstol',1e-12)
複製**
得到上述結果。
如果用的不是2009a,那麼你可以利用nit工具箱裡的quad2dggen函式。
那麼我們如果既沒有nit工具箱用的也不是2009a,怎麼辦呢?
答案是我們可以利用兩次quadl函式,注意到quadl函式要求積分表示式必須寫成向量化形式,所以我們構造的函式必須能接受向量輸入。見如下**
function
intdemo
function f1 = myfun1(x)
f1 = zeros(size(x));
for k =
1:length(x)
f1(k) = quadl(@(y)
x(k)*y,sin(x(k)),cos(x(k)));
endend
y =
quadl(@myfun1,1,2)
end
複製**
myfun1函式就是構造的原始被積函式對y積分後的函式,這時候是關於
x的函式,要能接受向量形式的x輸入,所以構造這個函式的時候考慮到x是向量的情況。
quadl(@(x)
arrayfun(@(xx) quadl(@(y)
xx*y,sin(xx),cos(xx)),x),1,2)
複製**
上面這行**體現了用matlab7.x求一般區域二重積分的一般方法。可以這麼理解這句**:
首先@(x)
arrayfun(@(xx) quadl(@(y)
xx*y,sin(xx),cos(xx)),x)
複製**
定義了乙個關於x的匿名函式,供quadl呼叫求最外重(x從1到2的)積分,這時候,x對於
arrayfun(@(xx)
quadl(@(y) xx*y,sin(xx),cos(xx)),x)
複製**
就是已知的了。
而@(xx) quadl(@(y)
xx*y,sin(xx),cos(xx))
複製**
定義的是對於給定的xx,求xx*y關於y的積分函式,這就相當於數學上積完第一重y的積分後得到乙個關於xx的函式
而arrayfun(@(xx)
quadl(@(y) xx*y,sin(xx),cos(xx)),x)
複製**
只是對@(xx) quadl(@(y)
xx*y,sin(xx),cos(xx))
複製**
加了乙個迴圈的殼,保證「積完第一重y的積分後得到乙個關於xx的函式」能夠接受向量化的xx的輸入,從而能夠被quadl呼叫。
有了這個模板,我們可以方便求其他一般積分區域(上下限是函式)形式的二重積分,例如被積函式
f = @(x,y)
exp(sin(x))*ln(y),y從5*x積分到x^2,x從10積分到20。
用quad2d,上面介紹的方法,還有dblquad幫助文件裡給的延拓函式的方法
tic,y1
= quad2d(@(x,y) exp(sin(x)).*log(y),10,20,@(x)5*x,@(x)x.^2),toc
tic,y2 =
quadl(@(x) arrayfun(@(x) quadl(@(y)
exp(sin(x)).*log(y),5*x,x.^2),x),10,20),toc
tic,y3 = dblquad(@(x,y)
exp(sin(x)).*log(y).*(y>=5*x & y<=x.^2),10,20,50,400),toc
y1 =
9.368671342614414e+003
elapsed time is 0.021152 seconds.
y2 =
9.368671342161189e+003
elapsed time is 0.276614 seconds.
y3 =
9.368671498376889e+003
elapsed time is 1.674544
seconds.
複製**
可見上述方法在2009a以外的版本中不失為一種方法,起碼效率高於dblquad幫助文件裡推薦的做法。更重要的是,這給我們求解一般區域三重積分提供了一種途徑。
matlab二重定積分 二重積分 matlab
第六章 用matlab 計算二重積分 由於二重積分可以化成二次積分來進行計算,因此只要確定出幾分區域,就可以反覆 使用int 命令來計算二重積分。例6.4.1 計算二重積分yd ixedxdy d是由直線 x 0,y 1,y x 所圍區域 解該積分可以寫成yy idyxe dx或yy idxxe d...
用python求一重積分和二重積分
首先是對一元函式求積分,使用scipy下的integrate函式 from scipy import integrate def g x return 1 x 2 0.5 用integrate下的quad函式可以同時求出積分結果和誤差 res,err integrate.quad g,1,1 1和1...
二重積分與多元函式結合的一題
分部積分的目的便在於轉換研究物件,在很久以前就總結過分部積分的各種使用場合,今天也翻出來看了看,在證明積分不等式中有非常完美的運用。在做真題時,遇到二重積分結合多元函式的大題,同樣使用一樣的思想 分部積分 產生 消去導數 這是2011真題大題,本題第一感覺 把f 塞到dxdy中,怎麼塞呢?我們只需要...