STRSEN  -  reorder  the real Schur factorization of a real
       matrix A = Q*T*Q**T, so that a selected cluster of  eigen­
       values appears in the leading diagonal blocks of the upper
       quasi-triangular matrix T,


SYNOPSIS

       SUBROUTINE STRSEN( JOB, COMPQ, SELECT, N, T, LDT, Q,  LDQ,
                          WR,  WI, M, S, SEP, WORK, LWORK, IWORK,
                          LIWORK, INFO )

           CHARACTER      COMPQ, JOB

           INTEGER        INFO, LDQ, LDT, LIWORK, LWORK, M, N

           REAL           S, SEP

           LOGICAL        SELECT( * )

           INTEGER        IWORK( * )

           REAL           Q( LDQ, * ), T( LDT,  *  ),  WI(  *  ),
                          WORK( * ), WR( * )


PURPOSE

       STRSEN  reorders  the  real  Schur factorization of a real
       matrix A = Q*T*Q**T, so that a selected cluster of  eigen­
       values appears in the leading diagonal blocks of the upper
       quasi-triangular matrix T, and the leading  columns  of  Q
       form  an  orthonormal  basis  of  the  corresponding right
       invariant subspace.

       Optionally the routine computes the  reciprocal  condition
       numbers of the cluster of eigenvalues and/or the invariant
       subspace.

       T must be in Schur canonical form (as returned by SHSEQR),
       that  is,  block  upper  triangular with 1-by-1 and 2-by-2
       diagonal blocks; each 2-by-2 diagonal block has its diago­
       nal  elemnts  equal and its off-diagonal elements of oppo­
       site sign.


ARGUMENTS

       JOB     (input) CHARACTER*1
               Specifies whether condition numbers  are  required
               for  the cluster of eigenvalues (S) or the invari­
               ant subspace (SEP):
               = 'N': none;
               = 'E': for eigenvalues only (S);
               = 'V': for invariant subspace only (SEP);
               = 'B': for both eigenvalues and invariant subspace
               (S and SEP).

               = 'V': update the matrix Q of Schur vectors;
               = 'N': do not update Q.

       SELECT  (input) LOGICAL array, dimension (N)
               SELECT  specifies  the eigenvalues in the selected
               cluster.  To  select  a  real   eigenvalue   w(j),
               SELECT(j)  must  be set to w(j) and w(j+1), corre­
               sponding  to  a  2-by-2  diagonal  block,   either
               SELECT(j)  or  SELECT(j+1)  or both must be set to
               either  both  included  in  the  cluster  or  both
               excluded.

       N       (input) INTEGER
               The order of the matrix T. N >= 0.

       T       (input/output) REAL array, dimension (LDT,N)
               On  entry, the upper quasi-triangular matrix T, in
               Schur canonical form.  On exit, T  is  overwritten
               by  the reordered matrix T, again in Schur canoni­
               cal form, with the  selected  eigenvalues  in  the
               leading diagonal blocks.

       LDT     (input) INTEGER
               The  leading  dimension  of  the  array  T. LDT >=
               max(1,N).

       Q       (input/output) REAL array, dimension (LDQ,N)
               On entry, if COMPQ = 'V', the matrix  Q  of  Schur
               vectors.   On  exit,  if  COMPQ  = 'V', Q has been
               postmultiplied by  the  orthogonal  transformation
               matrix  which reorders T; the leading M columns of
               Q form an  orthonormal  basis  for  the  specified
               invariant subspace.  If COMPQ = 'N', Q is not ref­
               erenced.

       LDQ     (input) INTEGER
               The leading dimension of the array Q.  LDQ  >=  1;
               and if COMPQ = 'V', LDQ >= N.

       WR      (output) REAL array, dimension (N)
               WI       (output)  REAL  array,  dimension (N) The
               real and imaginary  parts,  respectively,  of  the
               reordered  eigenvalues  of  T. The eigenvalues are
               stored in the same order as on the diagonal of  T,
               with  WR(i)  =  T(i,i) and, if T(i:i+1,i:i+1) is a
               2-by-2 diagonal block, WI(i) >  0  and  WI(i+1)  =
               -WI(i).  Note that if a complex eigenvalue is suf­
               ficiently ill-conditioned, then its value may dif­
               fer  significantly  from its value before reorder­
               ing.

       M       (output) INTEGER
               The dimension of the specified invariant subspace.


       S       (output) REAL
               If  JOB  =  'E'  or 'B', S is a lower bound on the
               reciprocal condition number for the selected clus­
               ter  of  eigenvalues.   S cannot underestimate the
               true reciprocal condition number by  more  than  a
               factor of sqrt(N). If M = 0 or N, S = 1.  If JOB =
               'N' or 'V', S is not referenced.

       SEP     (output) REAL
               If JOB = 'V' or 'B', SEP is the estimated recipro­
               cal  condition  number  of the specified invariant
               subspace. If M = 0 or N, SEP = norm(T).  If JOB  =
               'N' or 'E', SEP is not referenced.

       WORK    (workspace/output) REAL array, dimension (LWORK)
               On  exit, if INFO = 0, WORK(1) returns the optimal
               LWORK.

       LWORK   (input) INTEGER
               The dimension of the array WORK.  If  JOB  =  'N',
               LWORK >= max(1,N); if JOB = 'E', LWORK >= M*(N-M);
               if JOB = 'V' or 'B', LWORK >= 2*M*(N-M).

               If LWORK = -1, then a workspace query is  assumed;
               the  routine  only  calculates the optimal size of
               the WORK array, returns this value  as  the  first
               entry  of  the  WORK  array,  and no error message
               related to LWORK is issued by XERBLA.

       IWORK   (workspace) INTEGER array, dimension (LIWORK)
               IF JOB = 'N' or 'E', IWORK is not referenced.

       LIWORK  (input) INTEGER
               The dimension of the array IWORK.  If JOB = 'N' or
               'E',  LIWORK  >= 1; if JOB = 'V' or 'B', LIWORK >=
               M*(N-M).

               If LIWORK = -1, then a workspace query is assumed;
               the  routine  only  calculates the optimal size of
               the IWORK array, returns this value as  the  first
               entry  of  the  IWORK  array, and no error message
               related to LIWORK is issued by XERBLA.

       INFO    (output) INTEGER
               = 0: successful exit
               < 0: if INFO = -i, the i-th argument had an  ille­
               gal value
               = 1: reordering of T failed because some eigenval­
               ues are too close to separate (the problem is very
               ill-conditioned);   T   may  have  been  partially
               reordered, and WR and WI contain  the  eigenvalues
               in  the  same  order  as  in  T;  S  and  SEP  (if



FURTHER DETAILS

       STRSEN first collects the selected eigenvalues by  comput­
       ing an orthogonal transformation Z to move them to the top
       left corner of T.  In other words, the selected  eigenval­
       ues are the eigenvalues of T11 in:

                     Z'*T*Z = ( T11 T12 ) n1
                              (  0  T22 ) n2
                                 n1  n2

       where N = n1+n2 and Z' means the transpose of Z. The first
       n1 columns of Z span the specified invariant  subspace  of
       T.

       If  T  has been obtained from the real Schur factorization
       of a matrix A = Q*T*Q', then the reordered real Schur fac­
       torization of A is given by A = (Q*Z)*(Z'*T*Z)*(Q*Z)', and
       the first n1 columns of Q*Z span the corresponding invari­
       ant subspace of A.

       The  reciprocal  condition  number  of  the average of the
       eigenvalues of T11 may be returned in S. S lies between  0
       (very badly conditioned) and 1 (very well conditioned). It
       is computed as follows. First we compute R so that

                              P = ( I  R ) n1
                                  ( 0  0 ) n2
                                    n1 n2

       is the projector on the invariant subspace associated with
       T11.  R is the solution of the Sylvester equation:

                             T11*R - R*T22 = T12.

       Let F-norm(M) denote the Frobenius-norm of M and 2-norm(M)
       denote the two-norm of M. Then S is computed as the  lower
       bound

                           (1 + F-norm(R)**2)**(-1/2)

       on the reciprocal of 2-norm(P), the true reciprocal condi­
       tion number.  S cannot underestimate 1 / 2-norm(P) by more
       than a factor of sqrt(N).

       An approximate error bound for the computed average of the
       eigenvalues of T11 is

                              EPS * norm(T) / S

       where EPS is the machine precision.

       The reciprocal condition number  of  the  right  invariant

       is returned in SEP.  SEP is defined as the  separation  of
       T11 and T22:

                          sep( T11, T22 ) = sigma-min( C )

       where sigma-min(C) is the smallest singular value of the
       n1*n2-by-n1*n2 matrix

          C  = kprod( I(n2), T11 ) - kprod( transpose(T22), I(n1)
       )

       I(m) is an m by m identity matrix, and kprod  denotes  the
       Kronecker  product. We estimate sigma-min(C) by the recip­
       rocal of an estimate of the 1-norm of inverse(C). The true
       reciprocal  1-norm of inverse(C) cannot differ from sigma-
       min(C) by more than a factor of sqrt(n1*n2).

       When SEP is small, small changes  in  T  can  cause  large
       changes in the invariant subspace. An approximate bound on
       the maximum angular error in the computed right  invariant
       subspace is

                           EPS * norm(T) / SEP


Man(1) output converted with man2html