这太有趣了。我可以说的一件事是,这个问题是特定于np.matrix
子类。例如,以下工作正常:
h = np.array(hamiltonian)
unitary = [linalg.expm(-(1j)*t*h) for t in t_list]
更深入地研究回溯,异常是在_fragment_2_1
in scipy.sparse.linalg.matfuncs.py
, 具体来说这些行 https://github.com/scipy/scipy/blob/81c096001974f0b5efe29ec83b54f725cc681540/scipy/sparse/linalg/matfuncs.py#L748-L755:
n = X.shape[0]
diag_T = T.diagonal().copy()
# Replace diag(X) by exp(2^-s diag(T)).
scale = 2 ** -s
exp_diag = np.exp(scale * diag_T)
for k in range(n):
X[k, k] = exp_diag[k]
错误信息
X[k, k] = exp_diag[k]
TypeError: only length-1 arrays can be converted to Python scalars
向我建议exp_diag[k]
应该是一个标量,但是返回一个向量(并且你不能将向量分配给X[k, k]
,这是一个标量)。
设置断点并检查这些变量的形状证实了这一点:
ipdb> l
751 # Replace diag(X) by exp(2^-s diag(T)).
752 scale = 2 ** -s
753 exp_diag = np.exp(scale * diag_T)
754 for k in range(n):
755 import ipdb; ipdb.set_trace() # breakpoint e86ebbd4 //
--> 756 X[k, k] = exp_diag[k]
757
758 for i in range(s-1, -1, -1):
759 X = X.dot(X)
760
761 # Replace diag(X) by exp(2^-i diag(T)).
ipdb> exp_diag.shape
(1, 4)
ipdb> exp_diag[k].shape
(1, 4)
ipdb> X[k, k].shape
()
根本问题是exp_diag
假设为一维向量或列向量,但对角线np.matrix
对象是行向量。这突出了一个更普遍的观点:np.matrix
通常没有得到很好的支持np.ndarray
,所以大多数情况下最好使用后者。
一种可能的解决方案是使用np.ravel()
压扁diag_T
转化为一维np.ndarray
:
diag_T = np.ravel(T.diagonal().copy())
这似乎解决了您遇到的问题,尽管可能还有其他相关问题np.matrix
我还没有发现。
我已打开拉取请求here https://github.com/scipy/scipy/pull/5867.