If $n$ is small, LAPACK is your best bet. Instead of using cgees
, I'd use zgeev
, which is a routine that will calculate eigenvalues and optionally, the left and/or right eigenvectors of a general matrix. Compared to zgeev
, cgees
uses single-precision complex numbers, rather than double-precision complex numbers, and will return the Schur form of the matrix rather than the eigenvectors you would like to calculate . There exist other, similar routines that will take advantage of symmetry (or Hermitian symmetry) if it is available for your problem. LAPACK uses dense linear algebra, which is best for small matrices (say, less than 1000 by 1000 or so, as a rough estimate; you might be able to accommodate more or less depending on how much RAM you have available).
If $n$ is large, algorithms that combine a some variant of:
- a good estimate of an eigenvalue
- a shift-and-invert strategy
- and power iteration
are used to calculate eigenvalues and eigenvectors. (More advanced algorithms use Krylov subspace iteration, transformation to Hessenberg form, and other features; broadly speaking, these can be thought of as power iteration-like algorithms with some more desirable algorithmic properties.) If you want to find the smallest eigenvalue of a large matrix, you're best off using a package like SLEPc in concert with a specialized package for eigenvalue problems, like ARPACK, BLOPEX, etc.
Generally speaking, one does not calculate all eigenvectors and eigenvalues of a large, sparse matrix. Usually, eigenvalues at the extremes -- the eigenvalues with the largest and smallest magnitudes -- are easier to calculate accurately than eigenvalues in the middle of the spectrum. As alluded to earlier, these are also possible to calculate accurately if estimates of these eigenvalues are available (using shift-and-invert). It's also expensive to calculate all of these eigenvalues, because iterative methods calculate one eigenvalue and eigenvector at a time.
As for diagonalization, the Jordan canonical form is numerically unstable, so it is not calculated. The Schur decomposition is typically calculated instead, as alluded to in one of the comments above. If $A = QUQ^{*}$ is the Schur decomposition of $A$, then the main diagonal of $U$ gives the eigenvalues of $A$. Once these are known, it's possible to calculate the eigenvectors of $A$. If $A$ is a normal matrix (that is, $A^{*}A = AA^{*}$, where the asterisk denotes the conjugate transpose, or Hermitian transpose), then the Schur decomposition is also an eigendecomposition, and diagonalizes the matrix.
As many noted, it's possible for a matrix not to be diagonalizable over the reals. However, as you correctly point out, the fundamental theorem of algebra implies that for $n$ odd, an $n \times n$ matrix must have a real eigenvalue, thus a smallest eigenvalue must exist.