Why isn't `np.linalg.eig` returning a unitary matrix when I try to diagonalize a symmetric matrix?

Viewed 37

Using Python, I created a symmetric matrix X, and diagonalized it using the code Lambda, U = np.linalg.eig(X). My understanding is: since X is symmetric, U should be unitary, but I am finding that this is not the case, as U has eigenvalues that don't have absolute value 1, which also means that this is not a scaling issue either. My question is why is U not unitary?


Here is a minimum reproducible example to see what I am getting:

import numpy as np

#Symmetric matrix, which should be diagonalized by unitary
X = np.array([[-1.1918157 ,  0.        ,  0.        ,  0.        ,  0.09852097,
         0.        ,  0.        ,  0.        ],
       [ 0.        , -1.1918157 ,  0.        ,  0.        ,  0.        ,
         0.09852097,  0.        ,  0.        ],
       [ 0.        ,  0.        , -1.08529969,  0.        ,  0.        ,
         0.        ,  0.07826825,  0.        ],
       [ 0.        ,  0.        ,  0.        , -1.08529969,  0.        ,
         0.        ,  0.        ,  0.07826825],
       [ 0.09852097,  0.        ,  0.        ,  0.        , -1.00585682,
         0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.09852097,  0.        ,  0.        ,  0.        ,
        -1.00585682,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.07826825,  0.        ,  0.        ,
         0.        , -0.98276093,  0.        ],
       [ 0.        ,  0.        ,  0.        ,  0.07826825,  0.        ,
         0.        ,  0.        , -0.98276093]])

#check if X is symmetric
assert np.linalg.norm(X - X.conj().T) == 0

#diagonalize X
Lambda, U = np.linalg.eig(X)

#check if U is unitary
print(
    np.allclose(U @ U.conj().T, np.identity(U.shape[0])),
    np.allclose(U.conj().T,     np.linalg.inv(U))
)

#Since the above returns false, how non-unitary is it?
print(
    np.linalg.norm(U @ U.conj().T - np.identity(U.shape[0])),
    np.linalg.norm(U.conj().T - np.linalg.inv(U))
)

The output of the above code is:

False False
0.994793051566498 1.1641788210519939
1 Answers

It's never guaranteed to return a unitary matrix.

#diagonalize X
Lambda, U = np.linalg.eig(X)

From this line of code, you are finding eigenvalues of X (stored in Lambda) and its right eigenvectors (stored in matrix U). See the doc.

The only unitary things here are the eigenvectors, which are the columns of U.

In [4]: np.linalg.norm(U, axis=0)
Out[4]: array([1., 1., 1., 1., 1., 1., 1., 1.])

Read more about eigenvalues and eigenvectors: wiki, wolfram

Related