This is the mail archive of the gsl-discuss@sourceware.org mailing list for the GSL project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

error in Schur decomposition of a non-symmetric matrix ?


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


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]