How to make the eigenvalues and eigenvectors stay real instead of complex?

Question

I try to diagonalize a n*100*100 3d matrix K by numpy.linalg.eig and get the eigenvalues w and eigenvectors v. The matrix is 100*100, but I wanna do it with broadcasting, and that's the number n I set up. And the matrix is not hermitian.

w,v=np.linalg.eig(K)

At first, I tried n=1000, I get real eigenvalues and eigenvectors, i.e. xxxxxxxxxe+xx, but when I tried n=2000, the elements of w and v shows xxxxxxxxxe+xx+0.j. Due to +0.j, it gave complex numbers when using w and v do further calculation.

1. I thought it's because of the algorithm error for float number calculation but why?
2. Does numpy.linalg use LAPACK? Is that possible the error from LAPACK?
3. How can I get rid of +0.j?

Answer 1

According to documentation, numpy.linalg.eig uses (for real arguments) the LAPACK routine DGEEV which does not make any assumptions about the input matrix (apart from being real). If the matrix is within floating point precision sufficiently symmetric, the complex part of the returned eigenvalues will be zero (the output argument WI of DGEEV). However, due to finite precision, it might happen that you could get some spurious complex parts.

EDIT:

1. If you are sure that your matrices have only real eigenvalues, you could strip the complex part with numpy.real or use numpy.linalg.eigh specialized for symmetric matrices.

2. As for numpy.linalg.eig, the relevant part in numpy/linalg/linalg.py is:

   w, vt = _umath_linalg.eig(a, signature=signature, extobj=extobj)
   
   if not isComplexType(t) and all(w.imag == 0.0):
       w = w.real
       vt = vt.real
       result_t = _realType(result_t)
   else:
       result_t = _complexType(result_t)
   

So the test is a strict comparison all(w.imag == 0.0) and only then are the eigenvalues cast to real with w = w.real.