LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zgegv ( character  JOBVL,
character  JOBVR,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
complex*16, dimension( ldb, * )  B,
integer  LDB,
complex*16, dimension( * )  ALPHA,
complex*16, dimension( * )  BETA,
complex*16, dimension( ldvl, * )  VL,
integer  LDVL,
complex*16, dimension( ldvr, * )  VR,
integer  LDVR,
complex*16, dimension( * )  WORK,
integer  LWORK,
double precision, dimension( * )  RWORK,
integer  INFO 
)

ZGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

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

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

 ZGEGV computes the eigenvalues and, optionally, the left and/or right
 eigenvectors of a complex 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 COMPLEX*16 array, dimension (LDA, N)
          On entry, the matrix A.
          If JOBVL = 'V' or JOBVR = 'V', then on exit A
          contains the 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 elements
          of the Schur form will be correct.  See ZGGHRD and ZHGEQZ
          for details.
[in]LDA
          LDA is INTEGER
          The leading dimension of A.  LDA >= max(1,N).
[in,out]B
          B is COMPLEX*16 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 the diagonal
          elements of B will be correct.  See ZGGHRD and ZHGEQZ for
          details.
[in]LDB
          LDB is INTEGER
          The leading dimension of B.  LDB >= max(1,N).
[out]ALPHA
          ALPHA is COMPLEX*16 array, dimension (N)
          The complex scalars alpha that define the eigenvalues of
          GNEP.
[out]BETA
          BETA is COMPLEX*16 array, dimension (N)
          The complex scalars beta that define the eigenvalues of GNEP.
          
          Together, the quantities alpha = ALPHA(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 COMPLEX*16 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.
          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 COMPLEX*16 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.
          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 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 COMPLEX*16 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,2*N).
          For good performance, LWORK must generally be larger.
          To compute the optimal value of LWORK, call ILAENV to get
          blocksizes (for ZGEQRF, ZUNMQR, and ZUNGQR.)  Then compute:
          NB  -- MAX of the blocksizes for ZGEQRF, ZUNMQR, and ZUNGQR;
          The optimal LWORK is  MAX( 2*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]RWORK
          RWORK is DOUBLE PRECISION array, dimension (8*N)
[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 ALPHA(j) and BETA(j) should be
                correct for j=INFO+1,...,N.
          > N:  errors that usually indicate LAPACK problems:
                =N+1: error return from ZGGBAL
                =N+2: error return from ZGEQRF
                =N+3: error return from ZUNMQR
                =N+4: error return from ZUNGQR
                =N+5: error return from ZGGHRD
                =N+6: error return from ZHGEQZ (other than failed
                                               iteration)
                =N+7: error return from ZTGEVC
                =N+8: error return from ZGGBAK (computing VL)
                =N+9: error return from ZGGBAK (computing VR)
                =N+10: error return from ZLASCL (various calls)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
  Balancing
  ---------

  This driver calls ZGGBAL 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, ZGGBAK 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 complex Schur
  form[*] of the "balanced" versions of A and B.  If no eigenvectors
  are computed, then only the diagonal blocks will be correct.

  [*] In other words, upper triangular form.

Definition at line 284 of file zgegv.f.

284 *
285 * -- LAPACK driver routine (version 3.4.0) --
286 * -- LAPACK is a software package provided by Univ. of Tennessee, --
287 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
288 * November 2011
289 *
290 * .. Scalar Arguments ..
291  CHARACTER jobvl, jobvr
292  INTEGER info, lda, ldb, ldvl, ldvr, lwork, n
293 * ..
294 * .. Array Arguments ..
295  DOUBLE PRECISION rwork( * )
296  COMPLEX*16 a( lda, * ), alpha( * ), b( ldb, * ),
297  $ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
298  $ work( * )
299 * ..
300 *
301 * =====================================================================
302 *
303 * .. Parameters ..
304  DOUBLE PRECISION zero, one
305  parameter ( zero = 0.0d0, one = 1.0d0 )
306  COMPLEX*16 czero, cone
307  parameter ( czero = ( 0.0d0, 0.0d0 ),
308  $ cone = ( 1.0d0, 0.0d0 ) )
309 * ..
310 * .. Local Scalars ..
311  LOGICAL ilimit, ilv, ilvl, ilvr, lquery
312  CHARACTER chtemp
313  INTEGER icols, ihi, iinfo, ijobvl, ijobvr, ileft, ilo,
314  $ in, iright, irows, irwork, itau, iwork, jc, jr,
315  $ lopt, lwkmin, lwkopt, nb, nb1, nb2, nb3
316  DOUBLE PRECISION absai, absar, absb, anrm, anrm1, anrm2, bnrm,
317  $ bnrm1, bnrm2, eps, safmax, safmin, salfai,
318  $ salfar, sbeta, scale, temp
319  COMPLEX*16 x
320 * ..
321 * .. Local Arrays ..
322  LOGICAL ldumma( 1 )
323 * ..
324 * .. External Subroutines ..
325  EXTERNAL xerbla, zgeqrf, zggbak, zggbal, zgghrd, zhgeqz,
327 * ..
328 * .. External Functions ..
329  LOGICAL lsame
330  INTEGER ilaenv
331  DOUBLE PRECISION dlamch, zlange
332  EXTERNAL lsame, ilaenv, dlamch, zlange
333 * ..
334 * .. Intrinsic Functions ..
335  INTRINSIC abs, dble, dcmplx, dimag, int, max
336 * ..
337 * .. Statement Functions ..
338  DOUBLE PRECISION abs1
339 * ..
340 * .. Statement Function definitions ..
341  abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
342 * ..
343 * .. Executable Statements ..
344 *
345 * Decode the input arguments
346 *
347  IF( lsame( jobvl, 'N' ) ) THEN
348  ijobvl = 1
349  ilvl = .false.
350  ELSE IF( lsame( jobvl, 'V' ) ) THEN
351  ijobvl = 2
352  ilvl = .true.
353  ELSE
354  ijobvl = -1
355  ilvl = .false.
356  END IF
357 *
358  IF( lsame( jobvr, 'N' ) ) THEN
359  ijobvr = 1
360  ilvr = .false.
361  ELSE IF( lsame( jobvr, 'V' ) ) THEN
362  ijobvr = 2
363  ilvr = .true.
364  ELSE
365  ijobvr = -1
366  ilvr = .false.
367  END IF
368  ilv = ilvl .OR. ilvr
369 *
370 * Test the input arguments
371 *
372  lwkmin = max( 2*n, 1 )
373  lwkopt = lwkmin
374  work( 1 ) = lwkopt
375  lquery = ( lwork.EQ.-1 )
376  info = 0
377  IF( ijobvl.LE.0 ) THEN
378  info = -1
379  ELSE IF( ijobvr.LE.0 ) THEN
380  info = -2
381  ELSE IF( n.LT.0 ) THEN
382  info = -3
383  ELSE IF( lda.LT.max( 1, n ) ) THEN
384  info = -5
385  ELSE IF( ldb.LT.max( 1, n ) ) THEN
386  info = -7
387  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
388  info = -11
389  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
390  info = -13
391  ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
392  info = -15
393  END IF
394 *
395  IF( info.EQ.0 ) THEN
396  nb1 = ilaenv( 1, 'ZGEQRF', ' ', n, n, -1, -1 )
397  nb2 = ilaenv( 1, 'ZUNMQR', ' ', n, n, n, -1 )
398  nb3 = ilaenv( 1, 'ZUNGQR', ' ', n, n, n, -1 )
399  nb = max( nb1, nb2, nb3 )
400  lopt = max( 2*n, n*( nb+1 ) )
401  work( 1 ) = lopt
402  END IF
403 *
404  IF( info.NE.0 ) THEN
405  CALL xerbla( 'ZGEGV ', -info )
406  RETURN
407  ELSE IF( lquery ) THEN
408  RETURN
409  END IF
410 *
411 * Quick return if possible
412 *
413  IF( n.EQ.0 )
414  $ RETURN
415 *
416 * Get machine constants
417 *
418  eps = dlamch( 'E' )*dlamch( 'B' )
419  safmin = dlamch( 'S' )
420  safmin = safmin + safmin
421  safmax = one / safmin
422 *
423 * Scale A
424 *
425  anrm = zlange( 'M', n, n, a, lda, rwork )
426  anrm1 = anrm
427  anrm2 = one
428  IF( anrm.LT.one ) THEN
429  IF( safmax*anrm.LT.one ) THEN
430  anrm1 = safmin
431  anrm2 = safmax*anrm
432  END IF
433  END IF
434 *
435  IF( anrm.GT.zero ) THEN
436  CALL zlascl( 'G', -1, -1, anrm, one, n, n, a, lda, iinfo )
437  IF( iinfo.NE.0 ) THEN
438  info = n + 10
439  RETURN
440  END IF
441  END IF
442 *
443 * Scale B
444 *
445  bnrm = zlange( 'M', n, n, b, ldb, rwork )
446  bnrm1 = bnrm
447  bnrm2 = one
448  IF( bnrm.LT.one ) THEN
449  IF( safmax*bnrm.LT.one ) THEN
450  bnrm1 = safmin
451  bnrm2 = safmax*bnrm
452  END IF
453  END IF
454 *
455  IF( bnrm.GT.zero ) THEN
456  CALL zlascl( 'G', -1, -1, bnrm, one, n, n, b, ldb, iinfo )
457  IF( iinfo.NE.0 ) THEN
458  info = n + 10
459  RETURN
460  END IF
461  END IF
462 *
463 * Permute the matrix to make it more nearly triangular
464 * Also "balance" the matrix.
465 *
466  ileft = 1
467  iright = n + 1
468  irwork = iright + n
469  CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
470  $ rwork( iright ), rwork( irwork ), iinfo )
471  IF( iinfo.NE.0 ) THEN
472  info = n + 1
473  GO TO 80
474  END IF
475 *
476 * Reduce B to triangular form, and initialize VL and/or VR
477 *
478  irows = ihi + 1 - ilo
479  IF( ilv ) THEN
480  icols = n + 1 - ilo
481  ELSE
482  icols = irows
483  END IF
484  itau = 1
485  iwork = itau + irows
486  CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
487  $ work( iwork ), lwork+1-iwork, iinfo )
488  IF( iinfo.GE.0 )
489  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
490  IF( iinfo.NE.0 ) THEN
491  info = n + 2
492  GO TO 80
493  END IF
494 *
495  CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
496  $ work( itau ), a( ilo, ilo ), lda, work( iwork ),
497  $ lwork+1-iwork, iinfo )
498  IF( iinfo.GE.0 )
499  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
500  IF( iinfo.NE.0 ) THEN
501  info = n + 3
502  GO TO 80
503  END IF
504 *
505  IF( ilvl ) THEN
506  CALL zlaset( 'Full', n, n, czero, cone, vl, ldvl )
507  CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
508  $ vl( ilo+1, ilo ), ldvl )
509  CALL zungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
510  $ work( itau ), work( iwork ), lwork+1-iwork,
511  $ iinfo )
512  IF( iinfo.GE.0 )
513  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
514  IF( iinfo.NE.0 ) THEN
515  info = n + 4
516  GO TO 80
517  END IF
518  END IF
519 *
520  IF( ilvr )
521  $ CALL zlaset( 'Full', n, n, czero, cone, vr, ldvr )
522 *
523 * Reduce to generalized Hessenberg form
524 *
525  IF( ilv ) THEN
526 *
527 * Eigenvectors requested -- work on whole matrix.
528 *
529  CALL zgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
530  $ ldvl, vr, ldvr, iinfo )
531  ELSE
532  CALL zgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
533  $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, iinfo )
534  END IF
535  IF( iinfo.NE.0 ) THEN
536  info = n + 5
537  GO TO 80
538  END IF
539 *
540 * Perform QZ algorithm
541 *
542  iwork = itau
543  IF( ilv ) THEN
544  chtemp = 'S'
545  ELSE
546  chtemp = 'E'
547  END IF
548  CALL zhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
549  $ alpha, beta, vl, ldvl, vr, ldvr, work( iwork ),
550  $ lwork+1-iwork, rwork( irwork ), iinfo )
551  IF( iinfo.GE.0 )
552  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
553  IF( iinfo.NE.0 ) THEN
554  IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
555  info = iinfo
556  ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
557  info = iinfo - n
558  ELSE
559  info = n + 6
560  END IF
561  GO TO 80
562  END IF
563 *
564  IF( ilv ) THEN
565 *
566 * Compute Eigenvectors
567 *
568  IF( ilvl ) THEN
569  IF( ilvr ) THEN
570  chtemp = 'B'
571  ELSE
572  chtemp = 'L'
573  END IF
574  ELSE
575  chtemp = 'R'
576  END IF
577 *
578  CALL ztgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
579  $ vr, ldvr, n, in, work( iwork ), rwork( irwork ),
580  $ iinfo )
581  IF( iinfo.NE.0 ) THEN
582  info = n + 7
583  GO TO 80
584  END IF
585 *
586 * Undo balancing on VL and VR, rescale
587 *
588  IF( ilvl ) THEN
589  CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
590  $ rwork( iright ), n, vl, ldvl, iinfo )
591  IF( iinfo.NE.0 ) THEN
592  info = n + 8
593  GO TO 80
594  END IF
595  DO 30 jc = 1, n
596  temp = zero
597  DO 10 jr = 1, n
598  temp = max( temp, abs1( vl( jr, jc ) ) )
599  10 CONTINUE
600  IF( temp.LT.safmin )
601  $ GO TO 30
602  temp = one / temp
603  DO 20 jr = 1, n
604  vl( jr, jc ) = vl( jr, jc )*temp
605  20 CONTINUE
606  30 CONTINUE
607  END IF
608  IF( ilvr ) THEN
609  CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
610  $ rwork( iright ), n, vr, ldvr, iinfo )
611  IF( iinfo.NE.0 ) THEN
612  info = n + 9
613  GO TO 80
614  END IF
615  DO 60 jc = 1, n
616  temp = zero
617  DO 40 jr = 1, n
618  temp = max( temp, abs1( vr( jr, jc ) ) )
619  40 CONTINUE
620  IF( temp.LT.safmin )
621  $ GO TO 60
622  temp = one / temp
623  DO 50 jr = 1, n
624  vr( jr, jc ) = vr( jr, jc )*temp
625  50 CONTINUE
626  60 CONTINUE
627  END IF
628 *
629 * End of eigenvector calculation
630 *
631  END IF
632 *
633 * Undo scaling in alpha, beta
634 *
635 * Note: this does not give the alpha and beta for the unscaled
636 * problem.
637 *
638 * Un-scaling is limited to avoid underflow in alpha and beta
639 * if they are significant.
640 *
641  DO 70 jc = 1, n
642  absar = abs( dble( alpha( jc ) ) )
643  absai = abs( dimag( alpha( jc ) ) )
644  absb = abs( dble( beta( jc ) ) )
645  salfar = anrm*dble( alpha( jc ) )
646  salfai = anrm*dimag( alpha( jc ) )
647  sbeta = bnrm*dble( beta( jc ) )
648  ilimit = .false.
649  scale = one
650 *
651 * Check for significant underflow in imaginary part of ALPHA
652 *
653  IF( abs( salfai ).LT.safmin .AND. absai.GE.
654  $ max( safmin, eps*absar, eps*absb ) ) THEN
655  ilimit = .true.
656  scale = ( safmin / anrm1 ) / max( safmin, anrm2*absai )
657  END IF
658 *
659 * Check for significant underflow in real part of ALPHA
660 *
661  IF( abs( salfar ).LT.safmin .AND. absar.GE.
662  $ max( safmin, eps*absai, eps*absb ) ) THEN
663  ilimit = .true.
664  scale = max( scale, ( safmin / anrm1 ) /
665  $ max( safmin, anrm2*absar ) )
666  END IF
667 *
668 * Check for significant underflow in BETA
669 *
670  IF( abs( sbeta ).LT.safmin .AND. absb.GE.
671  $ max( safmin, eps*absar, eps*absai ) ) THEN
672  ilimit = .true.
673  scale = max( scale, ( safmin / bnrm1 ) /
674  $ max( safmin, bnrm2*absb ) )
675  END IF
676 *
677 * Check for possible overflow when limiting scaling
678 *
679  IF( ilimit ) THEN
680  temp = ( scale*safmin )*max( abs( salfar ), abs( salfai ),
681  $ abs( sbeta ) )
682  IF( temp.GT.one )
683  $ scale = scale / temp
684  IF( scale.LT.one )
685  $ ilimit = .false.
686  END IF
687 *
688 * Recompute un-scaled ALPHA, BETA if necessary.
689 *
690  IF( ilimit ) THEN
691  salfar = ( scale*dble( alpha( jc ) ) )*anrm
692  salfai = ( scale*dimag( alpha( jc ) ) )*anrm
693  sbeta = ( scale*beta( jc ) )*bnrm
694  END IF
695  alpha( jc ) = dcmplx( salfar, salfai )
696  beta( jc ) = sbeta
697  70 CONTINUE
698 *
699  80 CONTINUE
700  work( 1 ) = lwkopt
701 *
702  RETURN
703 *
704 * End of ZGEGV
705 *
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:105
subroutine zgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
ZGGHRD
Definition: zgghrd.f:206
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
subroutine zgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
ZGEQRF VARIANT: left-looking Level 3 BLAS of the algorithm.
Definition: zgeqrf.f:151
subroutine zggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
ZGGBAK
Definition: zggbak.f:150
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...
Definition: zlaset.f:108
subroutine zggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
ZGGBAL
Definition: zggbal.f:179
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMQR
Definition: zunmqr.f:169
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: zlange.f:117
subroutine zungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
ZUNGQR
Definition: zungqr.f:130
subroutine ztgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTGEVC
Definition: ztgevc.f:221
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine zhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
ZHGEQZ
Definition: zhgeqz.f:286
subroutine zlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: zlascl.f:145
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function: