我需要集成一个在 2Dims （z和 radius 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 循环？

[我已经阅读了函数 f 的输入 x 的形状。在我的情况下 d = 1 因为它是跨 r 的线积分。n = N_z 因为我想执行 N_z 这样的线积分，然后我将它们相加以获得 1 个单标量，即（整个）双积分的结果。

p = 1000 因为假设我想要 r 上的 1000 个积分点，对于 z 的每个值。

所以我需要在 N_z * 1000 点对函数进行采样。

函数 f 应返回一个形状为 (N_z, 1000) 的数组

这些参数的识别是否有任何帮助？]

谢谢！