我正在尝试用sympy进行一些符号矩阵计算。我的目标是获得一些矩阵计算结果的符号表示。我遇到过一些问题,我把它们归结为这个简单的例子,在这个例子中,我试图计算对一个指定矩阵取幂的结果,并将其乘以一个任意向量。
>>> import sympy
>>> v = sympy.MatrixSymbol('v', 2, 1)
>>> Z = sympy.zeros(2, 2) # create 2x2 zero matrix
>>> I = sympy.exp(Z) # exponentiate zero matrix to get identity matrix
>>> I * v
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "sympy/matrices/matrices.py", line 507, in __mul__
blst = B.T.tolist()
AttributeError: 'Transpose' object has no attribute 'tolist'
相反,如果我直接创建单位矩阵,然后乘以v,那么就没有问题了:
>>> I_ = sympy.eye(2) # directly create the identity matrix
>>> I_ == I # check the two matrices are equal
True
>>> I_ * v
v
我注意到的一件事是这两个单位矩阵属于不同的类:
>>> I.__class__
sympy.matrices.immutable.ImmutableMatrix
>>> I_.__class__
sympy.matrices.dense.MutableDenseMatrix
我还发现调用as_mutable()
方法提供了一个解决方案。
>>> I.as_mutable() * v
v
总是有必要把as_mutable()
调用在一个人的线性代数计算?我猜不是,相反,这些错误表明我使用了错误的策略来解决我的问题,但我不知道正确的策略是什么。有人有什么建议吗?
我已经阅读了关于不可变矩阵的文档页面,但我仍然可以使用一些帮助来理解它们与标准可变矩阵的差异在这里是重要的,以及为什么一些操作(例如sympy.exp)在这些不同的类之间转换。
我认为这是Sympy的一个bug:
在Python中,可以从两边重载乘法操作符。A*B
可以在内部通过调用A.__mul__(B)
或B.__rmul__(A)
来处理。Python首先调用A.__mul__
,如果此方法不存在或返回NotImplemented
,则Python自动尝试B.__rmul__
。SymPy使用一个名为call_highest_priority的装饰器来决定使用哪一个实现。它查找所涉及类的_op_priority
,并以更高的优先级调用实现的函数。在您的示例中,v
和I
的优先级为11,I_
的优先级为10.01,因此首选I
。而且,I
使用的__mul__
的基本实现缺少修饰符。
长话短说,I*v
最终总是调用I.__mul__
, __mul__
不能处理MatrixSymbol
s,但也不返回NotImplemented
。v.__rmul__(I)
工作正常
正确的修复方法是在matrices.py
中捕获AttributeError
并返回NotImplemented
,即
try:
blst = B.T.tolist()
except AttributeError:
return NotImplemented
Python会自动回退到__rmul__
。hack'ish修复将是调整_op_priority
。无论哪种方式,您都应该提交一个错误报告:如果错误是由设计引起的(也就是说,如果您不小心尝试了一些不应该工作的东西),那么错误消息将显示。