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