matlab數值積分中函式積分的4種方法
1. 採用inline行內函數
Matlab中可以有采用幾種不同的方式來指定被積函式。對於簡單的、長度不超過一行的公式採用inline命令比較方便。例如,
可用下面的語句進行計算
>> f=inline('1/sqrt(1+x^4)') %採用inline行內函數
f =
行內函數:
f(x) = 1/sqrt(1+x^4)
>> q=quadtx(f,0,1)
q =
0.9270
從matlab第七版開始,內聯(inline)物件被一種功能更強大的結構匿名函式(anonymousfunction)所替代。在matlab第七版內聯物件還允許被使用,但推薦用匿名函式,因為後者可以生成更高效率的程式程式碼。採用匿名函式,上面的例子變為
>> [email protected](x) 1/sqrt(1+x^4) %採用匿名函式@(x)
f =
@(x)1/sqrt(1+x^4)
>> q=quadtx(f,0,1)
q =
0.9270
2. 特殊點不可積函式,採用realmin
如果我們想要計算
可能使用下面的語句
>> f=inline('sin(x)/x')
f =
行內函數:
f(x) = sin(x)/x
>> q=quadtx(f,0,pi)
已達到最大遞迴限制 500。使用set(0,'RecursionLimit',N) 可更改此限制。請注意,超出可用堆疊空間可能會使 MATLAB 和/或計算機崩潰。
出錯 inlineeval
此時inline函式積分時好像出現了問題,因為在計算f(0)時,出現了除以0的情況,並且最終產生遞迴限制錯誤。一種補救的方法是,將積分的下限由0變為最小的正浮點數,realmin。
>> q=quadtx(f,realmin,pi)
q =
1.8519
這樣就可以避免遞迴錯誤情況的發生。
3. 使用m檔案寫函式
針對上面出現的問題,另一種方法是使用M檔案,而不是行內函數。建立包含下面程式的檔案sinc.m
function f=sinc(x)
%函式sin(x)/x
if x==0
f=1;
else
f=sin(x)/x;
end
然後使用函式控制代碼計算積分
>> Q=quadtx(@sinc,0,pi)
Q =
1.8519
4. 依賴於引數的積分
一個典型的例子是β函式,它定義為
matlab中已經實現了一個現成的β函式,但我們可以以它為例,說明如何處理積分中的引數。建立一個帶三個引數的行內函數
>> F=inline('t^(z-1)*(1-t)^(w-1)','t','z','w')
或者建立一個M檔案:
function f=betaf(t,z,w)
f= t^(z-1)*(1-t)^(w-1)
並將其命名為betaf.m。
就像任何函式一樣,引數的順序是很重要的。定義被積函式時,必須讓積分變數為其第一個引數。然後給出其他引數的值,作為傳遞給quadtx的附加引數。要計算β(8/3,10/3)
應該先設
>> z=8/3;
>> w=10/3;
>> tol=1e-6;
然後執行命令
>> Q=quadtx(F,0,1,tol,z,w)
Q =
0.0348
或Q=quadtx(@betaf,0,1,tol,z,w)