Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

ctgsna (3p)

Name

ctgsna - ues and/or eigenvectors of a matrix pair (A, B)

Synopsis

SUBROUTINE CTGSNA(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL, LDVL,
VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)

CHARACTER*1 JOB, HOWMNT
COMPLEX A(LDA,*), B(LDB,*), VL(LDVL,*), VR(LDVR,*), WORK(*)
INTEGER N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
INTEGER IWORK(*)
LOGICAL SELECT(*)
REAL S(*), DIF(*)

SUBROUTINE CTGSNA_64(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)

CHARACTER*1 JOB, HOWMNT
COMPLEX A(LDA,*), B(LDB,*), VL(LDVL,*), VR(LDVR,*), WORK(*)
INTEGER*8 N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
INTEGER*8 IWORK(*)
LOGICAL*8 SELECT(*)
REAL S(*), DIF(*)




F95 INTERFACE
SUBROUTINE TGSNA(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK,
INFO)

CHARACTER(LEN=1) :: JOB, HOWMNT
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, B, VL, VR
INTEGER :: N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
INTEGER, DIMENSION(:) :: IWORK
LOGICAL, DIMENSION(:) :: SELECT
REAL, DIMENSION(:) :: S, DIF

SUBROUTINE TGSNA_64(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK,
INFO)

CHARACTER(LEN=1) :: JOB, HOWMNT
COMPLEX, DIMENSION(:) :: WORK
COMPLEX, DIMENSION(:,:) :: A, B, VL, VR
INTEGER(8) :: N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
INTEGER(8), DIMENSION(:) :: IWORK
LOGICAL(8), DIMENSION(:) :: SELECT
REAL, DIMENSION(:) :: S, DIF




C INTERFACE
#include <sunperf.h>

void ctgsna(char job, char howmnt, int *select, int n, complex *a,  int
lda, complex *b, int ldb, complex *vl, int ldvl, complex *vr,
int ldvr, float *s, float *dif, int mm, int *m, int *info);

void ctgsna_64(char job, char howmnt, long *select, long n, complex *a,
long  lda, complex *b, long ldb, complex *vl, long ldvl, com-
plex *vr, long ldvr, float *s, float *dif, long mm, long  *m,
long *info);

Description

Oracle Solaris Studio Performance Library                           ctgsna(3P)



NAME
       ctgsna  - estimate reciprocal condition numbers for specified eigenval-
       ues and/or eigenvectors of a matrix pair (A, B)


SYNOPSIS
       SUBROUTINE CTGSNA(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL, LDVL,
             VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)

       CHARACTER*1 JOB, HOWMNT
       COMPLEX A(LDA,*), B(LDB,*), VL(LDVL,*), VR(LDVR,*), WORK(*)
       INTEGER N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
       INTEGER IWORK(*)
       LOGICAL SELECT(*)
       REAL S(*), DIF(*)

       SUBROUTINE CTGSNA_64(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
             LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)

       CHARACTER*1 JOB, HOWMNT
       COMPLEX A(LDA,*), B(LDB,*), VL(LDVL,*), VR(LDVR,*), WORK(*)
       INTEGER*8 N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
       INTEGER*8 IWORK(*)
       LOGICAL*8 SELECT(*)
       REAL S(*), DIF(*)




   F95 INTERFACE
       SUBROUTINE TGSNA(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
              LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK,
              INFO)

       CHARACTER(LEN=1) :: JOB, HOWMNT
       COMPLEX, DIMENSION(:) :: WORK
       COMPLEX, DIMENSION(:,:) :: A, B, VL, VR
       INTEGER :: N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
       INTEGER, DIMENSION(:) :: IWORK
       LOGICAL, DIMENSION(:) :: SELECT
       REAL, DIMENSION(:) :: S, DIF

       SUBROUTINE TGSNA_64(JOB, HOWMNT, SELECT, N, A, LDA, B, LDB, VL,
              LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK,
              INFO)

       CHARACTER(LEN=1) :: JOB, HOWMNT
       COMPLEX, DIMENSION(:) :: WORK
       COMPLEX, DIMENSION(:,:) :: A, B, VL, VR
       INTEGER(8) :: N, LDA, LDB, LDVL, LDVR, MM, M, LWORK, INFO
       INTEGER(8), DIMENSION(:) :: IWORK
       LOGICAL(8), DIMENSION(:) :: SELECT
       REAL, DIMENSION(:) :: S, DIF




   C INTERFACE
       #include <sunperf.h>

       void ctgsna(char job, char howmnt, int *select, int n, complex *a,  int
                 lda, complex *b, int ldb, complex *vl, int ldvl, complex *vr,
                 int ldvr, float *s, float *dif, int mm, int *m, int *info);

       void ctgsna_64(char job, char howmnt, long *select, long n, complex *a,
                 long  lda, complex *b, long ldb, complex *vl, long ldvl, com-
                 plex *vr, long ldvr, float *s, float *dif, long mm, long  *m,
                 long *info);



PURPOSE
       ctgsna estimates reciprocal condition numbers for specified eigenvalues
       and/or eigenvectors of a matrix pair (A, B).

       (A, B) must be in generalized Schur canonical form, that is,  A  and  B
       are both upper triangular.


ARGUMENTS
       JOB (input)
                 Specifies  whether  condition numbers are required for eigen-
                 values (S) or eigenvectors (DIF):
                 = 'E': for eigenvalues only (S);
                 = 'V': for eigenvectors only (DIF);
                 = 'B': for both eigenvalues and eigenvectors (S and DIF).


       HOWMNT (input)
                 = 'A': compute condition numbers for all eigenpairs;
                 = 'S': compute  condition  numbers  for  selected  eigenpairs
                 specified by the array SELECT.


       SELECT (input)
                 If  HOWMNT  =  'S', SELECT specifies the eigenpairs for which
                 condition numbers are required. To select  condition  numbers
                 for  the  corresponding  j-th  eigenvalue and/or eigenvector,
                 SELECT(j) must be set to .TRUE..  If HOWMNT = 'A', SELECT  is
                 not referenced.


       N (input) The order of the square matrix pair (A, B). N >= 0.


       A (input) The upper triangular matrix A in the pair (A,B).


       LDA (input)
                 The leading dimension of the array A. LDA >= max(1,N).


       B (input) The upper triangular matrix B in the pair (A, B).


       LDB (input)
                 The leading dimension of the array B. LDB >= max(1,N).


       VL (input)
                 If JOB = 'E' or 'B', VL must contain left eigenvectors of (A,
                 B), corresponding to the eigenpairs specified by  HOWMNT  and
                 SELECT.   The eigenvectors must be stored in consecutive col-
                 umns of VL, as returned by CTGEVC.  If JOB = 'V', VL  is  not
                 referenced.


       LDVL (input)
                 The  leading dimension of the array VL. LDVL >= 1; and If JOB
                 = 'E' or 'B', LDVL >= N.


       VR (input)
                 If JOB = 'E' or 'B', VR must contain  right  eigenvectors  of
                 (A,  B),  corresponding to the eigenpairs specified by HOWMNT
                 and SELECT.  The eigenvectors must be stored  in  consecutive
                 columns  of  VR,  as returned by CTGEVC.  If JOB = 'V', VR is
                 not referenced.


       LDVR (input)
                 The leading dimension of the array VR. LDVR >= 1;  If  JOB  =
                 'E' or 'B', LDVR >= N.


       S (output)
                 If  JOB = 'E' or 'B', the reciprocal condition numbers of the
                 selected eigenvalues, stored in consecutive elements  of  the
                 array.  If JOB = 'V', S is not referenced.


       DIF (output)
                 If  JOB = 'V' or 'B', the estimated reciprocal condition num-
                 bers of the selected eigenvectors, stored in consecutive ele-
                 ments  of  the array.  If the eigenvalues cannot be reordered
                 to compute DIF(j), DIF(j) is set to 0; this  can  only  occur
                 when the true value would be very small anyway.  For each ei-
                 genvalue/vector specified by SELECT, DIF stores  a  Frobenius
                 norm-based estimate of Difl.  If JOB = 'E', DIF is not refer-
                 enced.


       MM (input)
                 The number of elements in the arrays S and DIF. MM >= M.


       M (output)
                 The number of elements of the arrays S and DIF used to  store
                 the specified condition numbers; for each selected eigenvalue
                 one element is used. If HOWMNT = 'A', M is set to N.


       WORK (workspace)
                 If JOB = 'E', WORK is not referenced.  Otherwise, on exit, if
                 INFO = 0, WORK(1) returns the optimal LWORK.


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


       IWORK (workspace)
                 dimension(N+2) If JOB = 'E', IWORK is not referenced.


       INFO (output)
                 = 0: Successful exit
                 < 0: If INFO = -i, the i-th argument had an illegal value

FURTHER DETAILS
       The reciprocal of the condition number of the i-th  generalized  eigen-
       value w = (a, b) is defined as

               S(I) = (|v'Au|**2 + |v'Bu|**2)**(1/2) / (norm(u)*norm(v))

       where u and v are the right and left eigenvectors of (A, B) correspond-
       ing to w; |z| denotes the absolute value of  the  complex  number,  and
       norm(u) denotes the 2-norm of the vector u. The pair (a, b) corresponds
       to an eigenvalue w = a/b (= v'Au/v'Bu) of the matrix pair  (A,  B).  If
       both  a  and  b  equal  zero,  then  (A,B) is singular and S(I) = -1 is
       returned.

       An approximate error bound on the chordal  distance  between  the  i-th
       computed generalized eigenvalue w and the corresponding exact eigenval-
       ue lambda is

               chord(w, lambda) <=   EPS * norm(A, B) / S(I),

       where EPS is the machine precision.

       The reciprocal of the condition number of the right eigenvector  u  and
       left  eigenvector  v  corresponding  to the generalized eigenvalue w is
       defined as follows. Suppose

                        (A, B) = ( a   *  ) ( b  *  )  1
                                 ( 0  A22 ),( 0 B22 )  n-1
                                   1  n-1     1 n-1

       Then the reciprocal condition number DIF(I) is

               Difl[(a, b), (A22, B22)]  = sigma-min( Zl )

       where sigma-min(Zl) denotes the smallest singular value of

              Zl = [ kron(a, In-1) -kron(1, A22) ]
                   [ kron(b, In-1) -kron(1, B22) ].

       Here In-1 is the identity matrix of size n-1 and X'  is  the  conjugate
       transpose  of X. kron(X, Y) is the Kronecker product between the matri-
       ces X and Y.

       We approximate the smallest singular value of Zl with an  upper  bound.
       This is done by CLATDF.

       An approximate error bound for a computed eigenvector VL(i) or VR(i) is
       given by

                           EPS * norm(A, B) / DIF(i).

       See ref. [2-3] for more details and further references.

       Based on contributions by
          Bo Kagstrom and Peter Poromaa, Department of Computing Science,
          Umea University, S-901 87 Umea, Sweden.

       References
       ==========

       [1] B. Kagstrom; A Direct Method for Reordering Eigenvalues in the
           Generalized Real Schur Form of a Regular Matrix Pair (A, B), in
           M.S. Moonen et al (eds), Linear Algebra for Large Scale and
           Real-Time Applications, Kluwer Academic Publ. 1993, pp 195-218.

       [2] B. Kagstrom and P. Poromaa; Computing Eigenspaces with Specified
           Eigenvalues of a Regular Matrix Pair (A, B) and Condition
           Estimation: Theory, Algorithms and Software, Report
           UMINF - 94.04, Department of Computing Science, Umea University,
           S-901 87 Umea, Sweden, 1994. Also as LAPACK Working Note 87.
           To appear in Numerical Algorithms, 1996.

       [3] B. Kagstrom and P. Poromaa, LAPACK-Style Algorithms and Software
           for Solving the Generalized Sylvester Equation and Estimating the
           Separation between Regular Matrix Pairs, Report UMINF - 93.23,
           Department of Computing Science, Umea University, S-901 87 Umea,
           Sweden, December 1993, Revised April 1994, Also as LAPACK Working
           Note 75.
           To appear in ACM Trans. on Math. Software, Vol 22, No 1, 1996.




                                  7 Nov 2015                        ctgsna(3P)