罗兰在MATLAB的艺术

把想法变成MATLAB

请注意

罗兰在MATLAB的艺术已经退休,不会被更新。

双集成在MATLAB -方法和处理不连续,奇点,等等

在我们的最近的帖子迈克,何西阿书和我谈到调整的绝对和相对公差获得更准确的结果在计算二重积分。今天我们想谈谈选择的方法集成顺序以及选择的第一个维度(方向)的集成。

内容

创造了条件

integral2目前,两个不同的集成方法,“瓦”“迭代”,不包括“汽车”方法选择。

每个集成方法采用一种“分治”双集成方法但以非常不同的方式。的“瓦”方法是基于quad2d将该区域划分为象限的方法和近似积分每个象限二维正交法则。如果错误条件不满足一个矩形,矩形分为象限等等。的“迭代”方法,顾名思义,执行集成迭代一维积分,所以它的“各个击破”的策略是解决一维,然后转移到另一个。

那么,为什么有两个方法?每个人都有不同的优点。通常平铺的方法是越快,但它也有它的缺点。作为一个例子,考虑整合的方式在非矩形区域曾在MATLAB完成。因为dblquad只会对矩形积分,将“面具”输入函数零在一个矩形区域。例如,集成F在一个三角形的顶点(0,0),(1,0)和(1,1),你会写

F = @ (x, y) exp (10 * x)。* tan (y) Q = dblquad (@ (x, y) F (x, y)。* (y < = x)、0、1 0 1)
F = @ (x, y) exp (10 * x)。* tan (y) Q = 1072.6

现在你可以写

Q = integral2 (F, 0 1 0, @ x (x))
Q = 1072.6

这是一个多方便。告诉积分器边界在哪里,给了它一个巨大的优势。相反,掩盖了被积函数有毁灭性影响性能和精度。底层正交规则给最好的结果,当被积函数光滑,和被积函数的不连续急剧下降为零将导致自适应正交算法进一步细分和细分,直到误差估计在不连续是足够小。

这是特别的问题quad2dintegral2“瓦”方法。与“迭代”方法,任何曲线不连续被一个间隔覆盖在任何给定的集成,但在两个维度,细分每个瓷砖成4块,积分器最终覆盖大量的小矩形曲线路径。通常限制功能评估的数量将达到之前波兰了。

我们可以使用可视化quad2dFailurePlot选择和使用MaxFunEvals选项来控制多少工作完成之前放弃的代码。假设

F = @ (x, y) sin (x / y。* y + 1))。* < = x (y)
F = @ (x, y) sin (x / y (y。* + 1)) * < = x (y)。

然后

Q = quad2d (F, 0 1 0 1“MaxFunEvals”现年60岁的“FailurePlot”,真正的)
警告:达到的最大数量评估函数(60)。测试结果失败全球错误。Q = 0.26562

没有什么特别之处更大的矩形底部。quad2d优先攻击的矩形最大错误,所以当一个大矩形刚刚失败错误测试,它挂在如果积分器在其他地方也更大的问题。更大的问题是所指的小矩形。现在,如果我们增加允许的数量评估函数,我们看到quad2d与边界仍在苦苦挣扎。

Q = quad2d (F, 0 1 0 1“MaxFunEvals”,600,“FailurePlot”,真正的)
警告:达到的最大数量评估函数(600)。测试结果失败全球错误。Q = 0.26518

integral2“瓦”方法是基于quad2d的算法,它有同样的问题。

Q = integral2 (F, 0 1 0,1)
警告:达到的最大数量评估函数(10000)。测试结果失败全球错误。Q = 0.26514

虽然需要一定的时间,“迭代”方法可以处理这个问题。

抽搐Q = integral2 (F, 0 1 0 1“方法”,“迭代”)tq = toc
Q = 0.26513 tq = 5.0793

如果你发现“迭代”方法似乎缓慢,考虑放宽公差。integral2使用相同的公差subintegrations,事实证明,这是一个非常保守的事情。用默认的公差,前面的示例接管了4秒机器我们发布的博客。

抽搐Q = integral2 (F, 0 1 0 1“方法”,“迭代”,“AbsTol”1 e-6“RelTol”1 e - 3) tq (2) = toc
Q = 0.26514 tq = 5.0793 - 0.14705

这花了不到十分之一秒的时间在同一台机器上。

改变积分的顺序

不连续可能蛇通过集成的地区,但是偶尔他们可能是线性的和一个轴平行。在这种情况下,当使用积分的顺序问题“迭代”方法。由于每个外部集成执行许多内部集成,它往往是更有用的外部集成高效的如果可能的话,这意味着让内部集成照顾跳跃间断点。在这个例子中被积函数函数之间的交替z = 0z = y这取决于地板(x)是奇数还是偶数。

f = @ (x, y)国防部(地板(x), 2)。* y tic integral2 (20 f, 0, 0, 1,“方法”,“迭代”)目录
f = @ (x, y)国防部(地板(x), 2)。* y ans = 5运行时间是5.625905秒。

工作,但“内在”集成的y方向,对于一个给定的x这意味着函数等于零或一条直线。很容易集成。然而,这意味着外积分的被积函数不连续跳跃从0到0.5取决于地板(x)是奇数还是偶数。我们宁愿内积分处理任何不连续(或至少大部分)。改变积分的顺序是很容易的。

抽搐integral2 (@ (y, x) f (x, y), 0, - 1, 0, 20日“方法”,“迭代”)目录
ans = 5运行时间是0.264443秒。

奇异点

不像dblquad,integral2可以在边界积分奇异性如果他们不太严重。如果你的问题有一个奇点的内部地区的集成,你应该把积分成碎片,使所有奇点位于集成的边界地区。例如,

f = @ (x, y) 1. /√(abs (x。* y))
f = @ (x, y) 1. /√(abs (x。* y))

如果我们把这个广场1 < = x < = 1,1 < = y < = 1,我们有一个问题,因为f是沿着两个单数xy轴。很明显f对两轴是对称的,所以我们只需要积分0 < = x < = 1,0 < = y < = 1和乘以4,但让我们直接把这个放在一边,试着做它。

integral2 (f, 1, 1, 1, 1,“RelTol”1 e-12“AbsTol”1 e-14)
警告:达到的最大数量评估函数(10000)。测试结果失败全球错误。ans = 15.869

“迭代”方法不会拯救我们。

integral2 (f, 1, 1, 1, 1,“RelTol”1 e-12“AbsTol”1 e-14“方法”,“迭代”)
警告:到达间隔的最大数量限制在使用。近似误差是9.0 e-09上界。积分可能不存在,或者它可能很难近似数值精度要求。警告:集成是不成功的。ans =南

但是如果我们把集成分成四块,将边界上的奇点,绝对是程序的问题。

Q = 0;Q = Q + integral2 (f, 1 0 1 0“RelTol”1 e-12“AbsTol”1 e-14);Q = Q + integral2 (f, 1 0 0 1“RelTol”1 e-12“AbsTol”1 e-14);Q = Q + integral2 (f, 0 1 1 0“RelTol”1 e-12“AbsTol”1 e-14);Q = Q + integral2 (f, 0 1 0 1“RelTol”1 e-12“AbsTol”1 e-14)
Q = 16

你能利用这些新的集成程序吗?

你需要整合什么样的功能?会的灵活性integral2和其他的同伴帮你吗?让我们知道在这里




发表与MATLAB®R2013b

|
  • 打印
  • 发送电子邮件