Go to main content
Oracle Developer Studio 12.5 Man Pages

Exit Print View

Updated: June 2017
 
 

cggsvp (3p)

Name

cggsvp - compute unitary matrices

Synopsis

SUBROUTINE CGGSVP(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB, TOLA,
TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK, TAU, WORK,
INFO)

CHARACTER*1 JOBU, JOBV, JOBQ
COMPLEX  A(LDA,*),  B(LDB,*),  U(LDU,*),  V(LDV,*),  Q(LDQ,*),  TAU(*),
WORK(*)
INTEGER M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
INTEGER IWORK(*)
REAL TOLA, TOLB
REAL RWORK(*)

SUBROUTINE CGGSVP_64(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB, TOLA,
TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK, TAU, WORK,
INFO)

CHARACTER*1 JOBU, JOBV, JOBQ
COMPLEX  A(LDA,*),  B(LDB,*),  U(LDU,*),  V(LDV,*),  Q(LDQ,*),  TAU(*),
WORK(*)
INTEGER*8 M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
INTEGER*8 IWORK(*)
REAL TOLA, TOLB
REAL RWORK(*)




F95 INTERFACE
SUBROUTINE GGSVP(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB,
TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK,
TAU, WORK, INFO)

CHARACTER(LEN=1) :: JOBU, JOBV, JOBQ
COMPLEX, DIMENSION(:) :: TAU, WORK
COMPLEX, DIMENSION(:,:) :: A, B, U, V, Q
INTEGER :: M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
INTEGER, DIMENSION(:) :: IWORK
REAL :: TOLA, TOLB
REAL, DIMENSION(:) :: RWORK

SUBROUTINE GGSVP_64(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B,
LDB, TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK,
RWORK, TAU, WORK, INFO)

CHARACTER(LEN=1) :: JOBU, JOBV, JOBQ
COMPLEX, DIMENSION(:) :: TAU, WORK
COMPLEX, DIMENSION(:,:) :: A, B, U, V, Q
INTEGER(8) :: M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
INTEGER(8), DIMENSION(:) :: IWORK
REAL :: TOLA, TOLB
REAL, DIMENSION(:) :: RWORK




C INTERFACE
#include <sunperf.h>

void  cggsvp(char jobu, char jobv, char jobq, int m, int p, int n, com-
plex *a, int lda, complex *b,  int  ldb,  float  tola,  float
tolb,  int  *k,  int *l, complex *u, int ldu, complex *v, int
ldv, complex *q, int ldq, int *info);

void cggsvp_64(char jobu, char jobv, char jobq, long m, long p, long n,
complex *a, long lda, complex *b, long ldb, float tola, float
tolb, long *k, long *l, complex *u,  long  ldu,  complex  *v,
long ldv, complex *q, long ldq, long *info);

Description

Oracle Solaris Studio Performance Library                           cggsvp(3P)



NAME
       cggsvp - compute unitary matrices


SYNOPSIS
       SUBROUTINE CGGSVP(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB, TOLA,
             TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK, TAU, WORK,
             INFO)

       CHARACTER*1 JOBU, JOBV, JOBQ
       COMPLEX  A(LDA,*),  B(LDB,*),  U(LDU,*),  V(LDV,*),  Q(LDQ,*),  TAU(*),
       WORK(*)
       INTEGER M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
       INTEGER IWORK(*)
       REAL TOLA, TOLB
       REAL RWORK(*)

       SUBROUTINE CGGSVP_64(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB, TOLA,
             TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK, TAU, WORK,
             INFO)

       CHARACTER*1 JOBU, JOBV, JOBQ
       COMPLEX  A(LDA,*),  B(LDB,*),  U(LDU,*),  V(LDV,*),  Q(LDQ,*),  TAU(*),
       WORK(*)
       INTEGER*8 M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
       INTEGER*8 IWORK(*)
       REAL TOLA, TOLB
       REAL RWORK(*)




   F95 INTERFACE
       SUBROUTINE GGSVP(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB,
              TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK, RWORK,
              TAU, WORK, INFO)

       CHARACTER(LEN=1) :: JOBU, JOBV, JOBQ
       COMPLEX, DIMENSION(:) :: TAU, WORK
       COMPLEX, DIMENSION(:,:) :: A, B, U, V, Q
       INTEGER :: M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
       INTEGER, DIMENSION(:) :: IWORK
       REAL :: TOLA, TOLB
       REAL, DIMENSION(:) :: RWORK

       SUBROUTINE GGSVP_64(JOBU, JOBV, JOBQ, M, P, N, A, LDA, B,
              LDB, TOLA, TOLB, K, L, U, LDU, V, LDV, Q, LDQ, IWORK,
              RWORK, TAU, WORK, INFO)

       CHARACTER(LEN=1) :: JOBU, JOBV, JOBQ
       COMPLEX, DIMENSION(:) :: TAU, WORK
       COMPLEX, DIMENSION(:,:) :: A, B, U, V, Q
       INTEGER(8) :: M, P, N, LDA, LDB, K, L, LDU, LDV, LDQ, INFO
       INTEGER(8), DIMENSION(:) :: IWORK
       REAL :: TOLA, TOLB
       REAL, DIMENSION(:) :: RWORK




   C INTERFACE
       #include <sunperf.h>

       void  cggsvp(char jobu, char jobv, char jobq, int m, int p, int n, com-
                 plex *a, int lda, complex *b,  int  ldb,  float  tola,  float
                 tolb,  int  *k,  int *l, complex *u, int ldu, complex *v, int
                 ldv, complex *q, int ldq, int *info);

       void cggsvp_64(char jobu, char jobv, char jobq, long m, long p, long n,
                 complex *a, long lda, complex *b, long ldb, float tola, float
                 tolb, long *k, long *l, complex *u,  long  ldu,  complex  *v,
                 long ldv, complex *q, long ldq, long *info);



PURPOSE
       cggsvp computes unitary matrices U, V and Q such that
                     L ( 0     0   A23 )
                 M-K-L ( 0     0    0  )

                        N-K-L  K    L
               =     K ( 0    A12  A13 )  if M-K-L < 0;
                   M-K ( 0     0   A23 )

                      N-K-L  K    L
        V'*B*Q =   L ( 0     0   B13 )
                 P-L ( 0     0    0  )

       where the K-by-K matrix A12 and L-by-L matrix B13 are nonsingular upper
       triangular; A23 is L-by-L upper triangular if M-K-L >= 0, otherwise A23
       is  (M-K)-by-L upper trapezoidal. K+L = the effective numerical rank of
       the (M+P)-by-N matrix (A',B')'. Z' denotes the conjugate  transpose  of
       Z.

       This decomposition is the preprocessing step for computing the General-
       ized Singular Value Decomposition (GSVD), see subroutine CGGSVD.


ARGUMENTS
       JOBU (input)
                 = 'U':  Unitary matrix U is computed;
                 = 'N':  U is not computed.


       JOBV (input)
                 = 'V':  Unitary matrix V is computed;
                 = 'N':  V is not computed.


       JOBQ (input)
                 = 'Q':  Unitary matrix Q is computed;
                 = 'N':  Q is not computed.


       M (input) The number of rows of the matrix A.  M >= 0.


       P (input) The number of rows of the matrix B.  P >= 0.


       N (input) The number of columns of the matrices A and B.  N >= 0.


       A (input/output)
                 On entry, the M-by-N matrix A.  On exit, A contains the  tri-
                 angular (or trapezoidal) matrix described in the Purpose sec-
                 tion.


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


       B (input/output)
                 On entry, the P-by-N matrix B.
                 On exit, B contains the triangular matrix  described  in  the
                 Purpose section.


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


       TOLA (input)
                 TOLA  and  TOLB are the thresholds to determine the effective
                 numerical rank of matrix B and a subblock  of  A.  Generally,
                 they are set to
                 TOLA = MAX(M,N)*norm(A)*MACHEPS,
                 TOLB = MAX(P,N)*norm(B)*MACHEPS.
                 The  size  of  TOLA  and TOLB may affect the size of backward
                 errors of the decomposition.


       TOLB (input)
                 See description of TOLA.


       K (output)
                 On exit, K and L  specify  the  dimension  of  the  subblocks
                 described in Purpose section.
                 K + L = effective numerical rank of (A',B')'.


       L (output)
                 See the description of K.


       U (output)
                 If JOBU = 'U', U contains the unitary matrix U.
                 If JOBU = 'N', U is not referenced.


       LDU (input)
                 The leading dimension of the array U.
                 LDU >= max(1,M) if JOBU = 'U'; LDU >= 1 otherwise.


       V (output)
                 If JOBV = 'V', V contains the unitary matrix V.
                 If JOBV = 'N', V is not referenced.


       LDV (input)
                 The leading dimension of the array V.
                 LDV >= max(1,P) if JOBV = 'V'; LDV >= 1 otherwise.


       Q (output)
                 If JOBQ = 'Q', Q contains the unitary matrix Q.
                 If JOBQ = 'N', Q is not referenced.


       LDQ (input)
                 The leading dimension of the array Q.
                 LDQ >= max(1,N) if JOBQ = 'Q'; LDQ >= 1 otherwise.


       IWORK (workspace)
                 dimension(N)

       RWORK (workspace)
                 dimension(2*N)

       TAU (workspace)
                 dimension(N)

       WORK (workspace)
                 dimension(MAX(3*N,M,P))

       INFO (output)
                 = 0:  successful exit;
                 < 0:  if INFO = -i, the i-th argument had an illegal value.

FURTHER DETAILS
       The  subroutine  uses LAPACK subroutine CGEQPF for the QR factorization
       with column pivoting to detect the effective numerical rank  of  the  a
       matrix. It may be replaced by a better rank determination strategy.




                                  7 Nov 2015                        cggsvp(3P)