262 SUBROUTINE zggsvp( JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB,
263 $ tola, tolb, k, l, u, ldu, v, ldv, q, ldq,
264 $ iwork, rwork, tau, work, info )
272 CHARACTER JOBQ, JOBU, JOBV
273 INTEGER INFO, K, L, LDA, LDB, LDQ, LDU, LDV, M, N, P
274 DOUBLE PRECISION TOLA, TOLB
278 DOUBLE PRECISION RWORK( * )
279 COMPLEX*16 A( lda, * ), B( ldb, * ), Q( ldq, * ),
280 $ tau( * ), u( ldu, * ), v( ldv, * ), work( * )
286 COMPLEX*16 CZERO, CONE
287 parameter( czero = ( 0.0d+0, 0.0d+0 ),
288 $ cone = ( 1.0d+0, 0.0d+0 ) )
291 LOGICAL FORWRD, WANTQ, WANTU, WANTV
304 INTRINSIC abs, dble, dimag, max, min
307 DOUBLE PRECISION CABS1
310 cabs1( t ) = abs( dble( t ) ) + abs( dimag( t ) )
316 wantu = lsame( jobu,
'U' )
317 wantv = lsame( jobv,
'V' )
318 wantq = lsame( jobq,
'Q' )
322 IF( .NOT.( wantu .OR. lsame( jobu,
'N' ) ) )
THEN
324 ELSE IF( .NOT.( wantv .OR. lsame( jobv,
'N' ) ) )
THEN
326 ELSE IF( .NOT.( wantq .OR. lsame( jobq,
'N' ) ) )
THEN
328 ELSE IF( m.LT.0 )
THEN
330 ELSE IF( p.LT.0 )
THEN
332 ELSE IF( n.LT.0 )
THEN
334 ELSE IF( lda.LT.max( 1, m ) )
THEN
336 ELSE IF( ldb.LT.max( 1, p ) )
THEN
338 ELSE IF( ldu.LT.1 .OR. ( wantu .AND. ldu.LT.m ) )
THEN
340 ELSE IF( ldv.LT.1 .OR. ( wantv .AND. ldv.LT.p ) )
THEN
342 ELSE IF( ldq.LT.1 .OR. ( wantq .AND. ldq.LT.n ) )
THEN
346 CALL
xerbla(
'ZGGSVP', -info )
356 CALL
zgeqpf( p, n, b, ldb, iwork, tau, work, rwork, info )
360 CALL
zlapmt( forwrd, m, n, a, lda, iwork )
365 DO 20 i = 1, min( p, n )
366 IF( cabs1( b( i, i ) ).GT.tolb )
374 CALL
zlaset(
'Full', p, p, czero, czero, v, ldv )
376 $ CALL
zlacpy(
'Lower', p-1, n, b( 2, 1 ), ldb, v( 2, 1 ),
378 CALL
zung2r( p, p, min( p, n ), v, ldv, tau, work, info )
389 $ CALL
zlaset(
'Full', p-l, n, czero, czero, b( l+1, 1 ), ldb )
395 CALL
zlaset(
'Full', n, n, czero, cone, q, ldq )
396 CALL
zlapmt( forwrd, n, n, q, ldq, iwork )
399 IF( p.GE.l .AND. n.NE.l )
THEN
403 CALL
zgerq2( l, n, b, ldb, tau, work, info )
407 CALL
zunmr2(
'Right',
'Conjugate transpose', m, n, l, b, ldb,
408 $ tau, a, lda, work, info )
413 CALL
zunmr2(
'Right',
'Conjugate transpose', n, n, l, b,
414 $ ldb, tau, q, ldq, work, info )
419 CALL
zlaset(
'Full', l, n-l, czero, czero, b, ldb )
420 DO 60 j = n - l + 1, n
421 DO 50 i = j - n + l + 1, l
439 CALL
zgeqpf( m, n-l, a, lda, iwork, tau, work, rwork, info )
444 DO 80 i = 1, min( m, n-l )
445 IF( cabs1( a( i, i ) ).GT.tola )
451 CALL
zunm2r(
'Left',
'Conjugate transpose', m, l, min( m, n-l ),
452 $ a, lda, tau, a( 1, n-l+1 ), lda, work, info )
458 CALL
zlaset(
'Full', m, m, czero, czero, u, ldu )
460 $ CALL
zlacpy(
'Lower', m-1, n-l, a( 2, 1 ), lda, u( 2, 1 ),
462 CALL
zung2r( m, m, min( m, n-l ), u, ldu, tau, work, info )
469 CALL
zlapmt( forwrd, n, n-l, q, ldq, iwork )
481 $ CALL
zlaset(
'Full', m-k, n-l, czero, czero, a( k+1, 1 ), lda )
487 CALL
zgerq2( k, n-l, a, lda, tau, work, info )
493 CALL
zunmr2(
'Right',
'Conjugate transpose', n, n-l, k, a,
494 $ lda, tau, q, ldq, work, info )
499 CALL
zlaset(
'Full', k, n-l-k, czero, czero, a, lda )
500 DO 120 j = n - l - k + 1, n - l
501 DO 110 i = j - n + l + k + 1, k
512 CALL
zgeqr2( m-k, l, a( k+1, n-l+1 ), lda, tau, work, info )
518 CALL
zunm2r(
'Right',
'No transpose', m, m-k, min( m-k, l ),
519 $ a( k+1, n-l+1 ), lda, tau, u( 1, k+1 ), ldu,
525 DO 140 j = n - l + 1, n
526 DO 130 i = j - n + k + l + 1, m
subroutine zgeqr2(M, N, A, LDA, TAU, WORK, INFO)
ZGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm...
subroutine zlapmt(FORWRD, M, N, X, LDX, K)
ZLAPMT performs a forward or backward permutation of the columns of a matrix.
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
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
subroutine zunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
subroutine zgeqpf(M, N, A, LDA, JPVT, TAU, WORK, RWORK, INFO)
ZGEQPF
subroutine zunmr2(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNMR2 multiplies a general matrix by the unitary matrix from a RQ factorization determined by cgerqf...
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine zgerq2(M, N, A, LDA, TAU, WORK, INFO)
ZGERQ2 computes the RQ factorization of a general rectangular matrix using an unblocked algorithm...
subroutine zung2r(M, N, K, A, LDA, TAU, WORK, INFO)
ZUNG2R