LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ sgegv()

subroutine sgegv ( character  jobvl,
character  jobvr,
integer  n,
real, dimension( lda, * )  a,
integer  lda,
real, dimension( ldb, * )  b,
integer  ldb,
real, dimension( * )  alphar,
real, dimension( * )  alphai,
real, dimension( * )  beta,
real, dimension( ldvl, * )  vl,
integer  ldvl,
real, dimension( ldvr, * )  vr,
integer  ldvr,
real, dimension( * )  work,
integer  lwork,
integer  info 
)

SGEGV computes the eigenvalues and, optionally, the left and/or right eigenvectors of a real matrix pair (A,B).

Download SGEGV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 This routine is deprecated and has been replaced by routine SGGEV.

 SGEGV computes the eigenvalues and, optionally, the left and/or right
 eigenvectors of a real matrix pair (A,B).
 Given two square matrices A and B,
 the generalized nonsymmetric eigenvalue problem (GNEP) is to find the
 eigenvalues lambda and corresponding (non-zero) eigenvectors x such
 that

    A*x = lambda*B*x.

 An alternate form is to find the eigenvalues mu and corresponding
 eigenvectors y such that

    mu*A*y = B*y.

 These two forms are equivalent with mu = 1/lambda and x = y if
 neither lambda nor mu is zero.  In order to deal with the case that
 lambda or mu is zero or small, two values alpha and beta are returned
 for each eigenvalue, such that lambda = alpha/beta and
 mu = beta/alpha.

 The vectors x and y in the above equations are right eigenvectors of
 the matrix pair (A,B).  Vectors u and v satisfying

    u**H*A = lambda*u**H*B  or  mu*v**H*A = v**H*B

 are left eigenvectors of (A,B).

 Note: this routine performs "full balancing" on A and B
Parameters
[in]JOBVL
          JOBVL is CHARACTER*1
          = 'N':  do not compute the left generalized eigenvectors;
          = 'V':  compute the left generalized eigenvectors (returned
                  in VL).
[in]JOBVR
          JOBVR is CHARACTER*1
          = 'N':  do not compute the right generalized eigenvectors;
          = 'V':  compute the right generalized eigenvectors (returned
                  in VR).
[in]N
          N is INTEGER
          The order of the matrices A, B, VL, and VR.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the matrix A.
          If JOBVL = 'V' or JOBVR = 'V', then on exit A
          contains the real Schur form of A from the generalized Schur
          factorization of the pair (A,B) after balancing.
          If no eigenvectors were computed, then only the diagonal
          blocks from the Schur form will be correct.  See SGGHRD and
          SHGEQZ for details.
[in]LDA
          LDA is INTEGER
          The leading dimension of A.  LDA >= max(1,N).
[in,out]B
          B is REAL array, dimension (LDB, N)
          On entry, the matrix B.
          If JOBVL = 'V' or JOBVR = 'V', then on exit B contains the
          upper triangular matrix obtained from B in the generalized
          Schur factorization of the pair (A,B) after balancing.
          If no eigenvectors were computed, then only those elements of
          B corresponding to the diagonal blocks from the Schur form of
          A will be correct.  See SGGHRD and SHGEQZ for details.
[in]LDB
          LDB is INTEGER
          The leading dimension of B.  LDB >= max(1,N).
[out]ALPHAR
          ALPHAR is REAL array, dimension (N)
          The real parts of each scalar alpha defining an eigenvalue of
          GNEP.
[out]ALPHAI
          ALPHAI is REAL array, dimension (N)
          The imaginary parts of each scalar alpha defining an
          eigenvalue of GNEP.  If ALPHAI(j) is zero, then the j-th
          eigenvalue is real; if positive, then the j-th and
          (j+1)-st eigenvalues are a complex conjugate pair, with
          ALPHAI(j+1) = -ALPHAI(j).
[out]BETA
          BETA is REAL array, dimension (N)
          The scalars beta that define the eigenvalues of GNEP.

          Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and
          beta = BETA(j) represent the j-th eigenvalue of the matrix
          pair (A,B), in one of the forms lambda = alpha/beta or
          mu = beta/alpha.  Since either lambda or mu may overflow,
          they should not, in general, be computed.
[out]VL
          VL is REAL array, dimension (LDVL,N)
          If JOBVL = 'V', the left eigenvectors u(j) are stored
          in the columns of VL, in the same order as their eigenvalues.
          If the j-th eigenvalue is real, then u(j) = VL(:,j).
          If the j-th and (j+1)-st eigenvalues form a complex conjugate
          pair, then
             u(j) = VL(:,j) + i*VL(:,j+1)
          and
            u(j+1) = VL(:,j) - i*VL(:,j+1).

          Each eigenvector is scaled so that its largest component has
          abs(real part) + abs(imag. part) = 1, except for eigenvectors
          corresponding to an eigenvalue with alpha = beta = 0, which
          are set to zero.
          Not referenced if JOBVL = 'N'.
[in]LDVL
          LDVL is INTEGER
          The leading dimension of the matrix VL. LDVL >= 1, and
          if JOBVL = 'V', LDVL >= N.
[out]VR
          VR is REAL array, dimension (LDVR,N)
          If JOBVR = 'V', the right eigenvectors x(j) are stored
          in the columns of VR, in the same order as their eigenvalues.
          If the j-th eigenvalue is real, then x(j) = VR(:,j).
          If the j-th and (j+1)-st eigenvalues form a complex conjugate
          pair, then
            x(j) = VR(:,j) + i*VR(:,j+1)
          and
            x(j+1) = VR(:,j) - i*VR(:,j+1).

          Each eigenvector is scaled so that its largest component has
          abs(real part) + abs(imag. part) = 1, except for eigenvalues
          corresponding to an eigenvalue with alpha = beta = 0, which
          are set to zero.
          Not referenced if JOBVR = 'N'.
[in]LDVR
          LDVR is INTEGER
          The leading dimension of the matrix VR. LDVR >= 1, and
          if JOBVR = 'V', LDVR >= N.
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= max(1,8*N).
          For good performance, LWORK must generally be larger.
          To compute the optimal value of LWORK, call ILAENV to get
          blocksizes (for SGEQRF, SORMQR, and SORGQR.)  Then compute:
          NB  -- MAX of the blocksizes for SGEQRF, SORMQR, and SORGQR;
          The optimal LWORK is:
              2*N + MAX( 6*N, N*(NB+1) ).

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          = 1,...,N:
                The QZ iteration failed.  No eigenvectors have been
                calculated, but ALPHAR(j), ALPHAI(j), and BETA(j)
                should be correct for j=INFO+1,...,N.
          > N:  errors that usually indicate LAPACK problems:
                =N+1: error return from SGGBAL
                =N+2: error return from SGEQRF
                =N+3: error return from SORMQR
                =N+4: error return from SORGQR
                =N+5: error return from SGGHRD
                =N+6: error return from SHGEQZ (other than failed
                                                iteration)
                =N+7: error return from STGEVC
                =N+8: error return from SGGBAK (computing VL)
                =N+9: error return from SGGBAK (computing VR)
                =N+10: error return from SLASCL (various calls)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  Balancing
  ---------

  This driver calls SGGBAL to both permute and scale rows and columns
  of A and B.  The permutations PL and PR are chosen so that PL*A*PR
  and PL*B*R will be upper triangular except for the diagonal blocks
  A(i:j,i:j) and B(i:j,i:j), with i and j as close together as
  possible.  The diagonal scaling matrices DL and DR are chosen so
  that the pair  DL*PL*A*PR*DR, DL*PL*B*PR*DR have elements close to
  one (except for the elements that start out zero.)

  After the eigenvalues and eigenvectors of the balanced matrices
  have been computed, SGGBAK transforms the eigenvectors back to what
  they would have been (in perfect arithmetic) if they had not been
  balanced.

  Contents of A and B on Exit
  -------- -- - --- - -- ----

  If any eigenvectors are computed (either JOBVL='V' or JOBVR='V' or
  both), then on exit the arrays A and B will contain the real Schur
  form[*] of the "balanced" versions of A and B.  If no eigenvectors
  are computed, then only the diagonal blocks will be correct.

  [*] See SHGEQZ, SGEGS, or read the book "Matrix Computations",
      by Golub & van Loan, pub. by Johns Hopkins U. Press.

Definition at line 304 of file sgegv.f.

306*
307* -- LAPACK driver routine --
308* -- LAPACK is a software package provided by Univ. of Tennessee, --
309* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
310*
311* .. Scalar Arguments ..
312 CHARACTER JOBVL, JOBVR
313 INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
314* ..
315* .. Array Arguments ..
316 REAL A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
317 $ B( LDB, * ), BETA( * ), VL( LDVL, * ),
318 $ VR( LDVR, * ), WORK( * )
319* ..
320*
321* =====================================================================
322*
323* .. Parameters ..
324 REAL ZERO, ONE
325 parameter( zero = 0.0e0, one = 1.0e0 )
326* ..
327* .. Local Scalars ..
328 LOGICAL ILIMIT, ILV, ILVL, ILVR, LQUERY
329 CHARACTER CHTEMP
330 INTEGER ICOLS, IHI, IINFO, IJOBVL, IJOBVR, ILEFT, ILO,
331 $ IN, IRIGHT, IROWS, ITAU, IWORK, JC, JR, LOPT,
332 $ LWKMIN, LWKOPT, NB, NB1, NB2, NB3
333 REAL ABSAI, ABSAR, ABSB, ANRM, ANRM1, ANRM2, BNRM,
334 $ BNRM1, BNRM2, EPS, ONEPLS, SAFMAX, SAFMIN,
335 $ SALFAI, SALFAR, SBETA, SCALE, TEMP
336* ..
337* .. Local Arrays ..
338 LOGICAL LDUMMA( 1 )
339* ..
340* .. External Subroutines ..
341 EXTERNAL sgeqrf, sggbak, sggbal, sgghrd, shgeqz, slacpy,
343* ..
344* .. External Functions ..
345 LOGICAL LSAME
346 INTEGER ILAENV
347 REAL SLAMCH, SLANGE
348 EXTERNAL ilaenv, lsame, slamch, slange
349* ..
350* .. Intrinsic Functions ..
351 INTRINSIC abs, int, max
352* ..
353* .. Executable Statements ..
354*
355* Decode the input arguments
356*
357 IF( lsame( jobvl, 'N' ) ) THEN
358 ijobvl = 1
359 ilvl = .false.
360 ELSE IF( lsame( jobvl, 'V' ) ) THEN
361 ijobvl = 2
362 ilvl = .true.
363 ELSE
364 ijobvl = -1
365 ilvl = .false.
366 END IF
367*
368 IF( lsame( jobvr, 'N' ) ) THEN
369 ijobvr = 1
370 ilvr = .false.
371 ELSE IF( lsame( jobvr, 'V' ) ) THEN
372 ijobvr = 2
373 ilvr = .true.
374 ELSE
375 ijobvr = -1
376 ilvr = .false.
377 END IF
378 ilv = ilvl .OR. ilvr
379*
380* Test the input arguments
381*
382 lwkmin = max( 8*n, 1 )
383 lwkopt = lwkmin
384 work( 1 ) = lwkopt
385 lquery = ( lwork.EQ.-1 )
386 info = 0
387 IF( ijobvl.LE.0 ) THEN
388 info = -1
389 ELSE IF( ijobvr.LE.0 ) THEN
390 info = -2
391 ELSE IF( n.LT.0 ) THEN
392 info = -3
393 ELSE IF( lda.LT.max( 1, n ) ) THEN
394 info = -5
395 ELSE IF( ldb.LT.max( 1, n ) ) THEN
396 info = -7
397 ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
398 info = -12
399 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
400 info = -14
401 ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
402 info = -16
403 END IF
404*
405 IF( info.EQ.0 ) THEN
406 nb1 = ilaenv( 1, 'SGEQRF', ' ', n, n, -1, -1 )
407 nb2 = ilaenv( 1, 'SORMQR', ' ', n, n, n, -1 )
408 nb3 = ilaenv( 1, 'SORGQR', ' ', n, n, n, -1 )
409 nb = max( nb1, nb2, nb3 )
410 lopt = 2*n + max( 6*n, n*(nb+1) )
411 work( 1 ) = lopt
412 END IF
413*
414 IF( info.NE.0 ) THEN
415 CALL xerbla( 'SGEGV ', -info )
416 RETURN
417 ELSE IF( lquery ) THEN
418 RETURN
419 END IF
420*
421* Quick return if possible
422*
423 IF( n.EQ.0 )
424 $ RETURN
425*
426* Get machine constants
427*
428 eps = slamch( 'E' )*slamch( 'B' )
429 safmin = slamch( 'S' )
430 safmin = safmin + safmin
431 safmax = one / safmin
432 onepls = one + ( 4*eps )
433*
434* Scale A
435*
436 anrm = slange( 'M', n, n, a, lda, work )
437 anrm1 = anrm
438 anrm2 = one
439 IF( anrm.LT.one ) THEN
440 IF( safmax*anrm.LT.one ) THEN
441 anrm1 = safmin
442 anrm2 = safmax*anrm
443 END IF
444 END IF
445*
446 IF( anrm.GT.zero ) THEN
447 CALL slascl( 'G', -1, -1, anrm, one, n, n, a, lda, iinfo )
448 IF( iinfo.NE.0 ) THEN
449 info = n + 10
450 RETURN
451 END IF
452 END IF
453*
454* Scale B
455*
456 bnrm = slange( 'M', n, n, b, ldb, work )
457 bnrm1 = bnrm
458 bnrm2 = one
459 IF( bnrm.LT.one ) THEN
460 IF( safmax*bnrm.LT.one ) THEN
461 bnrm1 = safmin
462 bnrm2 = safmax*bnrm
463 END IF
464 END IF
465*
466 IF( bnrm.GT.zero ) THEN
467 CALL slascl( 'G', -1, -1, bnrm, one, n, n, b, ldb, iinfo )
468 IF( iinfo.NE.0 ) THEN
469 info = n + 10
470 RETURN
471 END IF
472 END IF
473*
474* Permute the matrix to make it more nearly triangular
475* Workspace layout: (8*N words -- "work" requires 6*N words)
476* left_permutation, right_permutation, work...
477*
478 ileft = 1
479 iright = n + 1
480 iwork = iright + n
481 CALL sggbal( 'P', n, a, lda, b, ldb, ilo, ihi, work( ileft ),
482 $ work( iright ), work( iwork ), iinfo )
483 IF( iinfo.NE.0 ) THEN
484 info = n + 1
485 GO TO 120
486 END IF
487*
488* Reduce B to triangular form, and initialize VL and/or VR
489* Workspace layout: ("work..." must have at least N words)
490* left_permutation, right_permutation, tau, work...
491*
492 irows = ihi + 1 - ilo
493 IF( ilv ) THEN
494 icols = n + 1 - ilo
495 ELSE
496 icols = irows
497 END IF
498 itau = iwork
499 iwork = itau + irows
500 CALL sgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
501 $ work( iwork ), lwork+1-iwork, iinfo )
502 IF( iinfo.GE.0 )
503 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
504 IF( iinfo.NE.0 ) THEN
505 info = n + 2
506 GO TO 120
507 END IF
508*
509 CALL sormqr( 'L', 'T', irows, icols, irows, b( ilo, ilo ), ldb,
510 $ work( itau ), a( ilo, ilo ), lda, work( iwork ),
511 $ lwork+1-iwork, iinfo )
512 IF( iinfo.GE.0 )
513 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
514 IF( iinfo.NE.0 ) THEN
515 info = n + 3
516 GO TO 120
517 END IF
518*
519 IF( ilvl ) THEN
520 CALL slaset( 'Full', n, n, zero, one, vl, ldvl )
521 CALL slacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
522 $ vl( ilo+1, ilo ), ldvl )
523 CALL sorgqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
524 $ work( itau ), work( iwork ), lwork+1-iwork,
525 $ iinfo )
526 IF( iinfo.GE.0 )
527 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
528 IF( iinfo.NE.0 ) THEN
529 info = n + 4
530 GO TO 120
531 END IF
532 END IF
533*
534 IF( ilvr )
535 $ CALL slaset( 'Full', n, n, zero, one, vr, ldvr )
536*
537* Reduce to generalized Hessenberg form
538*
539 IF( ilv ) THEN
540*
541* Eigenvectors requested -- work on whole matrix.
542*
543 CALL sgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
544 $ ldvl, vr, ldvr, iinfo )
545 ELSE
546 CALL sgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
547 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, iinfo )
548 END IF
549 IF( iinfo.NE.0 ) THEN
550 info = n + 5
551 GO TO 120
552 END IF
553*
554* Perform QZ algorithm
555* Workspace layout: ("work..." must have at least 1 word)
556* left_permutation, right_permutation, work...
557*
558 iwork = itau
559 IF( ilv ) THEN
560 chtemp = 'S'
561 ELSE
562 chtemp = 'E'
563 END IF
564 CALL shgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
565 $ alphar, alphai, beta, vl, ldvl, vr, ldvr,
566 $ work( iwork ), lwork+1-iwork, iinfo )
567 IF( iinfo.GE.0 )
568 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
569 IF( iinfo.NE.0 ) THEN
570 IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
571 info = iinfo
572 ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
573 info = iinfo - n
574 ELSE
575 info = n + 6
576 END IF
577 GO TO 120
578 END IF
579*
580 IF( ilv ) THEN
581*
582* Compute Eigenvectors (STGEVC requires 6*N words of workspace)
583*
584 IF( ilvl ) THEN
585 IF( ilvr ) THEN
586 chtemp = 'B'
587 ELSE
588 chtemp = 'L'
589 END IF
590 ELSE
591 chtemp = 'R'
592 END IF
593*
594 CALL stgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
595 $ vr, ldvr, n, in, work( iwork ), iinfo )
596 IF( iinfo.NE.0 ) THEN
597 info = n + 7
598 GO TO 120
599 END IF
600*
601* Undo balancing on VL and VR, rescale
602*
603 IF( ilvl ) THEN
604 CALL sggbak( 'P', 'L', n, ilo, ihi, work( ileft ),
605 $ work( iright ), n, vl, ldvl, iinfo )
606 IF( iinfo.NE.0 ) THEN
607 info = n + 8
608 GO TO 120
609 END IF
610 DO 50 jc = 1, n
611 IF( alphai( jc ).LT.zero )
612 $ GO TO 50
613 temp = zero
614 IF( alphai( jc ).EQ.zero ) THEN
615 DO 10 jr = 1, n
616 temp = max( temp, abs( vl( jr, jc ) ) )
617 10 CONTINUE
618 ELSE
619 DO 20 jr = 1, n
620 temp = max( temp, abs( vl( jr, jc ) )+
621 $ abs( vl( jr, jc+1 ) ) )
622 20 CONTINUE
623 END IF
624 IF( temp.LT.safmin )
625 $ GO TO 50
626 temp = one / temp
627 IF( alphai( jc ).EQ.zero ) THEN
628 DO 30 jr = 1, n
629 vl( jr, jc ) = vl( jr, jc )*temp
630 30 CONTINUE
631 ELSE
632 DO 40 jr = 1, n
633 vl( jr, jc ) = vl( jr, jc )*temp
634 vl( jr, jc+1 ) = vl( jr, jc+1 )*temp
635 40 CONTINUE
636 END IF
637 50 CONTINUE
638 END IF
639 IF( ilvr ) THEN
640 CALL sggbak( 'P', 'R', n, ilo, ihi, work( ileft ),
641 $ work( iright ), n, vr, ldvr, iinfo )
642 IF( iinfo.NE.0 ) THEN
643 info = n + 9
644 GO TO 120
645 END IF
646 DO 100 jc = 1, n
647 IF( alphai( jc ).LT.zero )
648 $ GO TO 100
649 temp = zero
650 IF( alphai( jc ).EQ.zero ) THEN
651 DO 60 jr = 1, n
652 temp = max( temp, abs( vr( jr, jc ) ) )
653 60 CONTINUE
654 ELSE
655 DO 70 jr = 1, n
656 temp = max( temp, abs( vr( jr, jc ) )+
657 $ abs( vr( jr, jc+1 ) ) )
658 70 CONTINUE
659 END IF
660 IF( temp.LT.safmin )
661 $ GO TO 100
662 temp = one / temp
663 IF( alphai( jc ).EQ.zero ) THEN
664 DO 80 jr = 1, n
665 vr( jr, jc ) = vr( jr, jc )*temp
666 80 CONTINUE
667 ELSE
668 DO 90 jr = 1, n
669 vr( jr, jc ) = vr( jr, jc )*temp
670 vr( jr, jc+1 ) = vr( jr, jc+1 )*temp
671 90 CONTINUE
672 END IF
673 100 CONTINUE
674 END IF
675*
676* End of eigenvector calculation
677*
678 END IF
679*
680* Undo scaling in alpha, beta
681*
682* Note: this does not give the alpha and beta for the unscaled
683* problem.
684*
685* Un-scaling is limited to avoid underflow in alpha and beta
686* if they are significant.
687*
688 DO 110 jc = 1, n
689 absar = abs( alphar( jc ) )
690 absai = abs( alphai( jc ) )
691 absb = abs( beta( jc ) )
692 salfar = anrm*alphar( jc )
693 salfai = anrm*alphai( jc )
694 sbeta = bnrm*beta( jc )
695 ilimit = .false.
696 scale = one
697*
698* Check for significant underflow in ALPHAI
699*
700 IF( abs( salfai ).LT.safmin .AND. absai.GE.
701 $ max( safmin, eps*absar, eps*absb ) ) THEN
702 ilimit = .true.
703 scale = ( onepls*safmin / anrm1 ) /
704 $ max( onepls*safmin, anrm2*absai )
705*
706 ELSE IF( salfai.EQ.zero ) THEN
707*
708* If insignificant underflow in ALPHAI, then make the
709* conjugate eigenvalue real.
710*
711 IF( alphai( jc ).LT.zero .AND. jc.GT.1 ) THEN
712 alphai( jc-1 ) = zero
713 ELSE IF( alphai( jc ).GT.zero .AND. jc.LT.n ) THEN
714 alphai( jc+1 ) = zero
715 END IF
716 END IF
717*
718* Check for significant underflow in ALPHAR
719*
720 IF( abs( salfar ).LT.safmin .AND. absar.GE.
721 $ max( safmin, eps*absai, eps*absb ) ) THEN
722 ilimit = .true.
723 scale = max( scale, ( onepls*safmin / anrm1 ) /
724 $ max( onepls*safmin, anrm2*absar ) )
725 END IF
726*
727* Check for significant underflow in BETA
728*
729 IF( abs( sbeta ).LT.safmin .AND. absb.GE.
730 $ max( safmin, eps*absar, eps*absai ) ) THEN
731 ilimit = .true.
732 scale = max( scale, ( onepls*safmin / bnrm1 ) /
733 $ max( onepls*safmin, bnrm2*absb ) )
734 END IF
735*
736* Check for possible overflow when limiting scaling
737*
738 IF( ilimit ) THEN
739 temp = ( scale*safmin )*max( abs( salfar ), abs( salfai ),
740 $ abs( sbeta ) )
741 IF( temp.GT.one )
742 $ scale = scale / temp
743 IF( scale.LT.one )
744 $ ilimit = .false.
745 END IF
746*
747* Recompute un-scaled ALPHAR, ALPHAI, BETA if necessary.
748*
749 IF( ilimit ) THEN
750 salfar = ( scale*alphar( jc ) )*anrm
751 salfai = ( scale*alphai( jc ) )*anrm
752 sbeta = ( scale*beta( jc ) )*bnrm
753 END IF
754 alphar( jc ) = salfar
755 alphai( jc ) = salfai
756 beta( jc ) = sbeta
757 110 CONTINUE
758*
759 120 CONTINUE
760 work( 1 ) = lwkopt
761*
762 RETURN
763*
764* End of SGEGV
765*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine sgeqrf(m, n, a, lda, tau, work, lwork, info)
SGEQRF
Definition sgeqrf.f:146
subroutine sggbak(job, side, n, ilo, ihi, lscale, rscale, m, v, ldv, info)
SGGBAK
Definition sggbak.f:147
subroutine sggbal(job, n, a, lda, b, ldb, ilo, ihi, lscale, rscale, work, info)
SGGBAL
Definition sggbal.f:177
subroutine sgghrd(compq, compz, n, ilo, ihi, a, lda, b, ldb, q, ldq, z, ldz, info)
SGGHRD
Definition sgghrd.f:207
subroutine shgeqz(job, compq, compz, n, ilo, ihi, h, ldh, t, ldt, alphar, alphai, beta, q, ldq, z, ldz, work, lwork, info)
SHGEQZ
Definition shgeqz.f:304
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
Definition slacpy.f:103
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
real function slange(norm, m, n, a, lda, work)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition slange.f:114
subroutine slascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition slascl.f:143
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.
Definition slaset.f:110
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine stgevc(side, howmny, select, n, s, lds, p, ldp, vl, ldvl, vr, ldvr, mm, m, work, info)
STGEVC
Definition stgevc.f:295
subroutine sorgqr(m, n, k, a, lda, tau, work, lwork, info)
SORGQR
Definition sorgqr.f:128
subroutine sormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
SORMQR
Definition sormqr.f:168
Here is the call graph for this function: