266 SUBROUTINE cgges3( JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B,
267 $ LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR,
268 $ WORK, LWORK, RWORK, BWORK, INFO )
275 CHARACTER JOBVSL, JOBVSR, SORT
276 INTEGER INFO, LDA, LDB, LDVSL, LDVSR, LWORK, N, SDIM
281 COMPLEX A( LDA, * ), ALPHA( * ), B( LDB, * ),
282 $ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
294 PARAMETER ( ZERO = 0.0e0, one = 1.0e0 )
296 parameter( czero = ( 0.0e0, 0.0e0 ),
297 $ cone = ( 1.0e0, 0.0e0 ) )
300 LOGICAL CURSL, ILASCL, ILBSCL, ILVSL, ILVSR, LASTSL,
302 INTEGER I, ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT,
303 $ ILO, IRIGHT, IROWS, IRWRK, ITAU, IWRK, LWKOPT
304 REAL ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS, PVSL,
319 EXTERNAL lsame, clange, slamch
328 IF( lsame( jobvsl,
'N' ) )
THEN
331 ELSE IF( lsame( jobvsl,
'V' ) )
THEN
339 IF( lsame( jobvsr,
'N' ) )
THEN
342 ELSE IF( lsame( jobvsr,
'V' ) )
THEN
350 wantst = lsame( sort,
'S' )
355 lquery = ( lwork.EQ.-1 )
356 IF( ijobvl.LE.0 )
THEN
358 ELSE IF( ijobvr.LE.0 )
THEN
360 ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort,
'N' ) ) )
THEN
362 ELSE IF( n.LT.0 )
THEN
364 ELSE IF( lda.LT.max( 1, n ) )
THEN
366 ELSE IF( ldb.LT.max( 1, n ) )
THEN
368 ELSE IF( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) )
THEN
370 ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) )
THEN
372 ELSE IF( lwork.LT.max( 1, 2*n ) .AND. .NOT.lquery )
THEN
379 CALL cgeqrf( n, n, b, ldb, work, work, -1, ierr )
380 lwkopt = max( 1, n + int( work( 1 ) ) )
381 CALL cunmqr(
'L',
'C', n, n, n, b, ldb, work, a, lda, work,
383 lwkopt = max( lwkopt, n + int( work( 1 ) ) )
385 CALL cungqr( n, n, n, vsl, ldvsl, work, work, -1,
387 lwkopt = max( lwkopt, n + int( work( 1 ) ) )
389 CALL cgghd3( jobvsl, jobvsr, n, 1, n, a, lda, b, ldb, vsl,
390 $ ldvsl, vsr, ldvsr, work, -1, ierr )
391 lwkopt = max( lwkopt, n + int( work( 1 ) ) )
392 CALL claqz0(
'S', jobvsl, jobvsr, n, 1, n, a, lda, b, ldb,
393 $ alpha, beta, vsl, ldvsl, vsr, ldvsr, work, -1,
395 lwkopt = max( lwkopt, int( work( 1 ) ) )
397 CALL ctgsen( 0, ilvsl, ilvsr, bwork, n, a, lda, b, ldb,
398 $ alpha, beta, vsl, ldvsl, vsr, ldvsr, sdim,
399 $ pvsl, pvsr, dif, work, -1, idum, 1, ierr )
400 lwkopt = max( lwkopt, int( work( 1 ) ) )
402 work( 1 ) = cmplx( lwkopt )
407 CALL xerbla(
'CGGES3 ', -info )
409 ELSE IF( lquery )
THEN
423 smlnum = slamch(
'S' )
424 bignum = one / smlnum
425 CALL slabad( smlnum, bignum )
426 smlnum = sqrt( smlnum ) / eps
427 bignum = one / smlnum
431 anrm = clange(
'M', n, n, a, lda, rwork )
433 IF( anrm.GT.zero .AND. anrm.LT.smlnum )
THEN
436 ELSE IF( anrm.GT.bignum )
THEN
442 $
CALL clascl(
'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
446 bnrm = clange(
'M', n, n, b, ldb, rwork )
448 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum )
THEN
451 ELSE IF( bnrm.GT.bignum )
THEN
457 $
CALL clascl(
'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
464 CALL cggbal(
'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
465 $ rwork( iright ), rwork( irwrk ), ierr )
469 irows = ihi + 1 - ilo
473 CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
474 $ work( iwrk ), lwork+1-iwrk, ierr )
478 CALL cunmqr(
'L',
'C', irows, icols, irows, b( ilo, ilo ), ldb,
479 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
480 $ lwork+1-iwrk, ierr )
485 CALL claset(
'Full', n, n, czero, cone, vsl, ldvsl )
486 IF( irows.GT.1 )
THEN
487 CALL clacpy(
'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
488 $ vsl( ilo+1, ilo ), ldvsl )
490 CALL cungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
491 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
497 $
CALL claset(
'Full', n, n, czero, cone, vsr, ldvsr )
501 CALL cgghd3( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
502 $ ldvsl, vsr, ldvsr, work( iwrk ), lwork+1-iwrk, ierr )
509 CALL claqz0(
'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
510 $ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwrk ),
511 $ lwork+1-iwrk, rwork( irwrk ), 0, ierr )
513 IF( ierr.GT.0 .AND. ierr.LE.n )
THEN
515 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n )
THEN
530 $
CALL clascl(
'G', 0, 0, anrm, anrmto, n, 1, alpha, n, ierr )
532 $
CALL clascl(
'G', 0, 0, bnrm, bnrmto, n, 1, beta, n, ierr )
537 bwork( i ) = selctg( alpha( i ), beta( i ) )
540 CALL ctgsen( 0, ilvsl, ilvsr, bwork, n, a, lda, b, ldb, alpha,
541 $ beta, vsl, ldvsl, vsr, ldvsr, sdim, pvsl, pvsr,
542 $ dif, work( iwrk ), lwork-iwrk+1, idum, 1, ierr )
551 $
CALL cggbak(
'P',
'L', n, ilo, ihi, rwork( ileft ),
552 $ rwork( iright ), n, vsl, ldvsl, ierr )
554 $
CALL cggbak(
'P',
'R', n, ilo, ihi, rwork( ileft ),
555 $ rwork( iright ), n, vsr, ldvsr, ierr )
560 CALL clascl(
'U', 0, 0, anrmto, anrm, n, n, a, lda, ierr )
561 CALL clascl(
'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
565 CALL clascl(
'U', 0, 0, bnrmto, bnrm, n, n, b, ldb, ierr )
566 CALL clascl(
'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
576 cursl = selctg( alpha( i ), beta( i ) )
579 IF( cursl .AND. .NOT.lastsl )
588 work( 1 ) = cmplx( lwkopt )
subroutine slabad(SMALL, LARGE)
SLABAD
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
recursive subroutine claqz0(WANTS, WANTQ, WANTZ, N, ILO, IHI, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, REC, INFO)
CLAQZ0
subroutine cgges3(JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK, BWORK, INFO)
CGGES3 computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE ...
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 clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine ctgsen(IJOB, WANTQ, WANTZ, SELECT, N, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, M, PL, PR, DIF, WORK, LWORK, IWORK, LIWORK, INFO)
CTGSEN
subroutine cgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
CGGHD3
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR