关于以下问题,我有点困惑:
我正在计算量子通道的固定点,这意味着我想计算特定矩阵的前导特征向量。矩阵使得其维数为n ^ 2 x n ^ 2并且以这样的方式定义:将具有形状n x n的矩阵重新形成的前导特征值是正矩阵(具有正特征值的自伴随)。
如果我用scipy.sparse.linalg.eigs
这样做,但我得到错误的结果。然而,精确的计算(使用scipy.linalg.eig
)工作正常。我尝试使用k
和ncv
为解算器的参数,但没有让我正常工作,除非我设置k=n**2
在这种情况下eigs
只是指eig
。但是,对于我实际上考虑到通道(下面脚本中的super_op
)实际编码为LinearOperator
的情况,这不起作用。所以我依靠使用eigs
:/
任何人都知道如何正确运行?
感谢大家提前!
import numpy as np
from numpy.random import rand
from numpy import tensordot as td
from scipy.sparse.linalg import eigs
from scipy.linalg import eig
n = 16
d = 3
kraus_op = .5 - rand(n, d, n) + 1j * (.5 - rand(n, d, n))
super_op = td(kraus_op, kraus_op.conj(), [[1], [1]]).transpose(0, 2, 1, 3)
########
# Sparse
########
vals, vecs = eigs(super_op.reshape(n**2, n**2), k=n*(n-1), which='LM')
rho = vecs[:,0].reshape(n, n)
print('is self adjoint: ', np.allclose(rho, rho.conj().T))
super_op_times_rho = td(super_op, rho, [[2, 3], [0, 1]])
print('super_op(rho) == lambda * rho :', np.allclose(rho, super_op_times_rho/vals[0]))
########
# Exact
########
vals, vecs = eig(super_op.reshape(n**2, n**2))
rho = vecs[:,0].reshape(n, n)
print('is self adjoint: ', np.allclose(rho, rho.conj().T))
super_op_times_rho = td(super_op, rho, [[2, 3], [0, 1]])
print('super_op(rho) == lambda * rho :', np.allclose(rho, super_op_times_rho/vals[0]))
结果是:
is self adjoint: False
super_op(rho) == lambda * rho : True
is self adjoint: True
super_op(rho) == lambda * rho : True
为了完整性:
Python 3.5.2 numpy 1.16.1 scipy 1.2.1
毕竟我在同事的帮助下找到了解决方案:
虽然eig
给出了特征向量(1)(按幅度排序)和(2)使得第一个条目是实数,eigs
有另一个排序,不能像在eig
中那样,也不能规范特征向量的复杂相位。通过将张量除以第一个条目的相位(对应于确保第一个对角线元素是真实的,以消除选择特征向量的复杂相位并使Hermiticity成为可能的自由......),可以轻松地校正相位。 。
因此,稀疏案例的更正代码段是:
vals, vecs = eigs(super_op.reshape(chi**2, chi**2), k=chi*(chi-1), which='LM')
# find the index corresponding to the largest eigenvalue
arg = np.argmax(np.abs(vals))
rho = vecs[:,arg].reshape(chi, chi)
# regularize the output array
rho *= np.abs(rho[0, 0])/rho[0, 0]