error in Schur decomposition of a non-symmetric matrix ?
Francesco Abbate
francesco.bbt@gmail.com
Tue Mar 23 15:33:00 GMT 2010
Hi all,
it seems that I've found an error in GSL manual in the section about
eigenvalues/vectors determination of non-symmetric real matrices. The
manual says:
14.3 Real Nonsymmetric Matrices
===============================
The solution of the real nonsymmetric eigensystem problem for a matrix
A involves computing the Schur decomposition
A = Z T Z^T
where Z is an orthogonal matrix of Schur vectors and T, the Schur
form, is quasi upper triangular with diagonal 1-by-1 blocks which are
real eigenvalues of A, and diagonal 2-by-2 blocks whose eigenvalues are
complex conjugate eigenvalues of A. The algorithm used is the
double-shift Francis method.
-----------------------
and after about the function gsl_eigen_nonsymmv and gsl_eigen_nonsymmv_Z :
-- Function: int gsl_eigen_nonsymmv (gsl_matrix * A,
gsl_vector_complex * EVAL, gsl_matrix_complex * EVEC,
gsl_eigen_nonsymmv_workspace * W)
This function computes eigenvalues and right eigenvectors of the
N-by-N real nonsymmetric matrix A. It first calls
`gsl_eigen_nonsymm' to compute the eigenvalues, Schur form T, and
Schur vectors. Then it finds eigenvectors of T and backtransforms
them using the Schur vectors. The Schur vectors are destroyed in
the process, but can be saved by using `gsl_eigen_nonsymmv_Z'. The
computed eigenvectors are normalized to have unit magnitude. On
output, the upper portion of A contains the Schur form T. If
`gsl_eigen_nonsymm' fails, no eigenvectors are computed, and an
error code is returned.
-- Function: int gsl_eigen_nonsymmv_Z (gsl_matrix * A,
gsl_vector_complex * EVAL, gsl_matrix_complex * EVEC,
gsl_matrix * Z, gsl_eigen_nonsymmv_workspace * W)
This function is identical to `gsl_eigen_nonsymmv' except that it
also saves the Schur vectors into Z.
------------------------------
The problems seems to be about the Schur decomposition formula. What
it is true is that:
A = Z T Z^(-1)
and not that: A = Z T Z^T as asserted in the documentation. I've found
that empirically.
Otherwise I'm perplexed about the statement in the gsl_eigen_nonsymmv:
"On output, the upper portion of A contains the Schur form T". I don't
understand this statement because before it was stated that: "T, the
Schur form, is quasi upper triangular with diagonal 1-by-1 blocks
which are real eigenvalues of A, and diagonal 2-by-2 blocks whose
eigenvalues are complex conjugate eigenvalues of A". So the T matrix
is not really upper triangular.
I've made some tests and it seems that the function gsl_eigen_nonsymmv
set the matrix A to the Schur form T and not only the upper part.
It would be nice if someone could check because what I've seen is
based on empirical tests using the functions.
Best regards,
Francesco
More information about the Gsl-discuss
mailing list