python - python中的数值积分与符号积分

Question

我想查看函数 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 的值并找到每个值的积分？

Answer 1

旁白:当您创建符号表达式时，请保持其符号性。不要混合使用真正的 float np.pi 和复杂的 float 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 的符号性质无关紧要。

直接回答您的问题:

Is there any way of performing numerical integration with symbolic variables?

没有，就像没有干涸的水。这是一种自相矛盾的说法。

Is the only other alternative defining the values of r a priori and finding the integral for each one of them?

是的。可以通过 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))