377 SUBROUTINE stgsja( JOBU, JOBV, JOBQ, M, P, N, K, L, A, LDA, B,
378 $ ldb, tola, tolb, alpha, beta, u, ldu, v, ldv,
379 $ q, ldq, work, ncycle, info )
387 CHARACTER JOBQ, JOBU, JOBV
388 INTEGER INFO, K, L, LDA, LDB, LDQ, LDU, LDV, M, N,
393 REAL A( lda, * ), ALPHA( * ), B( ldb, * ),
394 $ beta( * ), q( ldq, * ), u( ldu, * ),
395 $ v( ldv, * ), work( * )
402 parameter( maxit = 40 )
404 parameter( zero = 0.0e+0, one = 1.0e+0 )
408 LOGICAL INITQ, INITU, INITV, UPPER, WANTQ, WANTU, WANTV
410 REAL A1, A2, A3, B1, B2, B3, CSQ, CSU, CSV, ERROR,
411 $ gamma, rwk, snq, snu, snv, ssmin
422 INTRINSIC abs, max, min
428 initu = lsame( jobu,
'I' )
429 wantu = initu .OR. lsame( jobu,
'U' )
431 initv = lsame( jobv,
'I' )
432 wantv = initv .OR. lsame( jobv,
'V' )
434 initq = lsame( jobq,
'I' )
435 wantq = initq .OR. lsame( jobq,
'Q' )
438 IF( .NOT.( initu .OR. wantu .OR. lsame( jobu,
'N' ) ) )
THEN
440 ELSE IF( .NOT.( initv .OR. wantv .OR. lsame( jobv,
'N' ) ) )
THEN
442 ELSE IF( .NOT.( initq .OR. wantq .OR. lsame( jobq,
'N' ) ) )
THEN
444 ELSE IF( m.LT.0 )
THEN
446 ELSE IF( p.LT.0 )
THEN
448 ELSE IF( n.LT.0 )
THEN
450 ELSE IF( lda.LT.max( 1, m ) )
THEN
452 ELSE IF( ldb.LT.max( 1, p ) )
THEN
454 ELSE IF( ldu.LT.1 .OR. ( wantu .AND. ldu.LT.m ) )
THEN
456 ELSE IF( ldv.LT.1 .OR. ( wantv .AND. ldv.LT.p ) )
THEN
458 ELSE IF( ldq.LT.1 .OR. ( wantq .AND. ldq.LT.n ) )
THEN
462 CALL
xerbla(
'STGSJA', -info )
469 $ CALL
slaset(
'Full', m, m, zero, one, u, ldu )
471 $ CALL
slaset(
'Full', p, p, zero, one, v, ldv )
473 $ CALL
slaset(
'Full', n, n, zero, one, q, ldq )
478 DO 40 kcycle = 1, maxit
489 $ a1 = a( k+i, n-l+i )
491 $ a3 = a( k+j, n-l+j )
498 $ a2 = a( k+i, n-l+j )
502 $ a2 = a( k+j, n-l+i )
506 CALL
slags2( upper, a1, a2, a3, b1, b2, b3, csu, snu,
507 $ csv, snv, csq, snq )
512 $ CALL
srot( l, a( k+j, n-l+1 ), lda, a( k+i, n-l+1 ),
517 CALL
srot( l, b( j, n-l+1 ), ldb, b( i, n-l+1 ), ldb,
523 CALL
srot( min( k+l, m ), a( 1, n-l+j ), 1,
524 $ a( 1, n-l+i ), 1, csq, snq )
526 CALL
srot( l, b( 1, n-l+j ), 1, b( 1, n-l+i ), 1, csq,
531 $ a( k+i, n-l+j ) = zero
535 $ a( k+j, n-l+i ) = zero
541 IF( wantu .AND. k+j.LE.m )
542 $ CALL
srot( m, u( 1, k+j ), 1, u( 1, k+i ), 1, csu,
546 $ CALL
srot( p, v( 1, j ), 1, v( 1, i ), 1, csv, snv )
549 $ CALL
srot( n, q( 1, n-l+j ), 1, q( 1, n-l+i ), 1, csq,
555 IF( .NOT.upper )
THEN
564 DO 30 i = 1, min( l, m-k )
565 CALL
scopy( l-i+1, a( k+i, n-l+i ), lda, work, 1 )
566 CALL
scopy( l-i+1, b( i, n-l+i ), ldb, work( l+1 ), 1 )
567 CALL
slapll( l-i+1, work, 1, work( l+1 ), 1, ssmin )
568 error = max( error, ssmin )
571 IF( abs( error ).LE.min( tola, tolb ) )
595 DO 70 i = 1, min( l, m-k )
600 IF( a1.NE.zero )
THEN
605 IF( gamma.LT.zero )
THEN
606 CALL
sscal( l-i+1, -one, b( i, n-l+i ), ldb )
608 $ CALL
sscal( p, -one, v( 1, i ), 1 )
611 CALL
slartg( abs( gamma ), one, beta( k+i ), alpha( k+i ),
614 IF( alpha( k+i ).GE.beta( k+i ) )
THEN
615 CALL
sscal( l-i+1, one / alpha( k+i ), a( k+i, n-l+i ),
618 CALL
sscal( l-i+1, one / beta( k+i ), b( i, n-l+i ),
620 CALL
scopy( l-i+1, b( i, n-l+i ), ldb, a( k+i, n-l+i ),
628 CALL
scopy( l-i+1, b( i, n-l+i ), ldb, a( k+i, n-l+i ),
637 DO 80 i = m + 1, k + l
643 DO 90 i = k + l + 1, n
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine slags2(UPPER, A1, A2, A3, B1, B2, B3, CSU, SNU, CSV, SNV, CSQ, SNQ)
SLAGS2 computes 2-by-2 orthogonal matrices U, V, and Q, and applies them to matrices A and B such tha...
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine stgsja(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)
STGSJA
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
subroutine slapll(N, X, INCX, Y, INCY, SSMIN)
SLAPLL measures the linear dependence of two vectors.
subroutine slartg(F, G, CS, SN, R)
SLARTG generates a plane rotation with real cosine and real sine.
subroutine sscal(N, SA, SX, INCX)
SSCAL
subroutine srot(N, SX, INCX, SY, INCY, C, S)
SROT