# quad2d
计算二重数值积分 - tiled 方法
函数库: TyMath
# 语法
q,E = quad2d(fun,a,b,c,d)
q,E = quad2d(fun,a,b,c,d;Key=Value)
# 说明
q, E = quad2d(fun,a,b,c,d) 逼近 fun(x,y) 在平面区域 a≤x≤b 和 c(x)≤y≤d(x) 上的积分 q。边界 c 和 d 均可为标量或函数句柄。并返回绝对误差 E = |q - I| 的逼近上限,其中 I 是积分的确切值。示例
q, E = quad2d(fun,a,b,c,d,Key=Value) 指定具有一个或多个 Key=Value 关键字参数的其他选项。例如,您可以指定 AbsTol 和 RelTol 来调整算法必须满足的误差阈值。示例
# 示例
计算二重积分
对方程求积分:
using TyMath
fun = (x, y) -> y*sin(x) + x*cos(y)
Q, = quad2d(fun, pi, 2 * pi, 0, pi)
Q = -9.869604400888917
将结果与积分的真值
-pi^2
ans = -9.869604400888917
位于积分边界具有奇异性的被积函数
求函数的积分:
using TyMath
fun = (x, y) -> 1/(sqrt(x+y)*(1+x+y)^2)
ymin = 0
ymax = x -> 1 - x
Q, = quad2d(fun, 0, 1, ymin, ymax)
Q = 0.2853982593844486
积分的真值为
pi/4 - 0.5
ans = 0.2853981633974483
quad2d 的局限性
quad2d 首先将积分区域映射到矩形。如果区域没有四条边或它的某条边不能顺利映射到直线,对此区域求积分时必然出现问题。如果求积分失败,可以通过几种有用的方法补救:将 Singular 设置为其默认值 true,在笛卡尔坐标与极坐标之间切换,或将积分区域分为几个小块并将这些小块的积分结果相加。
例如:
using TyMath
fun = (x, y) -> abs(x^2 + y^2 - 0.25)
c = x -> -sqrt(1 - x^2)
d = x -> sqrt(1 - x^2)
quad2d(fun, -1, 1, c, d; AbsTol=1e-8, FailurePlot=true, Singular=false)
┌ Warning: Reached the maximum number of function evaluations (2000). The result fails the global error test.
失败绘图显示两个问题区域,一个在点 (-1,0) 和 (1,0) 附近,一个在圆
将 Singular 的值更改为 true 能解决 (-1,0) 和 (1,0) 处的几何奇异点。较大的着色区域可能需要细化,但问题不大。
Q, = quad2d(fun, -1, 1, c, d; AbsTol=1e-8, FailurePlot=true, Singular=true)
┌ Warning: Reached the maximum number of function evaluations (2000). The result passes the global error test.
从这里,您可以利用对称性:
Q = 4 * quad2d(fun, 0, 1, 0, d; AbsTol=1e-8, Singular=true, FailurePlot=true)[1]
Q = 0.9817476400794798
但代码在奇异点附近仍然应对困难。它可能无法提供更高的准确度:
Q = 4*quad2d(fun,0,1,0,d,AbsTol=1e-10,Singular=true,FailurePlot=true)[1]
┌ Warning: Reached the maximum number of function evaluations (2000). The result passes the global error test.
在更高准确度时,更改坐标可能更为奏效。
polarfun = (theta, r) -> fun(r*cos(theta), r*sin(theta))*r
Q = 4*quad2d(polarfun, 0, pi/2, 0, 1; AbsTol=1e-10)[1]
Q = 0.98174770424664
最好是通过将积分区域分为两个部分来将奇异点放在边界:
Q1 = 4*quad2d(polarfun, 0, pi/2, 0, 0.5; AbsTol=5e-11)[1]
Q2 = 4*quad2d(polarfun, 0, pi/2, 0.5, 1; AbsTol=5e-11)[1]
Q = Q1 + Q2
Q = 0.9817477042468746
# 输入参数
fun - 求积分的函数函数句柄
求积分的函数,指定为函数句柄。函数 Z = fun(X,Y) 必须接受标量 X 和 Y,并返回相应标量值 Z。
示例: (x, y)-> x^2 - y^2
数据类型: 句柄函数
a,b - 积分的 x 范围标量
积分的 x 范围,指定为标量。
数据类型: Int | Float
复数支持: 是
c,d - 积分的 y 范围标量 | 函数句柄
积分的 y 范围,指定为标量或函数句柄。每个范围都可以指定为标量或函数句柄。如果将范围指定为函数句柄,则它们是积分 ymin(x) = x->c(x) 和 ymax(x) =x->d(x) 的 x 范围的函数。函数句柄 ymin 和 ymax 必须接受并返回标量值。
数据类型: Int | Float | 函数句柄
复数支持: 是
# 关键字参数
指定可选的、以逗号分隔的 Key=Value 键值参数。Key 为参数名称,Value 为对应的值。写法为Key = Value。您可采用任意顺序指定多个键值对组参数,如 Key1=Value1,...,KeyN=ValueN 所示。
示例: quad2d((x,y) -> x*y^2, 0, 1, 0, 2; AbsTol = 1e-3) 将积分的绝对容差指定为 1e-3。
AbsTol - 绝对误差容限1e-5(默认) | 标量
绝对误差容限,指定为由 AbsTol 和非负实数组成的键值对组。
quad2d 尝试满足ERRBND <= max(atol,rtol*|Q|)。如果 |Q| 足够小,这是绝对误差控制;如果 |Q| 较大,这是相对误差控制。未指定容差时使用默认容差。AbsTol 的默认值是 1e-5。RelTol 的默认值是 100*eps()。这还是 RelTol 的最小值。较小的 RelTol 值会自动增大为默认值。
数据类型: Int | Float
RelTol - 相对误差容限100*eps()(默认) | 标量
相对误差容限,指定为由 RelTol 和非负实数组成的键值对组。
quad2d 尝试满足ERRBND <= max(atol,rtol*|Q|)。如果 |Q| 足够小,这是绝对误差控制;如果 |Q| 较大,这是相对误差控制。未指定容差时使用默认容差。AbsTol 的默认值是 1e-5。RelTol 的默认值是 100*eps()。这还是 RelTol 的最小值。较小的 RelTol 值会自动增大为默认值。
数据类型: Float
MaxFunEvals - fun 的最大计算次数2000(默认) | 标量
fun 的最大计算次数,指定为由 MaxFunEvals 和一个标量组成的键值对组。使用此选项可限制 quad2d 计算函数 fun 的次数。
数据类型: Int
FailurePlot - 开启或关闭生成失败绘图false(默认) | true
开启或关闭生成失败绘图,指定为由 FailurePlot 和一个逻辑值组成的键值队组。将 FailurePlot 设置为 true,可在达到 MaxFunEvals 时生成需要进一步细化的区域的图形表示形式。如果在达到 MaxFunEvals 前成功计算积分,则不生成绘图。失败绘图包含(通常是)在内部映射到矩形的 4 条边的区域。聚集的小区域表示积分中的问题区域。
数据类型: Bool
Singular - 开启或关闭转换边界奇异点true(默认) | false
开启或关闭转换边界奇异点,指定为以由 Singular 和一个逻辑值组成的键值对组。默认情况下,quad2d 将利用转换削弱边界奇异点以提高性能。将 Singular 设置为 false 可禁用这些转换,为一些平滑问题提供性能优势。
数据类型: Bool
# 输出参数
q - 计算的积分标量
计算的积分,以标量形式返回。
数据类型: Int | Float
E - 误差边界标量
误差边界,以标量形式返回。误差边界提供计算的积分 q 与积分 I 的精确值之间误差的上界,即
数据类型: Int | Float
# 参考文献
[1] L.F. Shampine, "MATLAB Program for Quadrature in 2D." Applied Mathematics and Computation. Vol. 202, Issue 1, 2008, pp. 266–274.