This chapter describes functions for computing eigenvalues and
eigenvectors of matrices. There are routines for real symmetric and
complex hermitian matrices, and eigenvalues can be computed with or
without eigenvectors. The algorithms used are symmetric
bidiagonalization followed by QR reduction.
These routines are intended for "small" systems where simple algorithms are
acceptable. Anyone interested finding eigenvalues and eigenvectors of
large matrices will want to use the sophisticated routines found in
LAPACK. The Fortran version of LAPACK is recommended as the
standard package for linear algebra.
The functions described in this chapter are declared in the header file
`gsl_eigen.h'.
This function frees the memory associated with the workspace w.
Function: int gsl_eigen_symm(gsl_matrix * A, gsl_vector * eval, gsl_eigen_symm_workspace * w)
This function computes the eigenvalues of the real symmetric matrix
A. Additional workspace of the appropriate size must be provided
in w. The diagonal and lower triangular part of A are
destroyed during the computation, but the strict upper triangular part
is not referenced. The eigenvalues are stored in the vector eval
and are unordered.
This function frees the memory associated with the workspace w.
Function: int gsl_eigen_symmv(gsl_matrix * A, gsl_vector * eval, gsl_matrix * evec, gsl_eigen_symmv_workspace * w)
This function computes the eigenvalues and eigenvectors of the real
symmetric matrix A. Additional workspace of the appropriate size
must be provided in w. The diagonal and lower triangular part of
A are destroyed during the computation, but the strict upper
triangular part is not referenced. The eigenvalues are stored in the
vector eval and are unordered. The corresponding eigenvectors are
stored in the columns of the matrix evec. For example, the
eigenvector in the first column corresponds to the first eigenvalue.
The eigenvectors are guaranteed to be mutually orthogonal and normalised
to unit magnitude.
This function frees the memory associated with the workspace w.
Function: int gsl_eigen_herm(gsl_matrix_complex * A, gsl_vector * eval, gsl_eigen_herm_workspace * w)
This function computes the eigenvalues of the complex hermitian matrix
A. Additional workspace of the appropriate size must be provided
in w. The diagonal and lower triangular part of A are
destroyed during the computation, but the strict upper triangular part
is not referenced. The imaginary parts of the diagonal are assumed to be
zero and are not referenced. The eigenvalues are stored in the vector
eval and are unordered.
This function allocates a workspace for computing eigenvalues and
eigenvectors of n-by-n complex hermitian matrices. The size of
the workspace is O(5n).
This function frees the memory associated with the workspace w.
Function: int gsl_eigen_hermv(gsl_matrix_complex * A, gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_hermv_workspace * w)
This function computes the eigenvalues and eigenvectors of the complex
hermitian matrix A. Additional workspace of the appropriate size
must be provided in w. The diagonal and lower triangular part of
A are destroyed during the computation, but the strict upper
triangular part is not referenced. The imaginary parts of the diagonal
are assumed to be zero and are not referenced. The eigenvalues are
stored in the vector eval and are unordered. The corresponding
complex eigenvectors are stored in the columns of the matrix evec.
For example, the eigenvector in the first column corresponds to the
first eigenvalue. The eigenvectors are guaranteed to be mutually
orthogonal and normalised to unit magnitude.
Function: int gsl_eigen_symmv_sort(gsl_vector * eval, gsl_matrix * evec, gsl_eigen_sort_t sort_type)
This function simultaneously sorts the eigenvalues stored in the vector
eval and the corresponding real eigenvectors stored in the columns
of the matrix evec into ascending or descending order according to
the value of the parameter sort_type,
GSL_EIGEN_SORT_VAL_ASC
ascending order in numerical value
GSL_EIGEN_SORT_VAL_DESC
descending order in numerical value
GSL_EIGEN_SORT_ABS_ASC
ascending order in magnitude
GSL_EIGEN_SORT_ABS_DESC
descending order in magnitude
Function: int gsl_eigen_hermv_sort(gsl_vector * eval, gsl_matrix_complex * evec, gsl_eigen_sort_t sort_type)
This function simultaneously sorts the eigenvalues stored in the vector
eval and the corresponding complex eigenvectors stored in the
columns of the matrix evec into ascending or descending order
according to the value of the parameter sort_type as shown above.