比较 Mathematica 和 Python 中的卷积

Question

我将 Python (using sympy's symbolic variables) 和 Mathematica 中的卷积结果与其Convolve 函数进行比较。

在 Python 中，我的 MWE 是

from numpy import linspace, pi
from numpy.random import randn
from scipy.signal import fftconvolve
import matplotlib.pyplot as plt
from sympy import symbols
from sympy.utilities.lambdify import lambdify

a = 0.43
b = 0.41
c = 0.65
d = 0.71

x = symbols('x')
f = 2*b / ((x-a)**2 + b**2)
g = 2*d / ((x-c)**2 + d**2)
fog = fftconvolve(f,g,mode='same')
fog_fun = lambdify(x,fog,'numpy') # returns a numpy-ready function
x = linspace(-20,20,int(1e3))
dx = x[1]-x[0]
fogS = fog_fun(x)

fogA = 4*pi*(b+d)/((x-a-c)**2+(b+d)**2) # correct analytic solution

plt.figure()
plt.plot(x,fogA,lw=2,label='analytic')
plt.plot(x,fogS,lw=2,label='sympy')
plt.grid()
plt.legend(loc='best')
plt.show()

使用符号变量x 计算卷积。结果函数（在lambdifying之前）是

fog = 1.1644/(((x - 0.65)**2 + 0.5041)*((x - 0.43)**2 + 0.1681))

analytic（fogA，Mathematica）和sympy（fogS，Python）之间没有一致：

我的 Mathematica 代码是：

a = 0.43; b = 0.41; c = 0.65; d = 0.71;
fogA = FullSimplify[Convolve[2*b/((t-a)^2+b^2),2*d/((t-c)^2+d^2), t, x]];
fogS = 1.1644/(((x - 0.65)^2 + 0.5041)*((x - 0.43)^2 + 0.1681));

在哪里

fogA = (17.683+x*(-30.4006+14.0743*x))/(3.04149+x*(-7.9428+x*(8.3428+x*(-4.32+1.*x))))

fogS 和 fogA 的图形与 Python 相同。

为什么analytic 和sympy 解决方案之间存在如此大的分歧？我怀疑问题出在sympy 上。另一种 Pythonic 方法是对两个 数组 进行卷积，这似乎与 analytic 解决方案一致。

f = 2*b / ((x-a)**2 + b**2)
g = 2*d / ((x-c)**2 + d**2)
fogN = fftconvolve(f,g,mode='same')*dx # numeric

（注意：这是一个 MWE。我想要卷积的实际 f 和 g 比本文中定义的洛伦兹要复杂得多。）

Answer 1

我不认为这是使用scipy + sympy 的合理方式。 我真的很惊讶你从lambdify 得到结果。

你应该做的是使用卷积的符号定义，而不是使用scipy.signal.fftconvolve()，例如：

from sympy import oo, Symbol, integrate

def convolve(f, g, t, lower=-oo, upper=oo):
    tau = Symbol('__very_unlikely_name__', real=True)
    return integrate(
        f.subs(t, tau) * g.subs(t, t - tau), (tau, lower, upper))

改编自here。