四联体中具有非常数极限的二重积分？

Question

我需要集成一个在 2Dim 中定义的函数（

z

和半径

r

），但我没有表达式。

我可以在任意位置查询函数

(z,r)

并获取返回值。

我有跨 z 的积分范围：

[-z_range, z_range]

，我将其划分为

N_z

点为：

z_is = -z_range + (np.arange(N_z) + 0.5) * (2.*z_range/N_z)

对于

z_is

中的每个值，跨

r

积分的范围为：

[0, r_thresh_at_this_z]

。

r_thresh_at_this_z

是从称为

z

的

z_i

的值获得的：

def get_r_thresh(z_i):
    return expression_of_z_i # returns a positive float.

因此径向积分的范围取决于

z

上的值。

我的函数 f 为：

def f(r,z):
    return interpolator_for_f(r,z)

我想以最有效的方式使用

quadpy

包，因为它被创建为能够以这种方式使用。

我正在考虑使用 for 循环来遍历

z_is

并针对

[0, r_thres_for_that_z]

中的每个值在

z_is

上执行高斯求积。

我可以使用：

results = np.zeros((N_z))
errs = np.zeros((N_z))
for i in range(N_z):
   def f(r):
      return interpolator_for_f(r,z_is[i])
   results[i], errs[i] = quadpy.quad(f, [0, r_thresh_at_this_z])

但我觉得for循环并不是使用quadpy最有效的方法。

你能告诉我在仅使用快速 numpy 数组进行积分时缺少什么吗，所以没有 for 循环？

[我已经阅读了函数

x

的输入

f

的形状。 就我而言

d = 1

因为它是 r 上的线积分。

n = N_z

因为我想执行

N_z

这样的线积分，然后我将其相加以获得 1 个单标量，即（整个）二重积分的结果。

p = 1000

因为假设对于

r

的每个值，我希望在

z

上有 1000 个积分点。

所以我需要在

N_z * 1000

点对函数进行采样。

函数

f

应返回一个数组形状

(N_z, 1000)

这些参数的识别有什么帮助吗？]

谢谢！

Answer 1

Quadpy 作者在这里。

我想以最有效的方式使用quadpy包，因为它被创建为能够以这种方式使用。

事实并非如此。 Quadpy 有助于集成各种规范域中的函数，但尚不具备在任意 2D 域上集成的功能。 doublequad/nquad 应该在quadpy 中，但现在还没有。