Sympy与Scipy结合出现积分异常，求错误原因及修复方法

编写计算导数与积分的代码时，使用sympy和scipy.integrate工具出现异常：前两步导数计算正确，但积分结果错误（正确值应为18，实际得到约8.975），且必须在积分函数中使用z[0]否则报错。复现代码如下：

from scipy.integrate import quadrature
import sympy
import math

def func(z):
    return math.e**(2*z**2)+3*z

z = symbols('z')
f = func(z)
dlogf_dz = sympy.Lambda(z, sympy.log(f).diff(z))
print(dlogf_dz)
print(dlogf_dz(10))

def integral(z):
    # 为什么必须用z[0]？
    d_ARF=dlogf_dz(z[0])
    return d_ARF

result = quadrature(integral, 0, 3)
print(result)

输出：

Lambda(z, (4.0*2.71828182845905**(2*z**2)*z + 3)/(2.71828182845905**(2*z**2) + 3*z))
40.0000000000000
(8.97457203290041, 0.00103422711649337)

但将导数表达式硬编码到积分函数中时，结果正确：

def integral(z):
    d_ARF=(4.0*2.71828182845905**(2*z**2)*z + 3)/(2.71828182845905**(2*z**2) + 3*z)
    return d_ARF

result = quadrature(integral, 0, 3)
print(result)

输出：

(18.000000063540558, 1.9408245677254854e-07)

错误原因 #

1. 符号运算与数值运算混淆：使用math.e（数值型常量）而非sympy.E（符号型常量），导致func(z)传入sympy符号z时，生成的是数值近似的表达式而非纯符号表达式，后续sympy生成的Lambda函数无法正确支持向量化计算。
2. sympy.Lambda不兼容scipy的向量化输入：scipy的quadrature会传入数组（而非单个标量）进行向量化采样计算，但sympy.Lambda生成的函数仅支持标量输入，因此必须用z[0]取单个元素，这会导致函数无法处理批量采样点，积分计算时采样错误，结果偏差。

修复步骤 #

1. 替换math.e为sympy.E，确保符号运算的纯粹性。
2. 使用sympy.lambdify替代sympy.Lambda，生成兼容numpy数组的向量化函数，适配scipy积分工具的输入要求。

修复后的完整代码 #

from scipy.integrate import quadrature
import sympy
import numpy as np

def func(z):
    return sympy.E**(2*z**2) + 3*z

# 定义符号变量
z_sym = sympy.symbols('z')
f_sym = func(z_sym)
# 计算对数导数的符号表达式
log_f_deriv = sympy.log(f_sym).diff(z_sym)
# 转换为支持numpy数组的函数
dlogf_dz = sympy.lambdify(z_sym, log_f_deriv, 'numpy')

# 测试导数计算
print(dlogf_dz(10))

# 积分函数直接使用dlogf_dz，无需z[0]
def integral(z):
    return dlogf_dz(z)

result = quadrature(integral, 0, 3)
print(result)

输出：

40.0
(18.000000063540558, 1.9408245677254854e-07)

额外说明 #

• sympy.lambdify可以将符号表达式转换为numpy兼容的函数，自动支持向量化计算，因此scipy传入数组时无需手动取单个元素。
• 使用sympy.E保证符号运算的精确性，避免数值近似导致的表达式变形，确保生成的导数函数与硬编码的数学表达式完全一致。

内容的提问来源于stack exchange，提问作者NeStack