Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

dggev (3p)

Name

dggev - N real nonsymmetric matrices (A,B)

Synopsis

SUBROUTINE DGGEV(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI,
BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

CHARACTER*1 JOBVL, JOBVR
INTEGER N, LDA, LDB, LDVL, LDVR, LWORK, INFO
DOUBLE  PRECISION  A(LDA,*),  B(LDB,*),  ALPHAR(*), ALPHAI(*), BETA(*),
VL(LDVL,*), VR(LDVR,*), WORK(*)

SUBROUTINE DGGEV_64(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI,
BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

CHARACTER*1 JOBVL, JOBVR
INTEGER*8 N, LDA, LDB, LDVL, LDVR, LWORK, INFO
DOUBLE PRECISION A(LDA,*),  B(LDB,*),  ALPHAR(*),  ALPHAI(*),  BETA(*),
VL(LDVL,*), VR(LDVR,*), WORK(*)




F95 INTERFACE
SUBROUTINE GGEV(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

CHARACTER(LEN=1) :: JOBVL, JOBVR
INTEGER :: N, LDA, LDB, LDVL, LDVR, LWORK, INFO
REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
REAL(8), DIMENSION(:,:) :: A, B, VL, VR

SUBROUTINE GGEV_64(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

CHARACTER(LEN=1) :: JOBVL, JOBVR
INTEGER(8) :: N, LDA, LDB, LDVL, LDVR, LWORK, INFO
REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
REAL(8), DIMENSION(:,:) :: A, B, VL, VR




C INTERFACE
#include <sunperf.h>

void  dggev(char  jobvl,  char jobvr, int n, double *a, int lda, double
*b, int ldb, double *alphar, double  *alphai,  double  *beta,
double *vl, int ldvl, double *vr, int ldvr, int *info);

void dggev_64(char jobvl, char jobvr, long n, double *a, long lda, dou-
ble *b, long ldb,  double  *alphar,  double  *alphai,  double
*beta,  double  *vl,  long  ldvl, double *vr, long ldvr, long
*info);

Description

Oracle Solaris Studio Performance Library                            dggev(3P)



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


SYNOPSIS
       SUBROUTINE DGGEV(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI,
             BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

       CHARACTER*1 JOBVL, JOBVR
       INTEGER N, LDA, LDB, LDVL, LDVR, LWORK, INFO
       DOUBLE  PRECISION  A(LDA,*),  B(LDB,*),  ALPHAR(*), ALPHAI(*), BETA(*),
       VL(LDVL,*), VR(LDVR,*), WORK(*)

       SUBROUTINE DGGEV_64(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI,
             BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

       CHARACTER*1 JOBVL, JOBVR
       INTEGER*8 N, LDA, LDB, LDVL, LDVR, LWORK, INFO
       DOUBLE PRECISION A(LDA,*),  B(LDB,*),  ALPHAR(*),  ALPHAI(*),  BETA(*),
       VL(LDVL,*), VR(LDVR,*), WORK(*)




   F95 INTERFACE
       SUBROUTINE GGEV(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
              ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

       CHARACTER(LEN=1) :: JOBVL, JOBVR
       INTEGER :: N, LDA, LDB, LDVL, LDVR, LWORK, INFO
       REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
       REAL(8), DIMENSION(:,:) :: A, B, VL, VR

       SUBROUTINE GGEV_64(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
              ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)

       CHARACTER(LEN=1) :: JOBVL, JOBVR
       INTEGER(8) :: N, LDA, LDB, LDVL, LDVR, LWORK, INFO
       REAL(8), DIMENSION(:) :: ALPHAR, ALPHAI, BETA, WORK
       REAL(8), DIMENSION(:,:) :: A, B, VL, VR




   C INTERFACE
       #include <sunperf.h>

       void  dggev(char  jobvl,  char jobvr, int n, double *a, int lda, double
                 *b, int ldb, double *alphar, double  *alphai,  double  *beta,
                 double *vl, int ldvl, double *vr, int ldvr, int *info);

       void dggev_64(char jobvl, char jobvr, long n, double *a, long lda, dou-
                 ble *b, long ldb,  double  *alphar,  double  *alphai,  double
                 *beta,  double  *vl,  long  ldvl, double *vr, long ldvr, long
                 *info);



PURPOSE
       dggev computes for a pair of N-by-N real  nonsymmetric  matrices  (A,B)
       the generalized eigenvalues, and optionally, the left and/or right gen-
       eralized eigenvectors.

       A generalized eigenvalue for a pair  of  matrices  (A,B)  is  a  scalar
       lambda or a ratio alpha/beta = lambda, such that A - lambda*B is singu-
       lar. It is usually represented as the pair (alpha,beta), as there is  a
       reasonable interpretation for beta=0, and even for both being zero.

       The right eigenvector v(j) corresponding to the eigenvalue lambda(j) of
       (A,B) satisfies

                        A * v(j) = lambda(j) * B * v(j).

       The left eigenvector u(j) corresponding to the eigenvalue lambda(j)  of
       (A,B) satisfies

                        u(j)**H * A  = lambda(j) * u(j)**H * B .

       where u(j)**H is the conjugate-transpose of u(j).


ARGUMENTS
       JOBVL (input)
                 = 'N':  do not compute the left generalized eigenvectors;
                 = 'V':  compute the left generalized eigenvectors.


       JOBVR (input)
                 = 'N':  do not compute the right generalized eigenvectors;
                 = 'V':  compute the right generalized eigenvectors.


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


       A (input/output)
                 On  entry,  the  matrix  A in the pair (A,B).  On exit, A has
                 been overwritten.


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


       B (input/output)
                 On entry, the matrix B in the pair (A,B).   On  exit,  B  has
                 been overwritten.


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


       ALPHAR (output)
                 On  exit,  (ALPHAR(j) + ALPHAI(j)*i)/BETA(j), j=1,...,N, will
                 be the generalized eigenvalues.  If ALPHAI(j) is  zero,  then
                 the  j-th  eigenvalue is real; if positive, 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
                 alpha/beta.  However, ALPHAR and ALPHAI will be  always  less
                 than  and  usually  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.


       VL (output)
                 If JOBVL = 'V', the left eigenvectors  u(j)  are  stored  one
                 after  another  in  the  columns  of VL, in the same order as
                 their eigenvalues. If the j-th eigenvalue is real, then  u(j)
                 =  VL(:,j),  the  j-th column of VL. If the j-th and (j+1)-th
                 eigenvalues form  a  complex  conjugate  pair,  then  u(j)  =
                 VL(:,j)+i*VL(:,j+1)  and  u(j+1) = VL(:,j)-i*VL(:,j+1).  Each
                 eigenvector will be scaled  so  the  largest  component  have
                 abs(real  part)+abs(imag. part)=1.  Not referenced if JOBVL =
                 'N'.


       LDVL (input)
                 The leading dimension of the matrix VL. LDVL  >=  1,  and  if
                 JOBVL = 'V', LDVL >= N.


       VR (output)
                 If  JOBVR  =  'V', the right eigenvectors v(j) are stored one
                 after another in the columns of VR,  in  the  same  order  as
                 their  eigenvalues. If the j-th eigenvalue is real, then v(j)
                 = VR(:,j), the j-th column of VR. If the  j-th  and  (j+1)-th
                 eigenvalues  form  a  complex  conjugate  pair,  then  v(j) =
                 VR(:,j)+i*VR(:,j+1) and v(j+1) =  VR(:,j)-i*VR(:,j+1).   Each
                 eigenvector  will  be  scaled  so  the largest component have
                 abs(real part)+abs(imag. part)=1.  Not referenced if JOBVR  =
                 'N'.


       LDVR (input)
                 The  leading  dimension  of  the matrix VR. LDVR >= 1, and if
                 JOBVR = 'V', LDVR >= N.


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


       LWORK (input)
                 The dimension of the array WORK.  LWORK >=  max(1,8*N).   For
                 good performance, LWORK must generally be larger.

                 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.


       INFO (output)
                 = 0:  successful exit
                 < 0:  if INFO = -i, the i-th argument had an illegal value.
                 = 1,...,N: The QZ iteration  failed.   No  eigenvectors  have
                 been calculated, 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: error return from DTGEVC.




                                  7 Nov 2015                         dggev(3P)