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

◆ zdrgev()

subroutine zdrgev ( integer  NSIZES,
integer, dimension( * )  NN,
integer  NTYPES,
logical, dimension( * )  DOTYPE,
integer, dimension( 4 )  ISEED,
double precision  THRESH,
integer  NOUNIT,
complex*16, dimension( lda, * )  A,
integer  LDA,
complex*16, dimension( lda, * )  B,
complex*16, dimension( lda, * )  S,
complex*16, dimension( lda, * )  T,
complex*16, dimension( ldq, * )  Q,
integer  LDQ,
complex*16, dimension( ldq, * )  Z,
complex*16, dimension( ldqe, * )  QE,
integer  LDQE,
complex*16, dimension( * )  ALPHA,
complex*16, dimension( * )  BETA,
complex*16, dimension( * )  ALPHA1,
complex*16, dimension( * )  BETA1,
complex*16, dimension( * )  WORK,
integer  LWORK,
double precision, dimension( * )  RWORK,
double precision, dimension( * )  RESULT,
integer  INFO 
)

ZDRGEV

Purpose:
 ZDRGEV checks the nonsymmetric generalized eigenvalue problem driver
 routine ZGGEV.

 ZGGEV computes for a pair of n-by-n nonsymmetric matrices (A,B) the
 generalized eigenvalues and, optionally, the left and right
 eigenvectors.

 A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
 or a ratio  alpha/beta = w, such that A - w*B is singular.  It is
 usually represented as the pair (alpha,beta), as there is reasonable
 interpretation for beta=0, and even for both being zero.

 A right generalized eigenvector corresponding to a generalized
 eigenvalue  w  for a pair of matrices (A,B) is a vector r  such that
 (A - wB) * r = 0.  A left generalized eigenvector is a vector l such
 that l**H * (A - wB) = 0, where l**H is the conjugate-transpose of l.

 When ZDRGEV is called, a number of matrix "sizes" ("n's") and a
 number of matrix "types" are specified.  For each size ("n")
 and each type of matrix, a pair of matrices (A, B) will be generated
 and used for testing.  For each matrix pair, the following tests
 will be performed and compared with the threshold THRESH.

 Results from ZGGEV:

 (1)  max over all left eigenvalue/-vector pairs (alpha/beta,l) of

      | VL**H * (beta A - alpha B) |/( ulp max(|beta A|, |alpha B|) )

      where VL**H is the conjugate-transpose of VL.

 (2)  | |VL(i)| - 1 | / ulp and whether largest component real

      VL(i) denotes the i-th column of VL.

 (3)  max over all left eigenvalue/-vector pairs (alpha/beta,r) of

      | (beta A - alpha B) * VR | / ( ulp max(|beta A|, |alpha B|) )

 (4)  | |VR(i)| - 1 | / ulp and whether largest component real

      VR(i) denotes the i-th column of VR.

 (5)  W(full) = W(partial)
      W(full) denotes the eigenvalues computed when both l and r
      are also computed, and W(partial) denotes the eigenvalues
      computed when only W, only W and r, or only W and l are
      computed.

 (6)  VL(full) = VL(partial)
      VL(full) denotes the left eigenvectors computed when both l
      and r are computed, and VL(partial) denotes the result
      when only l is computed.

 (7)  VR(full) = VR(partial)
      VR(full) denotes the right eigenvectors computed when both l
      and r are also computed, and VR(partial) denotes the result
      when only l is computed.


 Test Matrices
 ---- --------

 The sizes of the test matrices are specified by an array
 NN(1:NSIZES); the value of each element NN(j) specifies one size.
 The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if
 DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
 Currently, the list of possible types is:

 (1)  ( 0, 0 )         (a pair of zero matrices)

 (2)  ( I, 0 )         (an identity and a zero matrix)

 (3)  ( 0, I )         (an identity and a zero matrix)

 (4)  ( I, I )         (a pair of identity matrices)

         t   t
 (5)  ( J , J  )       (a pair of transposed Jordan blocks)

                                     t                ( I   0  )
 (6)  ( X, Y )         where  X = ( J   0  )  and Y = (      t )
                                  ( 0   I  )          ( 0   J  )
                       and I is a k x k identity and J a (k+1)x(k+1)
                       Jordan block; k=(N-1)/2

 (7)  ( D, I )         where D is diag( 0, 1,..., N-1 ) (a diagonal
                       matrix with those diagonal entries.)
 (8)  ( I, D )

 (9)  ( big*D, small*I ) where "big" is near overflow and small=1/big

 (10) ( small*D, big*I )

 (11) ( big*I, small*D )

 (12) ( small*I, big*D )

 (13) ( big*D, big*I )

 (14) ( small*D, small*I )

 (15) ( D1, D2 )        where D1 is diag( 0, 0, 1, ..., N-3, 0 ) and
                        D2 is diag( 0, N-3, N-4,..., 1, 0, 0 )
           t   t
 (16) Q ( J , J ) Z     where Q and Z are random orthogonal matrices.

 (17) Q ( T1, T2 ) Z    where T1 and T2 are upper triangular matrices
                        with random O(1) entries above the diagonal
                        and diagonal entries diag(T1) =
                        ( 0, 0, 1, ..., N-3, 0 ) and diag(T2) =
                        ( 0, N-3, N-4,..., 1, 0, 0 )

 (18) Q ( T1, T2 ) Z    diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
                        diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
                        s = machine precision.

 (19) Q ( T1, T2 ) Z    diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 )

                                                        N-5
 (20) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, 1, a, ..., a   =s, 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )

 (21) Q ( T1, T2 ) Z    diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 )
                        diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
                        where r1,..., r(N-4) are random.

 (22) Q ( big*T1, small*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (23) Q ( small*T1, big*T2 ) Z    diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (24) Q ( small*T1, small*T2 ) Z  diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (25) Q ( big*T1, big*T2 ) Z      diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
                                  diag(T2) = ( 0, 1, ..., 1, 0, 0 )

 (26) Q ( T1, T2 ) Z     where T1 and T2 are random upper-triangular
                         matrices.
Parameters
[in]NSIZES
          NSIZES is INTEGER
          The number of sizes of matrices to use.  If it is zero,
          ZDRGES does nothing.  NSIZES >= 0.
[in]NN
          NN is INTEGER array, dimension (NSIZES)
          An array containing the sizes to be used for the matrices.
          Zero values will be skipped.  NN >= 0.
[in]NTYPES
          NTYPES is INTEGER
          The number of elements in DOTYPE.   If it is zero, ZDRGEV
          does nothing.  It must be at least zero.  If it is MAXTYP+1
          and NSIZES is 1, then an additional type, MAXTYP+1 is
          defined, which is to use whatever matrix is in A.  This
          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
          DOTYPE(MAXTYP+1) is .TRUE. .
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          If DOTYPE(j) is .TRUE., then for each size in NN a
          matrix of that size and of type j will be generated.
          If NTYPES is smaller than the maximum number of types
          defined (PARAMETER MAXTYP), then types NTYPES+1 through
          MAXTYP will not be generated. If NTYPES is larger
          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
          will be ignored.
[in,out]ISEED
          ISEED is INTEGER array, dimension (4)
          On entry ISEED specifies the seed of the random number
          generator. The array elements should be between 0 and 4095;
          if not they will be reduced mod 4096. Also, ISEED(4) must
          be odd.  The random number generator uses a linear
          congruential sequence limited to small integers, and so
          should produce machine independent random numbers. The
          values of ISEED are changed on exit, and can be used in the
          next call to ZDRGES to continue the same random number
          sequence.
[in]THRESH
          THRESH is DOUBLE PRECISION
          A test will count as "failed" if the "error", computed as
          described above, exceeds THRESH.  Note that the error is
          scaled to be O(1), so THRESH should be a reasonably small
          multiple of 1, e.g., 10 or 100.  In particular, it should
          not depend on the precision (single vs. double) or the size
          of the matrix.  It must be at least zero.
[in]NOUNIT
          NOUNIT is INTEGER
          The FORTRAN unit number for printing out error messages
          (e.g., if a routine returns IERR not equal to 0.)
[in,out]A
          A is COMPLEX*16 array, dimension(LDA, max(NN))
          Used to hold the original A matrix.  Used as input only
          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
          DOTYPE(MAXTYP+1)=.TRUE.
[in]LDA
          LDA is INTEGER
          The leading dimension of A, B, S, and T.
          It must be at least 1 and at least max( NN ).
[in,out]B
          B is COMPLEX*16 array, dimension(LDA, max(NN))
          Used to hold the original B matrix.  Used as input only
          if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
          DOTYPE(MAXTYP+1)=.TRUE.
[out]S
          S is COMPLEX*16 array, dimension (LDA, max(NN))
          The Schur form matrix computed from A by ZGGEV.  On exit, S
          contains the Schur form matrix corresponding to the matrix
          in A.
[out]T
          T is COMPLEX*16 array, dimension (LDA, max(NN))
          The upper triangular matrix computed from B by ZGGEV.
[out]Q
          Q is COMPLEX*16 array, dimension (LDQ, max(NN))
          The (left) eigenvectors matrix computed by ZGGEV.
[in]LDQ
          LDQ is INTEGER
          The leading dimension of Q and Z. It must
          be at least 1 and at least max( NN ).
[out]Z
          Z is COMPLEX*16 array, dimension( LDQ, max(NN) )
          The (right) orthogonal matrix computed by ZGGEV.
[out]QE
          QE is COMPLEX*16 array, dimension( LDQ, max(NN) )
          QE holds the computed right or left eigenvectors.
[in]LDQE
          LDQE is INTEGER
          The leading dimension of QE. LDQE >= max(1,max(NN)).
[out]ALPHA
          ALPHA is COMPLEX*16 array, dimension (max(NN))
[out]BETA
          BETA is COMPLEX*16 array, dimension (max(NN))

          The generalized eigenvalues of (A,B) computed by ZGGEV.
          ( ALPHAR(k)+ALPHAI(k)*i ) / BETA(k) is the k-th
          generalized eigenvalue of A and B.
[out]ALPHA1
          ALPHA1 is COMPLEX*16 array, dimension (max(NN))
[out]BETA1
          BETA1 is COMPLEX*16 array, dimension (max(NN))

          Like ALPHAR, ALPHAI, BETA, these arrays contain the
          eigenvalues of A and B, but those computed when ZGGEV only
          computes a partial eigendecomposition, i.e. not the
          eigenvalues and left and right eigenvectors.
[out]WORK
          WORK is COMPLEX*16 array, dimension (LWORK)
[in]LWORK
          LWORK is INTEGER
          The number of entries in WORK.  LWORK >= N*(N+1)
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension (8*N)
          Real workspace.
[out]RESULT
          RESULT is DOUBLE PRECISION array, dimension (2)
          The values computed by the tests described above.
          The values are currently limited to 1/ulp, to avoid overflow.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          > 0:  A routine returned an error code.  INFO is the
                absolute value of the INFO value returned.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 395 of file zdrgev.f.

399*
400* -- LAPACK test routine --
401* -- LAPACK is a software package provided by Univ. of Tennessee, --
402* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
403*
404* .. Scalar Arguments ..
405 INTEGER INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
406 $ NTYPES
407 DOUBLE PRECISION THRESH
408* ..
409* .. Array Arguments ..
410 LOGICAL DOTYPE( * )
411 INTEGER ISEED( 4 ), NN( * )
412 DOUBLE PRECISION RESULT( * ), RWORK( * )
413 COMPLEX*16 A( LDA, * ), ALPHA( * ), ALPHA1( * ),
414 $ B( LDA, * ), BETA( * ), BETA1( * ),
415 $ Q( LDQ, * ), QE( LDQE, * ), S( LDA, * ),
416 $ T( LDA, * ), WORK( * ), Z( LDQ, * )
417* ..
418*
419* =====================================================================
420*
421* .. Parameters ..
422 DOUBLE PRECISION ZERO, ONE
423 parameter( zero = 0.0d+0, one = 1.0d+0 )
424 COMPLEX*16 CZERO, CONE
425 parameter( czero = ( 0.0d+0, 0.0d+0 ),
426 $ cone = ( 1.0d+0, 0.0d+0 ) )
427 INTEGER MAXTYP
428 parameter( maxtyp = 26 )
429* ..
430* .. Local Scalars ..
431 LOGICAL BADNN
432 INTEGER I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
433 $ MAXWRK, MINWRK, MTYPES, N, N1, NB, NERRS,
434 $ NMATS, NMAX, NTESTT
435 DOUBLE PRECISION SAFMAX, SAFMIN, ULP, ULPINV
436 COMPLEX*16 CTEMP
437* ..
438* .. Local Arrays ..
439 LOGICAL LASIGN( MAXTYP ), LBSIGN( MAXTYP )
440 INTEGER IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
441 $ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
442 $ KBMAGN( MAXTYP ), KBTYPE( MAXTYP ),
443 $ KBZERO( MAXTYP ), KCLASS( MAXTYP ),
444 $ KTRIAN( MAXTYP ), KZ1( 6 ), KZ2( 6 )
445 DOUBLE PRECISION RMAGN( 0: 3 )
446* ..
447* .. External Functions ..
448 INTEGER ILAENV
449 DOUBLE PRECISION DLAMCH
450 COMPLEX*16 ZLARND
451 EXTERNAL ilaenv, dlamch, zlarnd
452* ..
453* .. External Subroutines ..
454 EXTERNAL alasvm, dlabad, xerbla, zget52, zggev, zlacpy,
456* ..
457* .. Intrinsic Functions ..
458 INTRINSIC abs, dble, dconjg, max, min, sign
459* ..
460* .. Data statements ..
461 DATA kclass / 15*1, 10*2, 1*3 /
462 DATA kz1 / 0, 1, 2, 1, 3, 3 /
463 DATA kz2 / 0, 0, 1, 2, 1, 1 /
464 DATA kadd / 0, 0, 0, 0, 3, 2 /
465 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
466 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
467 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
468 $ 1, 1, -4, 2, -4, 8*8, 0 /
469 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
470 $ 4*5, 4*3, 1 /
471 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
472 $ 4*6, 4*4, 1 /
473 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
474 $ 2, 1 /
475 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
476 $ 2, 1 /
477 DATA ktrian / 16*0, 10*1 /
478 DATA lasign / 6*.false., .true., .false., 2*.true.,
479 $ 2*.false., 3*.true., .false., .true.,
480 $ 3*.false., 5*.true., .false. /
481 DATA lbsign / 7*.false., .true., 2*.false.,
482 $ 2*.true., 2*.false., .true., .false., .true.,
483 $ 9*.false. /
484* ..
485* .. Executable Statements ..
486*
487* Check for errors
488*
489 info = 0
490*
491 badnn = .false.
492 nmax = 1
493 DO 10 j = 1, nsizes
494 nmax = max( nmax, nn( j ) )
495 IF( nn( j ).LT.0 )
496 $ badnn = .true.
497 10 CONTINUE
498*
499 IF( nsizes.LT.0 ) THEN
500 info = -1
501 ELSE IF( badnn ) THEN
502 info = -2
503 ELSE IF( ntypes.LT.0 ) THEN
504 info = -3
505 ELSE IF( thresh.LT.zero ) THEN
506 info = -6
507 ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
508 info = -9
509 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax ) THEN
510 info = -14
511 ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax ) THEN
512 info = -17
513 END IF
514*
515* Compute workspace
516* (Note: Comments in the code beginning "Workspace:" describe the
517* minimal amount of workspace needed at that point in the code,
518* as well as the preferred amount for good performance.
519* NB refers to the optimal block size for the immediately
520* following subroutine, as returned by ILAENV.
521*
522 minwrk = 1
523 IF( info.EQ.0 .AND. lwork.GE.1 ) THEN
524 minwrk = nmax*( nmax+1 )
525 nb = max( 1, ilaenv( 1, 'ZGEQRF', ' ', nmax, nmax, -1, -1 ),
526 $ ilaenv( 1, 'ZUNMQR', 'LC', nmax, nmax, nmax, -1 ),
527 $ ilaenv( 1, 'ZUNGQR', ' ', nmax, nmax, nmax, -1 ) )
528 maxwrk = max( 2*nmax, nmax*( nb+1 ), nmax*( nmax+1 ) )
529 work( 1 ) = maxwrk
530 END IF
531*
532 IF( lwork.LT.minwrk )
533 $ info = -23
534*
535 IF( info.NE.0 ) THEN
536 CALL xerbla( 'ZDRGEV', -info )
537 RETURN
538 END IF
539*
540* Quick return if possible
541*
542 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
543 $ RETURN
544*
545 ulp = dlamch( 'Precision' )
546 safmin = dlamch( 'Safe minimum' )
547 safmin = safmin / ulp
548 safmax = one / safmin
549 CALL dlabad( safmin, safmax )
550 ulpinv = one / ulp
551*
552* The values RMAGN(2:3) depend on N, see below.
553*
554 rmagn( 0 ) = zero
555 rmagn( 1 ) = one
556*
557* Loop over sizes, types
558*
559 ntestt = 0
560 nerrs = 0
561 nmats = 0
562*
563 DO 220 jsize = 1, nsizes
564 n = nn( jsize )
565 n1 = max( 1, n )
566 rmagn( 2 ) = safmax*ulp / dble( n1 )
567 rmagn( 3 ) = safmin*ulpinv*n1
568*
569 IF( nsizes.NE.1 ) THEN
570 mtypes = min( maxtyp, ntypes )
571 ELSE
572 mtypes = min( maxtyp+1, ntypes )
573 END IF
574*
575 DO 210 jtype = 1, mtypes
576 IF( .NOT.dotype( jtype ) )
577 $ GO TO 210
578 nmats = nmats + 1
579*
580* Save ISEED in case of an error.
581*
582 DO 20 j = 1, 4
583 ioldsd( j ) = iseed( j )
584 20 CONTINUE
585*
586* Generate test matrices A and B
587*
588* Description of control parameters:
589*
590* KZLASS: =1 means w/o rotation, =2 means w/ rotation,
591* =3 means random.
592* KATYPE: the "type" to be passed to ZLATM4 for computing A.
593* KAZERO: the pattern of zeros on the diagonal for A:
594* =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
595* =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
596* =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
597* non-zero entries.)
598* KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
599* =2: large, =3: small.
600* LASIGN: .TRUE. if the diagonal elements of A are to be
601* multiplied by a random magnitude 1 number.
602* KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B.
603* KTRIAN: =0: don't fill in the upper triangle, =1: do.
604* KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
605* RMAGN: used to implement KAMAGN and KBMAGN.
606*
607 IF( mtypes.GT.maxtyp )
608 $ GO TO 100
609 ierr = 0
610 IF( kclass( jtype ).LT.3 ) THEN
611*
612* Generate A (w/o rotation)
613*
614 IF( abs( katype( jtype ) ).EQ.3 ) THEN
615 in = 2*( ( n-1 ) / 2 ) + 1
616 IF( in.NE.n )
617 $ CALL zlaset( 'Full', n, n, czero, czero, a, lda )
618 ELSE
619 in = n
620 END IF
621 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
622 $ kz2( kazero( jtype ) ), lasign( jtype ),
623 $ rmagn( kamagn( jtype ) ), ulp,
624 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
625 $ iseed, a, lda )
626 iadd = kadd( kazero( jtype ) )
627 IF( iadd.GT.0 .AND. iadd.LE.n )
628 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
629*
630* Generate B (w/o rotation)
631*
632 IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
633 in = 2*( ( n-1 ) / 2 ) + 1
634 IF( in.NE.n )
635 $ CALL zlaset( 'Full', n, n, czero, czero, b, lda )
636 ELSE
637 in = n
638 END IF
639 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
640 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
641 $ rmagn( kbmagn( jtype ) ), one,
642 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
643 $ iseed, b, lda )
644 iadd = kadd( kbzero( jtype ) )
645 IF( iadd.NE.0 .AND. iadd.LE.n )
646 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
647*
648 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
649*
650* Include rotations
651*
652* Generate Q, Z as Householder transformations times
653* a diagonal matrix.
654*
655 DO 40 jc = 1, n - 1
656 DO 30 jr = jc, n
657 q( jr, jc ) = zlarnd( 3, iseed )
658 z( jr, jc ) = zlarnd( 3, iseed )
659 30 CONTINUE
660 CALL zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
661 $ work( jc ) )
662 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
663 q( jc, jc ) = cone
664 CALL zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
665 $ work( n+jc ) )
666 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
667 z( jc, jc ) = cone
668 40 CONTINUE
669 ctemp = zlarnd( 3, iseed )
670 q( n, n ) = cone
671 work( n ) = czero
672 work( 3*n ) = ctemp / abs( ctemp )
673 ctemp = zlarnd( 3, iseed )
674 z( n, n ) = cone
675 work( 2*n ) = czero
676 work( 4*n ) = ctemp / abs( ctemp )
677*
678* Apply the diagonal matrices
679*
680 DO 60 jc = 1, n
681 DO 50 jr = 1, n
682 a( jr, jc ) = work( 2*n+jr )*
683 $ dconjg( work( 3*n+jc ) )*
684 $ a( jr, jc )
685 b( jr, jc ) = work( 2*n+jr )*
686 $ dconjg( work( 3*n+jc ) )*
687 $ b( jr, jc )
688 50 CONTINUE
689 60 CONTINUE
690 CALL zunm2r( 'L', 'N', n, n, n-1, q, ldq, work, a,
691 $ lda, work( 2*n+1 ), ierr )
692 IF( ierr.NE.0 )
693 $ GO TO 90
694 CALL zunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
695 $ a, lda, work( 2*n+1 ), ierr )
696 IF( ierr.NE.0 )
697 $ GO TO 90
698 CALL zunm2r( 'L', 'N', n, n, n-1, q, ldq, work, b,
699 $ lda, work( 2*n+1 ), ierr )
700 IF( ierr.NE.0 )
701 $ GO TO 90
702 CALL zunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
703 $ b, lda, work( 2*n+1 ), ierr )
704 IF( ierr.NE.0 )
705 $ GO TO 90
706 END IF
707 ELSE
708*
709* Random matrices
710*
711 DO 80 jc = 1, n
712 DO 70 jr = 1, n
713 a( jr, jc ) = rmagn( kamagn( jtype ) )*
714 $ zlarnd( 4, iseed )
715 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
716 $ zlarnd( 4, iseed )
717 70 CONTINUE
718 80 CONTINUE
719 END IF
720*
721 90 CONTINUE
722*
723 IF( ierr.NE.0 ) THEN
724 WRITE( nounit, fmt = 9999 )'Generator', ierr, n, jtype,
725 $ ioldsd
726 info = abs( ierr )
727 RETURN
728 END IF
729*
730 100 CONTINUE
731*
732 DO 110 i = 1, 7
733 result( i ) = -one
734 110 CONTINUE
735*
736* Call ZGGEV to compute eigenvalues and eigenvectors.
737*
738 CALL zlacpy( ' ', n, n, a, lda, s, lda )
739 CALL zlacpy( ' ', n, n, b, lda, t, lda )
740 CALL zggev( 'V', 'V', n, s, lda, t, lda, alpha, beta, q,
741 $ ldq, z, ldq, work, lwork, rwork, ierr )
742 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
743 result( 1 ) = ulpinv
744 WRITE( nounit, fmt = 9999 )'ZGGEV1', ierr, n, jtype,
745 $ ioldsd
746 info = abs( ierr )
747 GO TO 190
748 END IF
749*
750* Do the tests (1) and (2)
751*
752 CALL zget52( .true., n, a, lda, b, lda, q, ldq, alpha, beta,
753 $ work, rwork, result( 1 ) )
754 IF( result( 2 ).GT.thresh ) THEN
755 WRITE( nounit, fmt = 9998 )'Left', 'ZGGEV1',
756 $ result( 2 ), n, jtype, ioldsd
757 END IF
758*
759* Do the tests (3) and (4)
760*
761 CALL zget52( .false., n, a, lda, b, lda, z, ldq, alpha,
762 $ beta, work, rwork, result( 3 ) )
763 IF( result( 4 ).GT.thresh ) THEN
764 WRITE( nounit, fmt = 9998 )'Right', 'ZGGEV1',
765 $ result( 4 ), n, jtype, ioldsd
766 END IF
767*
768* Do test (5)
769*
770 CALL zlacpy( ' ', n, n, a, lda, s, lda )
771 CALL zlacpy( ' ', n, n, b, lda, t, lda )
772 CALL zggev( 'N', 'N', n, s, lda, t, lda, alpha1, beta1, q,
773 $ ldq, z, ldq, work, lwork, rwork, ierr )
774 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
775 result( 1 ) = ulpinv
776 WRITE( nounit, fmt = 9999 )'ZGGEV2', ierr, n, jtype,
777 $ ioldsd
778 info = abs( ierr )
779 GO TO 190
780 END IF
781*
782 DO 120 j = 1, n
783 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
784 $ beta1( j ) )result( 5 ) = ulpinv
785 120 CONTINUE
786*
787* Do test (6): Compute eigenvalues and left eigenvectors,
788* and test them
789*
790 CALL zlacpy( ' ', n, n, a, lda, s, lda )
791 CALL zlacpy( ' ', n, n, b, lda, t, lda )
792 CALL zggev( 'V', 'N', n, s, lda, t, lda, alpha1, beta1, qe,
793 $ ldqe, z, ldq, work, lwork, rwork, ierr )
794 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
795 result( 1 ) = ulpinv
796 WRITE( nounit, fmt = 9999 )'ZGGEV3', ierr, n, jtype,
797 $ ioldsd
798 info = abs( ierr )
799 GO TO 190
800 END IF
801*
802 DO 130 j = 1, n
803 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
804 $ beta1( j ) )result( 6 ) = ulpinv
805 130 CONTINUE
806*
807 DO 150 j = 1, n
808 DO 140 jc = 1, n
809 IF( q( j, jc ).NE.qe( j, jc ) )
810 $ result( 6 ) = ulpinv
811 140 CONTINUE
812 150 CONTINUE
813*
814* Do test (7): Compute eigenvalues and right eigenvectors,
815* and test them
816*
817 CALL zlacpy( ' ', n, n, a, lda, s, lda )
818 CALL zlacpy( ' ', n, n, b, lda, t, lda )
819 CALL zggev( 'N', 'V', n, s, lda, t, lda, alpha1, beta1, q,
820 $ ldq, qe, ldqe, work, lwork, rwork, ierr )
821 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
822 result( 1 ) = ulpinv
823 WRITE( nounit, fmt = 9999 )'ZGGEV4', ierr, n, jtype,
824 $ ioldsd
825 info = abs( ierr )
826 GO TO 190
827 END IF
828*
829 DO 160 j = 1, n
830 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
831 $ beta1( j ) )result( 7 ) = ulpinv
832 160 CONTINUE
833*
834 DO 180 j = 1, n
835 DO 170 jc = 1, n
836 IF( z( j, jc ).NE.qe( j, jc ) )
837 $ result( 7 ) = ulpinv
838 170 CONTINUE
839 180 CONTINUE
840*
841* End of Loop -- Check for RESULT(j) > THRESH
842*
843 190 CONTINUE
844*
845 ntestt = ntestt + 7
846*
847* Print out tests which fail.
848*
849 DO 200 jr = 1, 7
850 IF( result( jr ).GE.thresh ) THEN
851*
852* If this is the first test to fail,
853* print a header to the data file.
854*
855 IF( nerrs.EQ.0 ) THEN
856 WRITE( nounit, fmt = 9997 )'ZGV'
857*
858* Matrix types
859*
860 WRITE( nounit, fmt = 9996 )
861 WRITE( nounit, fmt = 9995 )
862 WRITE( nounit, fmt = 9994 )'Orthogonal'
863*
864* Tests performed
865*
866 WRITE( nounit, fmt = 9993 )
867*
868 END IF
869 nerrs = nerrs + 1
870 IF( result( jr ).LT.10000.0d0 ) THEN
871 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
872 $ result( jr )
873 ELSE
874 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
875 $ result( jr )
876 END IF
877 END IF
878 200 CONTINUE
879*
880 210 CONTINUE
881 220 CONTINUE
882*
883* Summary
884*
885 CALL alasvm( 'ZGV', nounit, nerrs, ntestt, 0 )
886*
887 work( 1 ) = maxwrk
888*
889 RETURN
890*
891 9999 FORMAT( ' ZDRGEV: ', a, ' returned INFO=', i6, '.', / 3x, 'N=',
892 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
893*
894 9998 FORMAT( ' ZDRGEV: ', a, ' Eigenvectors from ', a, ' incorrectly ',
895 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 3x,
896 $ 'N=', i4, ', JTYPE=', i3, ', ISEED=(', 3( i4, ',' ), i5,
897 $ ')' )
898*
899 9997 FORMAT( / 1x, a3, ' -- Complex Generalized eigenvalue problem ',
900 $ 'driver' )
901*
902 9996 FORMAT( ' Matrix types (see ZDRGEV for details): ' )
903*
904 9995 FORMAT( ' Special Matrices:', 23x,
905 $ '(J''=transposed Jordan block)',
906 $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
907 $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
908 $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
909 $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
910 $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
911 $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
912 9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
913 $ / ' 16=Transposed Jordan Blocks 19=geometric ',
914 $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
915 $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
916 $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
917 $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
918 $ '23=(small,large) 24=(small,small) 25=(large,large)',
919 $ / ' 26=random O(1) matrices.' )
920*
921 9993 FORMAT( / ' Tests performed: ',
922 $ / ' 1 = max | ( b A - a B )''*l | / const.,',
923 $ / ' 2 = | |VR(i)| - 1 | / ulp,',
924 $ / ' 3 = max | ( b A - a B )*r | / const.',
925 $ / ' 4 = | |VL(i)| - 1 | / ulp,',
926 $ / ' 5 = 0 if W same no matter if r or l computed,',
927 $ / ' 6 = 0 if l same no matter if l computed,',
928 $ / ' 7 = 0 if r same no matter if r computed,', / 1x )
929 9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
930 $ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
931 9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
932 $ 4( i4, ',' ), ' result ', i2, ' is', 1p, d10.3 )
933*
934* End of ZDRGEV
935*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:73
subroutine zget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA, WORK, RWORK, RESULT)
ZGET52
Definition: zget52.f:162
subroutine zlatm4(ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
ZLATM4
Definition: zlatm4.f:171
complex *16 function zlarnd(IDIST, ISEED)
ZLARND
Definition: zlarnd.f:75
subroutine zggev(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
ZGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition: zggev.f:217
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 zlarfg(N, ALPHA, X, INCX, TAU)
ZLARFG generates an elementary reflector (Householder matrix).
Definition: zlarfg.f:106
subroutine zunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
Definition: zunm2r.f:159
Here is the call graph for this function:
Here is the caller graph for this function: