MATLAB学习笔记：数值积分

当

（1）被积函数的原函数不能用初等函数表示。

（2）被积函数难以用公式表示，而是用图形或表格给出的。

就应该建立定积分的近似计算方法：数值积分方法。

梯形法：

z=trapz(x,y)

>> x=0:0.5:1;
>> y=exp(-x.^2);
>> z=trapz(x,y)

z =

    0.7314

>> x=0:0.05:1;
>> y=exp(-x.^2);
>> z=trapz(x,y)

z =

    0.7467

抛物线法：Simpson（辛普森）公式

自适应辛普森公式:

z=quad(f,a,b,tol)

a是积分下限，b是积分上限，tol是计算精度，缺省为0.001。

调用quad函数时，先要建立一个描述被积函数的函数文件或语句函数。

 

以下分别为：内联函数和函数句柄

>> g=inline('exp(-x.^2)');
>> z=quad(g,-1,1)

z =

    1.4936

>> f=@(x) exp(-x.^2);
>> z=quad(f,-1,1)

z =

    1.4936

分别使用左端点的矩形法，右端点的矩形法和梯形法

>> n=100;
>> x=linspace(0,1,n);
>> y=4./(1+x.^2);

>> jxl=sum(y(1:(n-1)))./n

jxl =

    3.1202

>> jxr=sum(y(2:n))./n

jxr =

    3.1002

>> tx=trapz(x,y)

tx =

    3.1416

还有两个常用的求数值积分的函数：

z=quadl(f,a,b,tol)  自适应复合Lobatto数值积分法

z=quadgk(f,a,c,tol)  自适应复合Gauss-Kronrod数值积分法，适用于高精度和震荡数值积分，以及广义数值积分

 

>> g=inline('4./(1+x.^2)');
>> quad(g,0,1)

ans =

    3.1416

>> quadl(g,0,1)

ans =

    3.1416

>> g=@(x) 4./(1+x.^2);
>> quadgk(g,0,1)

ans =

    3.1416

>> f1=@(x) 1./(x.^2);
>> z1=quadgk(f1,1,inf)

z1 =

     1

>> f2=@(x)exp(-x.^2);
>> z2=quadgk(f2,-inf,inf)

z2 =

    1.7725
>> f3=@(x)1./(x.*sqrt(1-log(x).^2));
>> z3=quadgk(f3,1,exp(1))

z3 =

    1.5708