Tridiagonal Decomposition of Hermitian Matrices
A hermitian matrix A can be factorized by similarity
transformations into the form,
A = U T U^T
where U is an unitary matrix and T is a real symmetric
tridiagonal matrix.
- Function: int gsl_linalg_hermtd_decomp (gsl_matrix_complex * A, gsl_vector_complex * tau)
-
This function factorizes the hermitian matrix A into the symmetric
tridiagonal decomposition U T U^T. On output the real parts of
the diagonal and subdiagonal part of the input matrix A contain
the tridiagonal matrix T. The remaining lower triangular part of
the input matrix contains the Householder vectors which, together with
the Householder coefficients tau, encode the orthogonal matrix
Q. This storage scheme is the same as used by LAPACK. The
upper triangular part of A and imaginary parts of the diagonal are
not referenced.
- Function: int gsl_linalg_hermtd_unpack (const gsl_matrix_complex * A, const gsl_vector_complex * tau, gsl_matrix_complex * Q, gsl_vector * diag, gsl_vector * subdiag)
-
This function unpacks the encoded tridiagonal decomposition (A,
tau) obtained from
gsl_linalg_hermtd_decomp
into the
unitary matrix U, the real vector of diagonal elements diag and
the real vector of subdiagonal elements subdiag.
- Function: int gsl_linalg_hermtd_unpack_T (const gsl_matrix_complex * A, gsl_vector * diag, gsl_vector * subdiag)
-
This function unpacks the diagonal and subdiagonal of the encoded
tridiagonal decomposition (A, tau) obtained from
gsl_linalg_hermtd_decomp
into the real vectors diag and
subdiag.