我想要查看一个函数 intensity(r) 在空间上的图像，因此 r 是具有辐射对称性的。
然而，我从 intensity(r) = integrate(integrand(r), (x,0,5)) 推导出我的强度函数，其中 integrand = exp(-x**2) * exp(np.pi*1j*(-x)) * besselj(0, r*x) * x。
上述所有语法都使用了 sympy 包，所以我首先定义了 x，y = symbols（'x r'）。
我使用符号变量，因为这样似乎更容易进行可视化处理，并将 r 作为变量保留到最后，在绘制图形并分配数值之前。
然而，使用符号变量进行积分计算似乎非常耗时间。
是否有任何方法可以使用符号变量进行数值积分？
唯一的替代方案是预先定义 r 的值并为每个值找到积分吗？
顺便提一句：当你创建一个符号表达式时，请将其保持为符号。不要混合使用实际的浮点数 np.pi 和复杂的浮点数 1j，而是使用 SymPy 的 pi 和 I。
from sympy import exp, pi, I, besselj, symbols
x, r = symbols('x r')
integrand = exp(-x**2) * exp(pi*I*(-x)) * besselj(0, r*x) * x
但是，看起来 SymPy 似乎无法对贝塞尔函数与exp(-x**2) * exp(pi*I*(-x))的乘积进行积分。这已经在将 r 替换为 1 时发生了，因此 r 的符号性质并不重要。
直接回答您的问题：
有没有一种方法可以在符号变量上执行数值积分？
没有，就像没有干水一样，这是个自相矛盾的说法。
唯一的其他选择是预先定义 r 的值并找到每个值的积分吗？
是的。可以通过 SymPy（它会调用 mpmath）完成：
>>> intensity = lambda r_: Integral(integrand.subs(r, r_), (x, 0, 5)).evalf()
>>> intensity(3)
0.0783849036516177 - 0.125648626220306*I
鉴于这个函数是复数值的，你打算如何绘制它还不太清楚。也许你想绘制强度的绝对值？
无论如何，使用 SymPy/mpmath（纯 Python）进行积分过于缓慢，不适合绘图。最好使用 SciPy 的 quad 进行积分。它不能处理复合积分，所以我将实部和虚部分别积分。
from scipy.integrate import quad
from scipy.special import jn
integrand = lambda x, r: np.exp(-x**2) * np.exp(np.pi*1j*(-x)) * jn(0, r*x) * x
intensity = lambda r: np.sqrt(quad(lambda x: np.real(integrand(x, r)), 0, 5)[0]**2 + quad(lambda x: np.imag(integrand(x, r)), 0, 5)[0]**2)
intensity(3)的计算速度比之前的版本要快得多。我们可以画出它的图像：import matplotlib.pyplot as plt
t = np.linspace(0, 3)
plt.plot(t, np.vectorize(intensity)(t))
