261 SUBROUTINE cggsvp( JOBU, JOBV, JOBQ, M, P, N, A, LDA, B, LDB,
262 $ tola, tolb, k, l, u, ldu, v, ldv, q, ldq,
263 $ iwork, rwork, tau, work, info )
271 CHARACTER JOBQ, JOBU, JOBV
272 INTEGER INFO, K, L, LDA, LDB, LDQ, LDU, LDV, M, N, P
278 COMPLEX A( lda, * ), B( ldb, * ), Q( ldq, * ),
279 $ tau( * ), u( ldu, * ), v( ldv, * ), work( * )
286 parameter ( czero = ( 0.0e+0, 0.0e+0 ),
287 $ cone = ( 1.0e+0, 0.0e+0 ) )
290 LOGICAL FORWRD, WANTQ, WANTU, WANTV
303 INTRINSIC abs, aimag, max, min, real
309 cabs1( t ) = abs(
REAL( T ) ) + abs( AIMAG( t ) )
315 wantu = lsame( jobu,
'U' )
316 wantv = lsame( jobv,
'V' )
317 wantq = lsame( jobq,
'Q' )
321 IF( .NOT.( wantu .OR. lsame( jobu,
'N' ) ) )
THEN
323 ELSE IF( .NOT.( wantv .OR. lsame( jobv,
'N' ) ) )
THEN
325 ELSE IF( .NOT.( wantq .OR. lsame( jobq,
'N' ) ) )
THEN
327 ELSE IF( m.LT.0 )
THEN
329 ELSE IF( p.LT.0 )
THEN
331 ELSE IF( n.LT.0 )
THEN
333 ELSE IF( lda.LT.max( 1, m ) )
THEN
335 ELSE IF( ldb.LT.max( 1, p ) )
THEN
337 ELSE IF( ldu.LT.1 .OR. ( wantu .AND. ldu.LT.m ) )
THEN
339 ELSE IF( ldv.LT.1 .OR. ( wantv .AND. ldv.LT.p ) )
THEN
341 ELSE IF( ldq.LT.1 .OR. ( wantq .AND. ldq.LT.n ) )
THEN
345 CALL xerbla(
'CGGSVP', -info )
355 CALL cgeqpf( p, n, b, ldb, iwork, tau, work, rwork, info )
359 CALL clapmt( forwrd, m, n, a, lda, iwork )
364 DO 20 i = 1, min( p, n )
365 IF( cabs1( b( i, i ) ).GT.tolb )
373 CALL claset(
'Full', p, p, czero, czero, v, ldv )
375 $
CALL clacpy(
'Lower', p-1, n, b( 2, 1 ), ldb, v( 2, 1 ),
377 CALL cung2r( p, p, min( p, n ), v, ldv, tau, work, info )
388 $
CALL claset(
'Full', p-l, n, czero, czero, b( l+1, 1 ), ldb )
394 CALL claset(
'Full', n, n, czero, cone, q, ldq )
395 CALL clapmt( forwrd, n, n, q, ldq, iwork )
398 IF( p.GE.l .AND. n.NE.l )
THEN
402 CALL cgerq2( l, n, b, ldb, tau, work, info )
406 CALL cunmr2(
'Right',
'Conjugate transpose', m, n, l, b, ldb,
407 $ tau, a, lda, work, info )
412 CALL cunmr2(
'Right',
'Conjugate transpose', n, n, l, b,
413 $ ldb, tau, q, ldq, work, info )
418 CALL claset(
'Full', l, n-l, czero, czero, b, ldb )
419 DO 60 j = n - l + 1, n
420 DO 50 i = j - n + l + 1, l
438 CALL cgeqpf( m, n-l, a, lda, iwork, tau, work, rwork, info )
443 DO 80 i = 1, min( m, n-l )
444 IF( cabs1( a( i, i ) ).GT.tola )
450 CALL cunm2r(
'Left',
'Conjugate transpose', m, l, min( m, n-l ),
451 $ a, lda, tau, a( 1, n-l+1 ), lda, work, info )
457 CALL claset(
'Full', m, m, czero, czero, u, ldu )
459 $
CALL clacpy(
'Lower', m-1, n-l, a( 2, 1 ), lda, u( 2, 1 ),
461 CALL cung2r( m, m, min( m, n-l ), u, ldu, tau, work, info )
468 CALL clapmt( forwrd, n, n-l, q, ldq, iwork )
480 $
CALL claset(
'Full', m-k, n-l, czero, czero, a( k+1, 1 ), lda )
486 CALL cgerq2( k, n-l, a, lda, tau, work, info )
492 CALL cunmr2(
'Right',
'Conjugate transpose', n, n-l, k, a,
493 $ lda, tau, q, ldq, work, info )
498 CALL claset(
'Full', k, n-l-k, czero, czero, a, lda )
499 DO 120 j = n - l - k + 1, n - l
500 DO 110 i = j - n + l + k + 1, k
511 CALL cgeqr2( m-k, l, a( k+1, n-l+1 ), lda, tau, work, info )
517 CALL cunm2r(
'Right',
'No transpose', m, m-k, min( m-k, l ),
518 $ a( k+1, n-l+1 ), lda, tau, u( 1, k+1 ), ldu,
524 DO 140 j = n - l + 1, n
525 DO 130 i = j - n + k + l + 1, m
subroutine cgeqr2(M, N, A, LDA, TAU, WORK, INFO)
CGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm...
subroutine cunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
CUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
subroutine cung2r(M, N, K, A, LDA, TAU, WORK, INFO)
CUNG2R
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cgeqpf(M, N, A, LDA, JPVT, TAU, WORK, RWORK, INFO)
CGEQPF
subroutine cunmr2(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
CUNMR2 multiplies a general matrix by the unitary matrix from a RQ factorization determined by cgerqf...
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine cgerq2(M, N, A, LDA, TAU, WORK, INFO)
CGERQ2 computes the RQ factorization of a general rectangular matrix using an unblocked algorithm...
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)
CGGSVP
subroutine clapmt(FORWRD, M, N, X, LDX, K)
CLAPMT performs a forward or backward permutation of the columns of a matrix.