Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

dgges (3p)

Name

dgges - N real nonsymmetric matrices (A,B),

Synopsis

SUBROUTINE DGGES(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK,
BWORK, INFO)

CHARACTER*1 JOBVSL, JOBVSR, SORT
INTEGER N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
LOGICAL DELCTG
LOGICAL BWORK(*)
DOUBLE  PRECISION  A(LDA,*),  B(LDB,*),  ALPHAR(*), ALPHAI(*), BETA(*),
VSL(LDVSL,*), VSR(LDVSR,*), WORK(*)

SUBROUTINE DGGES_64(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK,
BWORK, INFO)

CHARACTER*1 JOBVSL, JOBVSR, SORT
INTEGER*8 N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
LOGICAL*8 DELCTG
LOGICAL*8 BWORK(*)
DOUBLE PRECISION A(LDA,*),  B(LDB,*),  ALPHAR(*),  ALPHAI(*),  BETA(*),
VSL(LDVSL,*), VSR(LDVSR,*), WORK(*)




F95 INTERFACE
SUBROUTINE GGES(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK,
LWORK, BWORK, INFO)

CHARACTER(LEN=1) :: JOBVSL, JOBVSR, SORT
INTEGER :: N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
LOGICAL :: DELCTG
LOGICAL, DIMENSION(:) :: BWORK
REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
REAL(8), DIMENSION(:,:) :: A, B, VSL, VSR

SUBROUTINE GGES_64(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B,
LDB, SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR,
WORK, LWORK, BWORK, INFO)

CHARACTER(LEN=1) :: JOBVSL, JOBVSR, SORT
INTEGER(8) :: N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
LOGICAL(8) :: DELCTG
LOGICAL(8), DIMENSION(:) :: BWORK
REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
REAL(8), DIMENSION(:,:) :: A, B, VSL, VSR




C INTERFACE
#include <sunperf.h>

void  dgges(char  jobvsl,  char  jobvsr,  char  sort, int(*delctg)(dou-
ble,double,double), int n, double *a, int lda, double *b, int
ldb, int *sdim, double *alphar, double *alphai, double *beta,
double *vsl, int ldvsl, double *vsr, int ldvsr, int *info);

void dgges_64(char jobvsl, char jobvsr, char  sort,  long(*delctg)(dou-
ble,double,double),  long  n, double *a, long lda, double *b,
long ldb, long *sdim, double *alphar, double *alphai,  double
*beta, double *vsl, long ldvsl, double *vsr, long ldvsr, long
*info);

Description

Oracle Solaris Studio Performance Library                            dgges(3P)



NAME
       dgges - compute for a pair of N-by-N real nonsymmetric matrices (A,B),


SYNOPSIS
       SUBROUTINE DGGES(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
             SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK,
             BWORK, INFO)

       CHARACTER*1 JOBVSL, JOBVSR, SORT
       INTEGER N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
       LOGICAL DELCTG
       LOGICAL BWORK(*)
       DOUBLE  PRECISION  A(LDA,*),  B(LDB,*),  ALPHAR(*), ALPHAI(*), BETA(*),
       VSL(LDVSL,*), VSR(LDVSR,*), WORK(*)

       SUBROUTINE DGGES_64(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
             SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK,
             BWORK, INFO)

       CHARACTER*1 JOBVSL, JOBVSR, SORT
       INTEGER*8 N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
       LOGICAL*8 DELCTG
       LOGICAL*8 BWORK(*)
       DOUBLE PRECISION A(LDA,*),  B(LDB,*),  ALPHAR(*),  ALPHAI(*),  BETA(*),
       VSL(LDVSL,*), VSR(LDVSR,*), WORK(*)




   F95 INTERFACE
       SUBROUTINE GGES(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B, LDB,
              SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR, WORK,
              LWORK, BWORK, INFO)

       CHARACTER(LEN=1) :: JOBVSL, JOBVSR, SORT
       INTEGER :: N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
       LOGICAL :: DELCTG
       LOGICAL, DIMENSION(:) :: BWORK
       REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
       REAL(8), DIMENSION(:,:) :: A, B, VSL, VSR

       SUBROUTINE GGES_64(JOBVSL, JOBVSR, SORT, DELCTG, N, A, LDA, B,
              LDB, SDIM, ALPHAR, ALPHAI, BETA, VSL, LDVSL, VSR, LDVSR,
              WORK, LWORK, BWORK, INFO)

       CHARACTER(LEN=1) :: JOBVSL, JOBVSR, SORT
       INTEGER(8) :: N, LDA, LDB, SDIM, LDVSL, LDVSR, LWORK, INFO
       LOGICAL(8) :: DELCTG
       LOGICAL(8), DIMENSION(:) :: BWORK
       REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
       REAL(8), DIMENSION(:,:) :: A, B, VSL, VSR




   C INTERFACE
       #include <sunperf.h>

       void  dgges(char  jobvsl,  char  jobvsr,  char  sort, int(*delctg)(dou-
                 ble,double,double), int n, double *a, int lda, double *b, int
                 ldb, int *sdim, double *alphar, double *alphai, double *beta,
                 double *vsl, int ldvsl, double *vsr, int ldvsr, int *info);

       void dgges_64(char jobvsl, char jobvsr, char  sort,  long(*delctg)(dou-
                 ble,double,double),  long  n, double *a, long lda, double *b,
                 long ldb, long *sdim, double *alphar, double *alphai,  double
                 *beta, double *vsl, long ldvsl, double *vsr, long ldvsr, long
                 *info);



PURPOSE
       dgges computes for a pair of N-by-N real nonsymmetric  matrices  (A,B),
       the  generalized  eigenvalues,  the  generalized real Schur form (S,T),
       optionally, the left and/or right matrices of Schur  vectors  (VSL  and
       VSR). This gives the generalized Schur factorization

                (A,B) = ( (VSL)*S*(VSR)**T, (VSL)*T*(VSR)**T )

       Optionally,  it  also orders the eigenvalues so that a selected cluster
       of eigenvalues appears in the leading  diagonal  blocks  of  the  upper
       quasi-triangular matrix S and the upper triangular matrix T.The leading
       columns of VSL and VSR then form an orthonormal basis  for  the  corre-
       sponding left and right eigenspaces (deflating subspaces).

       (If  only  the generalized eigenvalues are needed, use the driver DGGEV
       instead, which is faster.)

       A generalized eigenvalue for a pair of matrices (A,B) is a scalar w  or
       a  ratio alpha/beta = w, such that  A - w*B is singular.  It is usually
       represented as the pair (alpha,beta), as there is a  reasonable  inter-
       pretation for beta=0 or both being zero.

       A  pair  of  matrices  (S,T)  is in generalized real Schur form if T is
       upper triangular with non-negative diagonal and S is block upper trian-
       gular  with 1-by-1 and 2-by-2 blocks.  1-by-1 blocks correspond to real
       generalized eigenvalues, while 2-by-2 blocks of S  will  be  "standard-
       ized" by making the corresponding elements of T have the form:
               [  a  0  ]
               [  0  b  ]

       and the pair of corresponding 2-by-2 blocks in S and T will have a com-
       plex conjugate pair of generalized eigenvalues.


ARGUMENTS
       JOBVSL (input)
                 = 'N':  do not compute the left Schur vectors;
                 = 'V':  compute the left Schur vectors.


       JOBVSR (input)
                 = 'N':  do not compute the right Schur vectors;
                 = 'V':  compute the right Schur vectors.


       SORT (input)
                 Specifies whether or not to  order  the  eigenvalues  on  the
                 diagonal  of the generalized Schur form.  = 'N':  Eigenvalues
                 are not ordered;
                 = 'S':  Eigenvalues are ordered (see DELCTG);


       DELCTG (input)
                 LOGICAL FUNCTION of three DOUBLE PRECISION  arguments  DELCTG
                 must be declared EXTERNAL in the calling subroutine.  If SORT
                 = 'N', DELCTG is not referenced.  If SORT =  'S',  DELCTG  is
                 used  to  select  eigenvalues  to sort to the top left of the
                 Schur form.  An eigenvalue  (ALPHAR(j)+ALPHAI(j))/BETA(j)  is
                 selected if DELCTG(ALPHAR(j),ALPHAI(j),BETA(j)) is true; i.e.
                 if either one of a complex conjugate pair of  eigenvalues  is
                 selected, then both complex eigenvalues are selected.

                 Note that in the ill-conditioned case, a selected complex ei-
                 genvalue may no  longer  satisfy  DELCTG(ALPHAR(j),ALPHAI(j),
                 BETA(j)) = .TRUE. after ordering. INFO is to be set to N+2 in
                 this case.


       N (input) The order of the matrices A, B, VSL, and VSR.  N >= 0.


       A (input/output)
                 DOUBLE PRECISION array, dimension(LDA,A) On entry, the  first
                 of  the pair of matrices.  On exit, A has been overwritten by
                 its generalized Schur form S.


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


       B (input/output)
                 DOUBLE PRECISION array, dimension (LDB, N) On entry, the sec-
                 ond of the pair of matrices.  On exit, B has been overwritten
                 by its generalized Schur form T.


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


       SDIM (output)
                 If SORT = 'N', SDIM = 0.  If SORT = 'S', SDIM = number of ei-
                 genvalues (after sorting) for which DELCTG is true.  (Complex
                 conjugate pairs for which DELCTG is true for either eigenval-
                 ue count as 2.)


       ALPHAR (output)
                 On  exit,  (ALPHAR(j) + ALPHAI(j)*i)/BETA(j), j=1,...,N, will
                 be the generalized eigenvalues.  ALPHAR(j) + ALPHAI(j)*i, and
                 BETA(j),j=1,...,N are the diagonals of the complex Schur form
                 (S,T) that would result if the 2-by-2 diagonal blocks of  the
                 real  Schur  form of (A,B) were further reduced to triangular
                 form  using  2-by-2  complex  unitary  transformations.    If
                 ALPHAI(j) is zero, then the j-th eigenvalue is real; if posi-
                 tive, then the j-th and (j+1)-st eigenvalues  are  a  complex
                 conjugate pair, with ALPHAI(j+1) negative.

                 Note:  the  quotients ALPHAR(j)/BETA(j) and ALPHAI(j)/BETA(j)
                 may easily over- or underflow, and BETA(j) may even be  zero.
                 Thus,  the  user  should  avoid  naively computing the ratio.
                 However, ALPHAR and ALPHAI will be always less than and  usu-
                 ally  comparable  with  norm(A) in magnitude, and BETA always
                 less than and usually comparable with norm(B).


       ALPHAI (output)
                 See the description for ALPHAR.


       BETA (output)
                 See the description for ALPHAR.


       VSL (output)
                 If JOBVSL = 'V', VSL will contain  the  left  Schur  vectors.
                 Not referenced if JOBVSL = 'N'.


       LDVSL (input)
                 The  leading  dimension  of the matrix VSL. LDVSL >=1, and if
                 JOBVSL = 'V', LDVSL >= N.


       VSR (output)
                 If JOBVSR = 'V', VSR will contain the  right  Schur  vectors.
                 Not referenced if JOBVSR = 'N'.


       LDVSR (input)
                 The  leading  dimension of the matrix VSR. LDVSR >= 1, and if
                 JOBVSR = 'V', LDVSR >= N.


       WORK (workspace)
                 On exit, if INFO = 0, WORK(1) returns the optimal LWORK.


       LWORK (input)
                 The dimension of the array WORK.  LWORK >= max(8*N,6*N+16).

                 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.


       BWORK (workspace)
                 dimension(N) Not referenced if SORT = 'N'.


       INFO (output)
                 = 0:  successful exit
                 < 0:  if INFO = -i, the i-th argument had an illegal value.
                 =  1,...,N:  The QZ iteration failed.  (A,B) are not in Schur
                 form, but ALPHAR(j), ALPHAI(j), and BETA(j) should be correct
                 for  j=INFO+1,...,N.   >  N:   =N+1:  other than QZ iteration
                 failed in DHGEQZ.
                 =N+2: after reordering, roundoff changed values of some  com-
                 plex  eigenvalues so that leading eigenvalues in the General-
                 ized Schur form no longer satisfy DELCTG=.TRUE.   This  could
                 also  be  caused  due to scaling.  =N+3: reordering failed in
                 DTGSEN.




                                  7 Nov 2015                         dgges(3P)