QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 4541|回复: 0
打印 上一主题 下一主题

告诉你!Matlab中几个数值积分函数的比较和优缺点

[复制链接]
字体大小: 正常 放大

398

主题

13

听众

1342

积分

  • TA的每日心情
    慵懒
    2015-12-12 14:33
  • 签到天数: 81 天

    [LV.6]常住居民II

    跳转到指定楼层
    1#
    发表于 2015-10-15 23:50 |只看该作者 |倒序浏览
    |招呼Ta 关注Ta

    一、Z = trapz(X,Y,dim)
    梯形数值积分,通过已知参数x,y按dim维使用梯形公式进行积分

    例1 计算int(sin(x),0,pi)

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>x=0:pi/100:2*pi;
    >>y=sin(x);
    >>z=trapz(x,y)%或者说使用z = pi/100*trapz(y)
    z =

    1.0300e-017

    >>z = pi/100*trapz(y)

    二、[q,fcnt]= quad(fun,a,b,tol,trace,p1,p2...)
    自适应simpson公式数值积分,适用于精度要求低,被积函数平滑性较差的数值积分

    注意事项:
    1.被积函数fun必须是函数句柄
    2.积分限[a,b]必须是有限的,因此不能为inf
    3.p1为其他需要传递的参数,一般是数值

    可能警告:
    1.'Minimum step size reached'
    意味着子区间的长度与计算机舍入误差相当,无法继续计算了。原因可能是有不可积的奇点
    2.'Maximum function count exceeded'
    意味着积分递归计算超过了10000次。原因可能是有不可积的奇点
    3.'Infinite or Not-a-Number function value encountered'
    意味着在积分计算时,区间内出现了浮点数溢出或者被零除。

    2 计算积分1/(x^3-2*x-p),其中参数p=5,积分区间为[0,2]

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F = @(x,n)1./(x.^3-2*x-n);
    >>Q = quad(@(x)F(x,5),0,2)%或者使用 quad(F,0,2,[],[],5)效果是一样的,只是前者使用的函数嵌套

    Q =

    -0.4605

    >>quad(F,0,2,[],[],5)

    ans =

    -0.4605

    三、[q,fcnt] = quadl(fun,a,b,tol,trace,p1,p2...)
    自适应Lobatto数值积分,适用于精度要求高,被积函数曲线比较平滑的数值积分

    注意事项:
    同quad

    可能警告:
    同quad

    3 计算积分1/(x^3-2*x-p),其中参数p=5,积分区间为[0,2]

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F=@(x,p)1./(x.^3-2*x-p);
    >>Q = quadl(F,0,2,[],[],5)%或者Q = quadl(@(x)F(x,5),0,2)

    Q =

    -0.4605

    四、[q,errbnd] = quadgk(fun,a,b,param1,val1,param2,val2,...)
    自适应Gauss-Kronrod数值积分,适用于高精度和震荡数值积分,支持无穷区间,并且能够处理端点包含奇点的情况,同时还支持沿着不连续函数积分,复数域线性路径的围道积分法

    注意事项:
    1.积分限[a,b]可以是[-inf,inf],但必须快速衰减
    2.被积函数在端点可以有奇点,如果区间内部有奇点,将以奇点区间划分成多个,也就是说奇点只能出现在端点上
    3.被积函数可以剧烈震荡
    4.可以计算不连续积分,此时需要用到'Waypoints'参数,'Waypoints'中的点必须严格单调
    5.可以计算围道积分,此时需要用到'Waypoints'参数,并且为复数,各点之间使用直线连接
    6.param,val为函数的其它控制参数,比如上面的'waypoints'就是,具体看帮助

    出现错误:
    1.'Reached the limit on the maximum number of intervals in use'
    2.'Infinite or Not-a-Number function value encountered'

    4 计算有奇点积分int(exp(x)*log(x),0,1)

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F=@(x)exp(x).*log(x);%奇点必须在端点上,否则请先进行区间划分
    >>Q = quadgk(F,0,1)

    Q =

    -1.3179

    例5 计算半无限震荡积分int(x^5*exp(-x)*sin(x),0,inf)

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F=@(x)x.^5.*exp(-x).*sin(x);
    >>fplot(F,[0,100])%绘图,看看函数的图形
    >>[q,errbnd] = quadgk(F,0,inf,'RelTol',1e-8,'AbsTol',1e-12)%积分限中可以有inf,但必须快速收敛

    q =

    -15.0000


    errbnd =

    9.4386e-009

    例6 计算不连续积分,积分函数为f(x)=x^5*exp(-x)*sin(x),但是人为定义f(2)=1000,f(5)=-100,积分区间为[1 10]

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F=@(x)x.^5.*exp(-x).*sin(x);
    >>[q,errbnd] = quadgk(F,1,10,'Waypoints',[2 5])%显然2,5为间断点

    q =

    -10.9408


    errbnd =

    3.2296e-014

    例7 计算围道积分,在复数域内,积分函数1/(2*z-1),积分路径为由[-1-i 1-i 1+i -1+i -1-i]围成的矩形边框

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>Waypoints=[-1-i 1-i 1+i -1+i -1-i];
    >>plot(Waypoints);%绘制积分路径
    >>xlabel('Real axis');ylabel('Image axis');axis([-1.5 1.5 -1.5 1.5]);grid on;
    >>Q = quadgk(@(z)1./(2*z - 1),-1-i,-1-i,'Waypoints',[1-i,1+i,-1+i])%注意各点间使用直线连接

    ans =

    0.0000 + 3.1416i

    >> quadgk(@(z)1./(2*z - 1),-1-i,-1-i,'Waypoints',Waypoints)%使用这个的效果也是一样的,就是说始末点可以随便包不包含在Waypoints中

    ans =

    0.0000 + 3.1416i

    五、[Q,fcnt] = quadv(fun,a,b,tol,trace)
    矢量化自适应simpson数值积分

    注意事项:
    1.该函将quad函数矢量化了,就是一次可以计算多个积分
    2.所有的要求完全与quad相同

    8 计算下面积分,分别计算n=1,2...,5时的5个积分值,被积函数1/(n+x),积分限为[0,1]

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>for k = 1:5, Qs(k) = quadv(@(x)1/(k+x),0,1);end;Qs

    Qs =

    0.6931 0.4055 0.2877 0.2231 0.1823

    >>F=@(x,n)1./((1:n)+x);%定义被积函数
    >>quadv(@(x)F(x,5),0,1)%我们可以完全使用quadv函数替换上面循环语句的,建议使用后者

    ans =

    0.6931 0.4055 0.2877 0.2231 0.1823

    六、q = dblquad(fun,xmin,xmax,ymin,ymax,tol,method)
    矩形区域二重数值积分,一般区域二重积分参见NIT(数值积分工具箱)的quad2dggen函数

    9 计算下面二重积分

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F = @(x,y)y*sin(x)+x*cos(y);
    >Q = dblquad(F,pi,2*pi,0,pi)

    Q =

    -9.8696

    七、q=triplequad(fun,xmin,xmax,ymin,ymax,zmin,zmax,tol,method)
    长方体区域三重数值积分,注意此时没有一般区域的三重积分

    例10 计算下面三重积分

    %by dynamic
    %all rights reserved by www.matlabsky.com
    >>F = @(x,y,z)y*sin(x)+z*cos(x);
    >>Q = triplequad(F,0,pi,0,1,-1,1)

    Q =

    2.0000

    八、超维长方体区域多重积分
    quadndg:NIT工具箱函数,可以解决多重超维长方体边界的定积分问题,但没有现成的一般积分区域求解函数


    下面总结下:
    (1)quad:采用自适应变步长simpson方法,速度和精度都是最差的,建议不要使用
    (2)quad8:使用8阶Newton-Cotes算法,精度和速度均优于quad,但在目前版本下已被取消
    (3)quadl:采用lobbato算法,精度和速度均较好,建议全部使用该函数
    (4)quadg:NIT(数值积分)工具箱函数,效率最高,但该工具箱需要另外下载
    (5)quadv:quad的矢量化函数,可以同时计算多个积分
    (6)quadgk:很有用的函数,功能在Matlab中最强大
    (7)quad2dggen:一般区域二重积分,效率很好,需要NIT支持
    (8)dblquad:长方形区域二重积分
    (9)triplequadL:长方体区域三重积分
    (10)quadndg:超维长方体区域积分,需要NIT支持

    NIT数值积分工具箱下载参见这里http://www.matlabsky.com/thread-225-1-2.html



    zan
    转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信
    您需要登录后才可以回帖 登录 | 注册地址

    qq
    收缩
    • 电话咨询

    • 04714969085
    fastpost

    关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

    手机版|Archiver| |繁體中文 手机客户端  

    蒙公网安备 15010502000194号

    Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

    GMT+8, 2025-7-21 08:33 , Processed in 0.278733 second(s), 49 queries .

    回顶部