大范数正半定矩阵的取幂和特征值计算中的差异

Question

我正在用埃尔米特（复）正半定矩阵进行一些计算。计算涉及取它们的指数。矩阵

A

的指数（以 2 为底），我将其表示为

2^A

，可以通过对特征值应用幂函数来获得。

因此

2^A

的每个特征值必须是2**（

A

的特征值）。然而，当我在 numpy/scipy 中实现这一点时，具有大范数（迹范数/核范数 >= 60）的矩阵似乎存在差异。例如，这两组特征值不一致，并且这两个特征值向量之间的欧几里德距离可能非常大，并且随着

A

变大（以迹范数测量）而增长。

这是我的代码

import numpy as np
import scipy as sp
import scipy.linalg as sla

d = 2
X = np.random.randn(d, d) + 1j*(np.random.randn(d, d))  # Generate an arbitrary complex matrix
P = [email protected]().T                         # Multiply by its transpose to obtain PSD matrix
P = P/np.trace(P)                        # Normalise to obtain density matrix (PSD matrix with unit trace)

for t in range(20, 75, 5):
  Q = t*P                                 # Scaling the matrix
  D, V = sla.eigh(Q)                      # Eigendecomposition 
  M = V @ np.diag(2**D) @ V.conj().T      # Constructing M = 2**Q

  #Comparing the eigenvalues of M vs 2**(EigVals(Q)). This should be equal.
  print(t, norm(sla.eigvalsh(M) - 2**sla.eigvalsh(Q))) 

这是

X

随机初始化的输出。 第一列是

t

，它是迹范数，第二列是 使用这两种方法的特征值向量。

20 3.654601116809922e-12
25 6.374134689249895e-11
30 9.33093987260751e-10
35 3.728611604386753e-09
40 5.962181140973812e-08
45 9.573880441485143e-07
50 1.4168248818502024e-05
55 4.866518429480493e-05
60 0.0024326798970430987
65 0.025923866144586645
70 0.12695862169940617

差异随着 t 的增加而增大。可能是什么原因？我的实现有错误吗？

Answer 1

您的实施是正确的。然而，当比较具有大范数的数组时，您应该考虑对结果进行标准化。不可能比

machine_precision * norm(a)

更好

如果将最后一行替换为

print(t, lg.norm(lg.eigvalsh(M) - 2**lg.eigvalsh(Q)) / lg.norm(M))

，您会发现结果在机器精度范围内保持正确