LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ zgegv()

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.
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 280 of file zgegv.f.

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