为什么对称矩阵的 Numpy 特征向量不能构造原始矩阵

Question

对于一个对称实矩阵A，可以分解为A=Q'UQ，其中Q是特征向量， U是特征值矩阵，Q'是Q的转置矩阵。但是，当我使用numpy.linalg.eig() 计算特征值和特征向量时， 在某些情况下，结果是正确的，而在另一些情况下，结果是错误的。例如：

1. A = [[3, -1, -1, -1], [-1, 3, -1, -1], [-1, -1, 3, -1], [-1, -1, -1, 3]]
2. A = [[1, 0, -1, 0], [0, 1, -1, 0], [-1, -1, 3, -1], [0, 0, -1, 1]]

在Case1中，原矩阵A可以重构成功，但在Case2中重构失败。对于第二种情况矩阵，我手动计算特征值和特征向量。结果是对的，如下图。我真的很想知道为什么？！

实验代码如下：

import numpy as np
import scipy.linalg as spl

N = 4
# case 1
# A = np.array([[3, -1, -1, -1], [-1, 3, -1, -1], [-1, -1, 3, -1], [-1, -1, -1, 3]])
# case 2
A = np.array([[1, 0, -1, 0], [0, 1, -1, 0], [-1, -1, 3, -1], [0, 0, -1, 1]])
lam, vec = np.linalg.eig(A)

# calculate the orthonormal eigenvectors matrix Q
vec = spl.orth(vec)

# orthonormal eigenvectors matrix Q calculated by hand in case 2
# vec = np.array([[np.sqrt(12)/12, np.sqrt(12)/12, -3*np.sqrt(12)/12, np.sqrt(12)/12], [np.sqrt(4)/4, np.sqrt(4)/4, np.sqrt(4)/4, np.sqrt(4)/4], [-np.sqrt(2)/2, np.sqrt(2)/2, 0, 0], [-np.sqrt(6)/6, -np.sqrt(6)/6, 0, 2*np.sqrt(6)/6]]).T

# calculate eigenvalues matrix U
lam_matrix = np.zeros((N,N))
i_0 = [i for i in range(N)]
j_0 = [i for i in range(N)]
lam_matrix[i_0, j_0] = lam

# print the experimental result
print('#### Result ####')
print('eigenvalues')
print(lam)
print('eigenvectors')
print(vec)
print('orthogonality of eigenvectors')
print(vec.T.dot(vec))
print('reconstruct the orginal matix')
print(vec.dot(lam_matrix).dot(vec.T))

Answer 1

scipy.linalg.orth 使用 SVD 为输入范围构造一个正交基，它并不总是承诺返回矩阵 A 的正交特征向量。

要计算正交特征值分解，请改用eigh。

import numpy as np
import scipy.linalg as spl

N = 4
# case 1
# A = np.array([[3, -1, -1, -1], [-1, 3, -1, -1], [-1, -1, 3, -1], [-1, -1, -1, 3]])
# case 2
A = np.array([[1, 0, -1, 0], [0, 1, -1, 0], [-1, -1, 3, -1], [0, 0, -1, 1]])
lam, vec = np.linalg.eigh(A)

# calculate the orthonormal eigenvectors matrix Q
#vec = spl.orth(vec)

# orthonormal eigenvectors matrix Q calculated by hand in case 2
# vec = np.array([[np.sqrt(12)/12, np.sqrt(12)/12, -3*np.sqrt(12)/12, np.sqrt(12)/12], [np.sqrt(4)/4, np.sqrt(4)/4, np.sqrt(4)/4, np.sqrt(4)/4], [-np.sqrt(2)/2, np.sqrt(2)/2, 0, 0], [-np.sqrt(6)/6, -np.sqrt(6)/6, 0, 2*np.sqrt(6)/6]]).T

# calculate eigenvalues matrix U
lam_matrix = np.zeros((N,N))
i_0 = [i for i in range(N)]
j_0 = [i for i in range(N)]
lam_matrix[i_0, j_0] = lam


# print the experimental result
print('#### Result ####')
print('eigenvalues')
print(lam)
print('eigenvectors')
print(vec)
print('orthogonality of eigenvectors')
print(vec.T.dot(vec))
print('reconstruct the orginal matix')
print(vec.dot(lam_matrix).dot(vec.T))

返回

#### Result ####
eigenvalues
[-2.29037709e-16  1.00000000e+00  1.00000000e+00  4.00000000e+00]
eigenvectors
[[-5.00000000e-01  2.26548862e-01 -7.84437556e-01  2.88675135e-01]
 [-5.00000000e-01 -7.92617282e-01  1.96021708e-01  2.88675135e-01]
 [-5.00000000e-01  1.11022302e-16 -5.55111512e-17 -8.66025404e-01]
 [-5.00000000e-01  5.66068420e-01  5.88415848e-01  2.88675135e-01]]
orthogonality of eigenvectors
[[ 1.00000000e+00 -2.40880415e-17  1.99197095e-16  2.65824870e-16]
 [-2.40880415e-17  1.00000000e+00 -1.37886642e-17  2.08372994e-16]
 [ 1.99197095e-16 -1.37886642e-17  1.00000000e+00 -1.85512875e-16]
 [ 2.65824870e-16  2.08372994e-16 -1.85512875e-16  1.00000000e+00]]
reconstruct the orginal matix
[[ 1.00000000e+00 -4.68812634e-16 -1.00000000e+00 -2.12417609e-16]
 [-4.68812634e-16  1.00000000e+00 -1.00000000e+00 -5.88422541e-16]
 [-1.00000000e+00 -1.00000000e+00  3.00000000e+00 -1.00000000e+00]
 [-2.12417609e-16 -5.32911390e-16 -1.00000000e+00  1.00000000e+00]]