在 Sympy 中使用lambdify 将具有复指数的表达式转换为函数

Question

我有一个 Sympy 形式的表达式

exp =  0.5*y1**2 + 2*y1*y2 + 2.5*x5**2 + 5.0*x6**2 - 25.0*exp(2.0*I*x4)*exp(1.0*I*x5) - 25.0*exp(2.0*I*x4) + 100.0 - 25.0*exp(-2.0*I*x4) - 25.0*exp(-2.0*I*x4)*exp(-1.0*I*x5)

我想将上述表达式转换为使用变量作为参数的函数。此外，参数将稀疏 scipy csc 矩阵作为输入代替变量。所以，我用下面的把它转换成一个函数。

func = lambdify((y1,y2,x4,x5,x6),exp,"scipy")

虽然，此函数生成时没有任何错误，但该函数不适用于 csc 稀疏矩阵。它给出了以下错误：

from scipy import sparse
a = sparse.identity(12000,format="csc")
func(a,a,a,a,a) # just to test the output
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
~/.local/lib/python3.9/site-packages/scipy/sparse/base.py in __getattr__(self, attr)
    686         else:
--> 687             raise AttributeError(attr + " not found")
    688 

AttributeError: exp not found

The above exception was the direct cause of the following exception:

TypeError                                 Traceback (most recent call last)
<ipython-input-321-ee4555e1ca2e> in <module>
      1 a = sparse.identity(21,format="csc")
----> 2 func(a,a)

<lambdifygenerated-23> in _lambdifygenerated(x2, n2)
      1 def _lambdifygenerated(x2, n2):
----> 2     return (0.125*n2**2 - 25.0*exp(2.0*1j*x2) - 25.0*exp(-2.0*1j*x2))

~/.local/lib/python3.9/site-packages/scipy/_lib/deprecation.py in call(*args, **kwargs)
     18             warnings.warn(msg, category=DeprecationWarning,
     19                           stacklevel=stacklevel)
---> 20             return fun(*args, **kwargs)
     21         call.__doc__ = msg
     22         return call

TypeError: loop of ufunc does not support argument 0 of type csc_matrix which has no callable exp method

这表明稀疏矩阵的指数未定义。因此，我尝试将 exp(I*x1) 替换为单独的变量，然后将其用作参数。但是，当我使用 subs 时，我得到以下结果。

$exp.subs(exp(I*x4),e4)
0.5*y1**2 + 2*y1*y2 + 2.5*x5**2 + 5.0*x6**2 - 25.0*e4**2.0*exp(1.0*I*x5) - 25.0*e4**2.0 + 100.0 - 25.0/e4**2.0 - 25.0/e4**2.0*exp(-1.0*I*x5)

上述结果的问题是指数的负幂被分母中的变量e4代替了。因此，相同的lambdify版本不适用于矩阵，因为矩阵的除法运算是未定义的。我尝试使用 replace 分别替换exp(1.0*I*x5) 和exp(-1.0*I*x5) 形式的术语，但这不适用于exp(2.0*I*x5) 等术语。

有什么办法可以解决这种情况吗？任何帮助将不胜感激。

编辑：修改上述答案的某些部分以包含可重现的示例。另外，我使用的是复数方阵，其大小可以大得多，例如 12000 x 12000。

谢谢！

Answer 1

np.exp 无法处理稀疏矩阵：

In [64]: from scipy import sparse
In [65]: M = sparse.random(10,10,.2, 'csc')
In [66]: np.exp(M)
Traceback (most recent call last):
  File "/usr/local/lib/python3.8/dist-packages/scipy/sparse/base.py", line 687, in __getattr__
    raise AttributeError(attr + " not found")
AttributeError: exp not found

The above exception was the direct cause of the following exception:
Traceback (most recent call last):
  File "<ipython-input-66-1b6c717c2aff>", line 1, in <module>
    np.exp(M)
TypeError: loop of ufunc does not support argument 0 of type csc_matrix which has no callable exp method

np.exp 和许多numpy 函数一样，首先将非数组参数转换为数组：

In [67]: np.array(M)
Out[67]: 
array(<10x10 sparse matrix of type '<class 'numpy.float64'>'
    with 20 stored elements in Compressed Sparse Column format>, dtype=object)

但是这种简单的转换只会产生一个 0d 对象 dtpe 数组，而不会改变稀疏矩阵。反过来，它会尝试使用不存在的exp 方法。

将稀疏矩阵正确转换为（密集）数组确实有效：

np.exp(M.A)

产生一个数组结果。切换到密集数组不会丢失任何东西，因为 exp(0) 是 1。

使用help(f)（或ipython f?）查看函数的文档（和代码），例如：

Created with lambdify. Signature:

func(x, y)

Expression:

0.5*x**2 + 25.0*exp(2*I*x)

Source code:

def _lambdifygenerated(x, y):
    return (0.5*x**2 + 25.0*exp(2*1j*x))

np.expm1(M) 也有效，因为 M.expm1() 已定义。

由于 exp(0) 为 1，在稀疏矩阵上执行 exp 将产生一个数组，其中原始数组为 0，即非稀疏数组。另一方面，expm1 保留了稀疏性，实际上它可以将函数应用于所有非零元素，忽略（默认）0。

总而言之，将您的 x4 和 x5 更改为密集数组，它应该可以运行。

Answer 2

我能够通过编辑 lambdify 函数使用的全局字典来使其工作。此函数描述了替换 sympy 表达式中使用的数学函数的函数。由于这里的问题在于 exp 的运算符，因此可以使用以下命令从 scipy 将其替换为 expm。结果函数现在可以将 csc 稀疏矩阵作为输入。

f = lambdify((symbols("v2,x2")),exp.expand(), [{'exp':scipy.linalg.expm},"scipy"])

希望这可以帮助其他面临类似问题的人！