我想评估以下积分：

$$\int_{0}^{60} \ \left(\int_{0}^{2z} 0.5\cdot t \left(\mathrm{erf}(t-a) -1 \right)J_{0}(qt)\mathrm{d}t \right)^2 \mathrm{exp}\left(-\frac{(z-a)^2}{2s^2}\right)\mathrm{d}z$$

其中 a、q 和 s 是常数。如果没有第一个积分之外的正方形，可以使用dblquadfrom轻松评估scipy：

import math
from scipy.integrate import dblquad
from scipy.special import erf, jv


def h(t, z):
    return f(t) * g(z)


def f(t):
    return 0.5 * t * (erf(t - a) - 1) * jv(0, q * t)


def g(z):
    return math.exp(-((z - a) ** 2)/(2 * (s ** 2)))


if __name__ == '__main__':
    a, q, s = 0, 1, 2  # set the constants
    result, abserr = dblquad(h, 0, 60, lambda z: 0, lambda z: 2 * z)
    print(f'result: {result}, abserr: {abserr}')

我查看了积分的平方，如图所示，积分的平方等价于二重积分：

$$\left(\int_a^bf(x)\text{d}x\right)^2 = \int_a^b \int_a^b f(x) f(y) \text{d}x\text{d}y$$

这使问题变得更加复杂，因为现在我必须评估三重积分，这会使评估速度变慢。有没有更好的方法来解决这个问题？任何帮助将不胜感激。