LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
 subroutine zggsvd ( character JOBU, character JOBV, character JOBQ, integer M, integer N, integer P, integer K, integer L, complex*16, dimension( lda, * ) A, integer LDA, complex*16, dimension( ldb, * ) B, integer LDB, double precision, dimension( * ) ALPHA, double precision, dimension( * ) BETA, complex*16, dimension( ldu, * ) U, integer LDU, complex*16, dimension( ldv, * ) V, integer LDV, complex*16, dimension( ldq, * ) Q, integer LDQ, complex*16, dimension( * ) WORK, double precision, dimension( * ) RWORK, integer, dimension( * ) IWORK, integer INFO )

ZGGSVD computes the singular value decomposition (SVD) for OTHER matrices

Purpose:
``` This routine is deprecated and has been replaced by routine ZGGSVD3.

ZGGSVD computes the generalized singular value decomposition (GSVD)
of an M-by-N complex matrix A and P-by-N complex matrix B:

U**H*A*Q = D1*( 0 R ),    V**H*B*Q = D2*( 0 R )

where U, V and Q are unitary matrices.
Let K+L = the effective numerical rank of the
matrix (A**H,B**H)**H, then R is a (K+L)-by-(K+L) nonsingular upper
triangular matrix, D1 and D2 are M-by-(K+L) and P-by-(K+L) "diagonal"
matrices and of the following structures, respectively:

If M-K-L >= 0,

K  L
D1 =     K ( I  0 )
L ( 0  C )
M-K-L ( 0  0 )

K  L
D2 =   L ( 0  S )
P-L ( 0  0 )

N-K-L  K    L
( 0 R ) = K (  0   R11  R12 )
L (  0    0   R22 )
where

C = diag( ALPHA(K+1), ... , ALPHA(K+L) ),
S = diag( BETA(K+1),  ... , BETA(K+L) ),
C**2 + S**2 = I.

R is stored in A(1:K+L,N-K-L+1:N) on exit.

If M-K-L < 0,

K M-K K+L-M
D1 =   K ( I  0    0   )
M-K ( 0  C    0   )

K M-K K+L-M
D2 =   M-K ( 0  S    0  )
K+L-M ( 0  0    I  )
P-L ( 0  0    0  )

N-K-L  K   M-K  K+L-M
( 0 R ) =     K ( 0    R11  R12  R13  )
M-K ( 0     0   R22  R23  )
K+L-M ( 0     0    0   R33  )

where

C = diag( ALPHA(K+1), ... , ALPHA(M) ),
S = diag( BETA(K+1),  ... , BETA(M) ),
C**2 + S**2 = I.

(R11 R12 R13 ) is stored in A(1:M, N-K-L+1:N), and R33 is stored
( 0  R22 R23 )
in B(M-K+1:L,N+M-K-L+1:N) on exit.

The routine computes C, S, R, and optionally the unitary
transformation matrices U, V and Q.

In particular, if B is an N-by-N nonsingular matrix, then the GSVD of
A and B implicitly gives the SVD of A*inv(B):
A*inv(B) = U*(D1*inv(D2))*V**H.
If ( A**H,B**H)**H has orthnormal columns, then the GSVD of A and B is also
equal to the CS decomposition of A and B. Furthermore, the GSVD can
be used to derive the solution of the eigenvalue problem:
A**H*A x = lambda* B**H*B x.
In some literature, the GSVD of A and B is presented in the form
U**H*A*X = ( 0 D1 ),   V**H*B*X = ( 0 D2 )
where U and V are orthogonal and X is nonsingular, and D1 and D2 are
``diagonal''.  The former GSVD form can be converted to the latter
form by taking the nonsingular matrix X as

X = Q*(  I   0    )
(  0 inv(R) )```
Parameters
 [in] JOBU ``` JOBU is CHARACTER*1 = 'U': Unitary matrix U is computed; = 'N': U is not computed.``` [in] JOBV ``` JOBV is CHARACTER*1 = 'V': Unitary matrix V is computed; = 'N': V is not computed.``` [in] JOBQ ``` JOBQ is CHARACTER*1 = 'Q': Unitary matrix Q is computed; = 'N': Q is not computed.``` [in] M ``` M is INTEGER The number of rows of the matrix A. M >= 0.``` [in] N ``` N is INTEGER The number of columns of the matrices A and B. N >= 0.``` [in] P ``` P is INTEGER The number of rows of the matrix B. P >= 0.``` [out] K ` K is INTEGER` [out] L ``` L is INTEGER On exit, K and L specify the dimension of the subblocks described in Purpose. K + L = effective numerical rank of (A**H,B**H)**H.``` [in,out] A ``` A is COMPLEX*16 array, dimension (LDA,N) On entry, the M-by-N matrix A. On exit, A contains the triangular matrix R, or part of R. See Purpose for details.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,M).``` [in,out] B ``` B is COMPLEX*16 array, dimension (LDB,N) On entry, the P-by-N matrix B. On exit, B contains part of the triangular matrix R if M-K-L < 0. See Purpose for details.``` [in] LDB ``` LDB is INTEGER The leading dimension of the array B. LDB >= max(1,P).``` [out] ALPHA ` ALPHA is DOUBLE PRECISION array, dimension (N)` [out] BETA ``` BETA is DOUBLE PRECISION array, dimension (N) On exit, ALPHA and BETA contain the generalized singular value pairs of A and B; ALPHA(1:K) = 1, BETA(1:K) = 0, and if M-K-L >= 0, ALPHA(K+1:K+L) = C, BETA(K+1:K+L) = S, or if M-K-L < 0, ALPHA(K+1:M)=C, ALPHA(M+1:K+L)=0 BETA(K+1:M) =S, BETA(M+1:K+L) =1 and ALPHA(K+L+1:N) = 0 BETA(K+L+1:N) = 0``` [out] U ``` U is COMPLEX*16 array, dimension (LDU,M) If JOBU = 'U', U contains the M-by-M unitary matrix U. If JOBU = 'N', U is not referenced.``` [in] LDU ``` LDU is INTEGER The leading dimension of the array U. LDU >= max(1,M) if JOBU = 'U'; LDU >= 1 otherwise.``` [out] V ``` V is COMPLEX*16 array, dimension (LDV,P) If JOBV = 'V', V contains the P-by-P unitary matrix V. If JOBV = 'N', V is not referenced.``` [in] LDV ``` LDV is INTEGER The leading dimension of the array V. LDV >= max(1,P) if JOBV = 'V'; LDV >= 1 otherwise.``` [out] Q ``` Q is COMPLEX*16 array, dimension (LDQ,N) If JOBQ = 'Q', Q contains the N-by-N unitary matrix Q. If JOBQ = 'N', Q is not referenced.``` [in] LDQ ``` LDQ is INTEGER The leading dimension of the array Q. LDQ >= max(1,N) if JOBQ = 'Q'; LDQ >= 1 otherwise.``` [out] WORK ` WORK is COMPLEX*16 array, dimension (max(3*N,M,P)+N)` [out] RWORK ` RWORK is DOUBLE PRECISION array, dimension (2*N)` [out] IWORK ``` IWORK is INTEGER array, dimension (N) On exit, IWORK stores the sorting information. More precisely, the following loop will sort ALPHA for I = K+1, min(M,K+L) swap ALPHA(I) and ALPHA(IWORK(I)) endfor such that ALPHA(1) >= ALPHA(2) >= ... >= ALPHA(N).``` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value. > 0: if INFO = 1, the Jacobi-type procedure failed to converge. For further details, see subroutine ZTGSJA.```
Internal Parameters:
```  TOLA    DOUBLE PRECISION
TOLB    DOUBLE PRECISION
TOLA and TOLB are the thresholds to determine the effective
rank of (A**H,B**H)**H. Generally, they are set to
TOLA = MAX(M,N)*norm(A)*MAZHEPS,
TOLB = MAX(P,N)*norm(B)*MAZHEPS.
The size of TOLA and TOLB may affect the size of backward
errors of the decomposition.```
Date
November 2011
Contributors:
Ming Gu and Huan Ren, Computer Science Division, University of California at Berkeley, USA

Definition at line 339 of file zggsvd.f.

339 *
340 * -- LAPACK driver routine (version 3.4.0) --
341 * -- LAPACK is a software package provided by Univ. of Tennessee, --
342 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
343 * November 2011
344 *
345 * .. Scalar Arguments ..
346  CHARACTER jobq, jobu, jobv
347  INTEGER info, k, l, lda, ldb, ldq, ldu, ldv, m, n, p
348 * ..
349 * .. Array Arguments ..
350  INTEGER iwork( * )
351  DOUBLE PRECISION alpha( * ), beta( * ), rwork( * )
352  COMPLEX*16 a( lda, * ), b( ldb, * ), q( ldq, * ),
353  \$ u( ldu, * ), v( ldv, * ), work( * )
354 * ..
355 *
356 * =====================================================================
357 *
358 * .. Local Scalars ..
359  LOGICAL wantq, wantu, wantv
360  INTEGER i, ibnd, isub, j, ncycle
361  DOUBLE PRECISION anorm, bnorm, smax, temp, tola, tolb, ulp, unfl
362 * ..
363 * .. External Functions ..
364  LOGICAL lsame
365  DOUBLE PRECISION dlamch, zlange
366  EXTERNAL lsame, dlamch, zlange
367 * ..
368 * .. External Subroutines ..
369  EXTERNAL dcopy, xerbla, zggsvp, ztgsja
370 * ..
371 * .. Intrinsic Functions ..
372  INTRINSIC max, min
373 * ..
374 * .. Executable Statements ..
375 *
376 * Decode and test the input parameters
377 *
378  wantu = lsame( jobu, 'U' )
379  wantv = lsame( jobv, 'V' )
380  wantq = lsame( jobq, 'Q' )
381 *
382  info = 0
383  IF( .NOT.( wantu .OR. lsame( jobu, 'N' ) ) ) THEN
384  info = -1
385  ELSE IF( .NOT.( wantv .OR. lsame( jobv, 'N' ) ) ) THEN
386  info = -2
387  ELSE IF( .NOT.( wantq .OR. lsame( jobq, 'N' ) ) ) THEN
388  info = -3
389  ELSE IF( m.LT.0 ) THEN
390  info = -4
391  ELSE IF( n.LT.0 ) THEN
392  info = -5
393  ELSE IF( p.LT.0 ) THEN
394  info = -6
395  ELSE IF( lda.LT.max( 1, m ) ) THEN
396  info = -10
397  ELSE IF( ldb.LT.max( 1, p ) ) THEN
398  info = -12
399  ELSE IF( ldu.LT.1 .OR. ( wantu .AND. ldu.LT.m ) ) THEN
400  info = -16
401  ELSE IF( ldv.LT.1 .OR. ( wantv .AND. ldv.LT.p ) ) THEN
402  info = -18
403  ELSE IF( ldq.LT.1 .OR. ( wantq .AND. ldq.LT.n ) ) THEN
404  info = -20
405  END IF
406  IF( info.NE.0 ) THEN
407  CALL xerbla( 'ZGGSVD', -info )
408  RETURN
409  END IF
410 *
411 * Compute the Frobenius norm of matrices A and B
412 *
413  anorm = zlange( '1', m, n, a, lda, rwork )
414  bnorm = zlange( '1', p, n, b, ldb, rwork )
415 *
416 * Get machine precision and set up threshold for determining
417 * the effective numerical rank of the matrices A and B.
418 *
419  ulp = dlamch( 'Precision' )
420  unfl = dlamch( 'Safe Minimum' )
421  tola = max( m, n )*max( anorm, unfl )*ulp
422  tolb = max( p, n )*max( bnorm, unfl )*ulp
423 *
424  CALL zggsvp( jobu, jobv, jobq, m, p, n, a, lda, b, ldb, tola,
425  \$ tolb, k, l, u, ldu, v, ldv, q, ldq, iwork, rwork,
426  \$ work, work( n+1 ), info )
427 *
428 * Compute the GSVD of two upper "triangular" matrices
429 *
430  CALL ztgsja( jobu, jobv, jobq, m, p, n, k, l, a, lda, b, ldb,
431  \$ tola, tolb, alpha, beta, u, ldu, v, ldv, q, ldq,
432  \$ work, ncycle, info )
433 *
434 * Sort the singular values and store the pivot indices in IWORK
435 * Copy ALPHA to RWORK, then sort ALPHA in RWORK
436 *
437  CALL dcopy( n, alpha, 1, rwork, 1 )
438  ibnd = min( l, m-k )
439  DO 20 i = 1, ibnd
440 *
441 * Scan for largest ALPHA(K+I)
442 *
443  isub = i
444  smax = rwork( k+i )
445  DO 10 j = i + 1, ibnd
446  temp = rwork( k+j )
447  IF( temp.GT.smax ) THEN
448  isub = j
449  smax = temp
450  END IF
451  10 CONTINUE
452  IF( isub.NE.i ) THEN
453  rwork( k+isub ) = rwork( k+i )
454  rwork( k+i ) = smax
455  iwork( k+i ) = k + isub
456  ELSE
457  iwork( k+i ) = k + i
458  END IF
459  20 CONTINUE
460 *
461  RETURN
462 *
463 * End of ZGGSVD
464 *
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zggsvp(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)
ZGGSVP
Definition: zggsvp.f:267
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: zlange.f:117
subroutine ztgsja(JOBU, JOBV, JOBQ, M, P, N, K, L, A, LDA, B, LDB, TOLA, TOLB, ALPHA, BETA, U, LDU, V, LDV, Q, LDQ, WORK, NCYCLE, INFO)
ZTGSJA
Definition: ztgsja.f:381
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function: