一般区域n重积分MATLAB计算方法[有代码]
这里说的n理论上讲不限大小，但是鉴于n重积分的特点，如果积分重数高，同时积分区域范围较大，数值积分计算量十分巨大，如果各重积分原函数都能显式表示，采用符号积分是第一选择。当只能数值积分的时候，很多时候4重5重积分计算时间都较长。所以虽说不限制n的大小，但是能实际应用中n一般不会很大。下面介绍下这个函数：f = nIntegrate(fun,low,up)
f，函数返回值是n重积分积分结果。fun是被积函数字符串形式，不同的变量依次以x1,x2,...xn表示，（注意，必须以x1,x2,...xn这种形式，其余像y1,y2,...yn或是其他表示方法都不行），low和up都是长度为n的cell数组，low存储从外到内各重积分的积分下限函数，up存储从外到内各重积分的积分上限函数，（都是字符串形式）。举例说明：
计算如下积分：x1*x2*x3，x1从1到2，x2从x1到2*x1，x3从x1*x2到2*x1*x2，构造输入变量
fun = 'x1*x2*x3';
up = {'2','2*x1','2*x1*x2'};
low = {'1','x1','x1*x2'};
然后按f = nIntegrate(fun,low,up)调用即可。
需要说明的是这次没打算上传源m文件，而是上传的加密后的p文件，但是不影响函数使用。该函数可运行7.0以后的MATLAB版本，其中又以2009a为界设置了不同的求解方法。如果用的是2009a,那么速度明显会比09a以前的版本快不少。譬如下面附件里比较复杂的问题：在我的电脑（CPU 酷睿2 T8100 2G内存）上用2009a算得到结果要560多秒，而用08a要5000多秒。需要说明的是由于函数是调用现成的MATLAB函数来编写n重积分，所以函数运行效率也不是很高。
[*]
[*]fun5 = 'sin(x1*exp(x2*sqrt(x3)))+x4^x5'
[*]up5 = {'1','exp(x1)','x1+sin(x2)','x1+x3','2*x4'}
[*]low5 = {'0','exp(x1)/2','x1/2+sin(x2)/2','x1/2+x3/2','x4'}
[*]tic;f = nIntegrate(fun5,low5,up5),toc
[*]fun5 =
[*]sin(x1*exp(x2*sqrt(x3)))+x4^x5
[*]up5 =
[*] '1' 'exp(x1)' 'x1+sin(x2)' 'x1+x3' '2*x4'
[*]low5 =
[*]Columns 1 through 4
[*] '0' 'exp(x1)/2' 'x1/2+sin(x2)/2' 'x1/2+x3/2'
[*]Column 5
[*] 'x4'
[*]f =
[*] 5.7981
[*]Elapsed time is 564.504780 seconds.
1＃的积分用Forcal求解：
高斯算法（通常耗时较长）：
[*]mvar:
[*]f(x1,x2,x3,x4,x5)=sin{x1*exp}+x4^x5;
[*]yx(j,x1,x2,x3,x4,x5,y0,y1)=
[*]{
[*] which{
[*] j==0:,
[*] j==1:,
[*] j==2:/2, y1=x1+sin(x2)],
[*] j==3:,
[*] j==4:
[*] }
[*]};
[*]t0=sys::clock();
[*]XSLSF::igaus;//一般，区间分隔数10,10,10,10,10越大越精确，但耗时越长
[*]/1000;
[*]
复制代码
结果：
5.798135338576838
耗时113.875秒
用连分式法混合了IMSL算法（精度取1e-6）：
[*]mvar:
[*]Fx4x5(x4,x5)=sin{x1*exp}+x4^x5;
[*]x4x5(x4,y0,y1)=
[*]{
[*] y0=x4,
[*] y1=2*x4
[*]};
[*]Fx2x3(_x2,_x3)= x2=_x2, x3=_x3, XSLSF::pqg2;
[*]Gx2(x2)=/2;
[*]Hx2(x2)=x1+sin(x2);
[*]Fx1(_x1)= x1=_x1, IMSL::TWODQ;
[*]t0=sys::clock();
[*]IMSL::QDAGS;
[*]/1000;
[*]
复制代码
结果：
5.798135342841069 8＃的积分用Forcal求解：
[*]mvar: a=1, omiga=2, c1=3, c2=4;
[*]Ft(tao)=exp(-tao*tao);
[*]F(t)=cos(omiga*t)*{c1*exp(-t*t)-c2*XSLSF::fpqg};
[*]IMSL::QDAGS;
结果：-0.1417543457105889
完全用IMSL函数，且相对精度取1e-6：
[*]mvar: a=1, omiga=2, c1=3, c2=4;
[*]Ft(tao)=exp(-tao*tao);
[*]F(t)=cos(omiga*t)*{c1*exp(-t*t)-c2*IMSL::QDAGS};
[*]IMSL::QDAGS;
结果：
-0.1417543562201364
完全使用IMSL的Forcal解法
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
用单重积分构造5重积分（相对精度取1e-1）：
[*]mvar:
[*]Fx5(x5)=sin{x1*exp}+x4^x5;
[*]Fx4(_x4)= x4=_x4, IMSL::QDAGS;
[*]Fx3(_x3)= x3=_x3, IMSL::QDAGS;
[*]Fx2(_x2)= x2=_x2, IMSL::QDAGS/2,x1+sin(x2),0,1e-1,0];
[*]Fx1(_x1)= x1=_x1, IMSL::QDAGS;
[*]t0=sys::clock();
[*]IMSL::QDAGS;
[*]/1000;
结果：
5.798135465850306
时间2.953秒
如果相对精度取1e-2或更小，计算时间延长，但结果越来越不准确，这个还不知道为什么？看来用单重构造多重不是很好的方法。
＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝
用单重＋二重＋二重积分构造5重积分（相对精度取1e-3）：
[*]mvar:
[*]Fx4x5(x4,x5)=sin{x1*exp}+x4^x5;
[*]Gx4(x4)=x4;
[*]Hx4(x4)=2*x4;
[*]Fx2x3(_x2,_x3)= x2=_x2, x3=_x3, IMSL::TWODQ;
[*]Gx2(x2)=/2;
[*]Hx2(x2)=x1+sin(x2);
[*]Fx1(_x1)= x1=_x1, IMSL::TWODQ;
[*]t0=sys::clock();
[*]IMSL::QDAGS;
[*]/1000;
结果：
5.798135365426686
耗时6.141秒。
相对精度取1e-2或1e-3，计算结果是一样的，时间也几乎相等；但相对精度取1e-4或更小时，计算时间延长，且结果越来越不准确，这个还不知道为什么？
用单重积分构造5重积分，似乎不如用单重＋二重＋二重积分构造5重积分精确。
页:
[1]