LAPACK 3.12.1
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  () and a
!> number of matrix  are specified.  For each size ()
!> 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  are specified by a logical array DOTYPE( 1:NTYPES ); if
!> DOTYPE(j) is .TRUE., then matrix type  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  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  if the , 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, 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 ulpinv = one / ulp
550*
551* The values RMAGN(2:3) depend on N, see below.
552*
553 rmagn( 0 ) = zero
554 rmagn( 1 ) = one
555*
556* Loop over sizes, types
557*
558 ntestt = 0
559 nerrs = 0
560 nmats = 0
561*
562 DO 220 jsize = 1, nsizes
563 n = nn( jsize )
564 n1 = max( 1, n )
565 rmagn( 2 ) = safmax*ulp / dble( n1 )
566 rmagn( 3 ) = safmin*ulpinv*n1
567*
568 IF( nsizes.NE.1 ) THEN
569 mtypes = min( maxtyp, ntypes )
570 ELSE
571 mtypes = min( maxtyp+1, ntypes )
572 END IF
573*
574 DO 210 jtype = 1, mtypes
575 IF( .NOT.dotype( jtype ) )
576 $ GO TO 210
577 nmats = nmats + 1
578*
579* Save ISEED in case of an error.
580*
581 DO 20 j = 1, 4
582 ioldsd( j ) = iseed( j )
583 20 CONTINUE
584*
585* Generate test matrices A and B
586*
587* Description of control parameters:
588*
589* KZLASS: =1 means w/o rotation, =2 means w/ rotation,
590* =3 means random.
591* KATYPE: the "type" to be passed to ZLATM4 for computing A.
592* KAZERO: the pattern of zeros on the diagonal for A:
593* =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
594* =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
595* =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
596* non-zero entries.)
597* KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
598* =2: large, =3: small.
599* LASIGN: .TRUE. if the diagonal elements of A are to be
600* multiplied by a random magnitude 1 number.
601* KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B.
602* KTRIAN: =0: don't fill in the upper triangle, =1: do.
603* KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
604* RMAGN: used to implement KAMAGN and KBMAGN.
605*
606 IF( mtypes.GT.maxtyp )
607 $ GO TO 100
608 ierr = 0
609 IF( kclass( jtype ).LT.3 ) THEN
610*
611* Generate A (w/o rotation)
612*
613 IF( abs( katype( jtype ) ).EQ.3 ) THEN
614 in = 2*( ( n-1 ) / 2 ) + 1
615 IF( in.NE.n )
616 $ CALL zlaset( 'Full', n, n, czero, czero, a, lda )
617 ELSE
618 in = n
619 END IF
620 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
621 $ kz2( kazero( jtype ) ), lasign( jtype ),
622 $ rmagn( kamagn( jtype ) ), ulp,
623 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
624 $ iseed, a, lda )
625 iadd = kadd( kazero( jtype ) )
626 IF( iadd.GT.0 .AND. iadd.LE.n )
627 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
628*
629* Generate B (w/o rotation)
630*
631 IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
632 in = 2*( ( n-1 ) / 2 ) + 1
633 IF( in.NE.n )
634 $ CALL zlaset( 'Full', n, n, czero, czero, b, lda )
635 ELSE
636 in = n
637 END IF
638 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
639 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
640 $ rmagn( kbmagn( jtype ) ), one,
641 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
642 $ iseed, b, lda )
643 iadd = kadd( kbzero( jtype ) )
644 IF( iadd.NE.0 .AND. iadd.LE.n )
645 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
646*
647 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
648*
649* Include rotations
650*
651* Generate Q, Z as Householder transformations times
652* a diagonal matrix.
653*
654 DO 40 jc = 1, n - 1
655 DO 30 jr = jc, n
656 q( jr, jc ) = zlarnd( 3, iseed )
657 z( jr, jc ) = zlarnd( 3, iseed )
658 30 CONTINUE
659 CALL zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
660 $ work( jc ) )
661 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
662 q( jc, jc ) = cone
663 CALL zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
664 $ work( n+jc ) )
665 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
666 z( jc, jc ) = cone
667 40 CONTINUE
668 ctemp = zlarnd( 3, iseed )
669 q( n, n ) = cone
670 work( n ) = czero
671 work( 3*n ) = ctemp / abs( ctemp )
672 ctemp = zlarnd( 3, iseed )
673 z( n, n ) = cone
674 work( 2*n ) = czero
675 work( 4*n ) = ctemp / abs( ctemp )
676*
677* Apply the diagonal matrices
678*
679 DO 60 jc = 1, n
680 DO 50 jr = 1, n
681 a( jr, jc ) = work( 2*n+jr )*
682 $ dconjg( work( 3*n+jc ) )*
683 $ a( jr, jc )
684 b( jr, jc ) = work( 2*n+jr )*
685 $ dconjg( work( 3*n+jc ) )*
686 $ b( jr, jc )
687 50 CONTINUE
688 60 CONTINUE
689 CALL zunm2r( 'L', 'N', n, n, n-1, q, ldq, work, a,
690 $ lda, work( 2*n+1 ), ierr )
691 IF( ierr.NE.0 )
692 $ GO TO 90
693 CALL zunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
694 $ a, lda, work( 2*n+1 ), ierr )
695 IF( ierr.NE.0 )
696 $ GO TO 90
697 CALL zunm2r( 'L', 'N', n, n, n-1, q, ldq, work, b,
698 $ lda, work( 2*n+1 ), ierr )
699 IF( ierr.NE.0 )
700 $ GO TO 90
701 CALL zunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
702 $ b, lda, work( 2*n+1 ), ierr )
703 IF( ierr.NE.0 )
704 $ GO TO 90
705 END IF
706 ELSE
707*
708* Random matrices
709*
710 DO 80 jc = 1, n
711 DO 70 jr = 1, n
712 a( jr, jc ) = rmagn( kamagn( jtype ) )*
713 $ zlarnd( 4, iseed )
714 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
715 $ zlarnd( 4, iseed )
716 70 CONTINUE
717 80 CONTINUE
718 END IF
719*
720 90 CONTINUE
721*
722 IF( ierr.NE.0 ) THEN
723 WRITE( nounit, fmt = 9999 )'Generator', ierr, n, jtype,
724 $ ioldsd
725 info = abs( ierr )
726 RETURN
727 END IF
728*
729 100 CONTINUE
730*
731 DO 110 i = 1, 7
732 result( i ) = -one
733 110 CONTINUE
734*
735* Call ZGGEV to compute eigenvalues and eigenvectors.
736*
737 CALL zlacpy( ' ', n, n, a, lda, s, lda )
738 CALL zlacpy( ' ', n, n, b, lda, t, lda )
739 CALL zggev( 'V', 'V', n, s, lda, t, lda, alpha, beta, q,
740 $ ldq, z, ldq, work, lwork, rwork, ierr )
741 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
742 result( 1 ) = ulpinv
743 WRITE( nounit, fmt = 9999 )'ZGGEV1', ierr, n, jtype,
744 $ ioldsd
745 info = abs( ierr )
746 GO TO 190
747 END IF
748*
749* Do the tests (1) and (2)
750*
751 CALL zget52( .true., n, a, lda, b, lda, q, ldq, alpha, beta,
752 $ work, rwork, result( 1 ) )
753 IF( result( 2 ).GT.thresh ) THEN
754 WRITE( nounit, fmt = 9998 )'Left', 'ZGGEV1',
755 $ result( 2 ), n, jtype, ioldsd
756 END IF
757*
758* Do the tests (3) and (4)
759*
760 CALL zget52( .false., n, a, lda, b, lda, z, ldq, alpha,
761 $ beta, work, rwork, result( 3 ) )
762 IF( result( 4 ).GT.thresh ) THEN
763 WRITE( nounit, fmt = 9998 )'Right', 'ZGGEV1',
764 $ result( 4 ), n, jtype, ioldsd
765 END IF
766*
767* Do test (5)
768*
769 CALL zlacpy( ' ', n, n, a, lda, s, lda )
770 CALL zlacpy( ' ', n, n, b, lda, t, lda )
771 CALL zggev( 'N', 'N', n, s, lda, t, lda, alpha1, beta1, q,
772 $ ldq, z, ldq, work, lwork, rwork, ierr )
773 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
774 result( 1 ) = ulpinv
775 WRITE( nounit, fmt = 9999 )'ZGGEV2', ierr, n, jtype,
776 $ ioldsd
777 info = abs( ierr )
778 GO TO 190
779 END IF
780*
781 DO 120 j = 1, n
782 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
783 $ beta1( j ) )result( 5 ) = ulpinv
784 120 CONTINUE
785*
786* Do test (6): Compute eigenvalues and left eigenvectors,
787* and test them
788*
789 CALL zlacpy( ' ', n, n, a, lda, s, lda )
790 CALL zlacpy( ' ', n, n, b, lda, t, lda )
791 CALL zggev( 'V', 'N', n, s, lda, t, lda, alpha1, beta1, qe,
792 $ ldqe, z, ldq, work, lwork, rwork, ierr )
793 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
794 result( 1 ) = ulpinv
795 WRITE( nounit, fmt = 9999 )'ZGGEV3', ierr, n, jtype,
796 $ ioldsd
797 info = abs( ierr )
798 GO TO 190
799 END IF
800*
801 DO 130 j = 1, n
802 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
803 $ beta1( j ) )result( 6 ) = ulpinv
804 130 CONTINUE
805*
806 DO 150 j = 1, n
807 DO 140 jc = 1, n
808 IF( q( j, jc ).NE.qe( j, jc ) )
809 $ result( 6 ) = ulpinv
810 140 CONTINUE
811 150 CONTINUE
812*
813* Do test (7): Compute eigenvalues and right eigenvectors,
814* and test them
815*
816 CALL zlacpy( ' ', n, n, a, lda, s, lda )
817 CALL zlacpy( ' ', n, n, b, lda, t, lda )
818 CALL zggev( 'N', 'V', n, s, lda, t, lda, alpha1, beta1, q,
819 $ ldq, qe, ldqe, work, lwork, rwork, ierr )
820 IF( ierr.NE.0 .AND. ierr.NE.n+1 ) THEN
821 result( 1 ) = ulpinv
822 WRITE( nounit, fmt = 9999 )'ZGGEV4', ierr, n, jtype,
823 $ ioldsd
824 info = abs( ierr )
825 GO TO 190
826 END IF
827*
828 DO 160 j = 1, n
829 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
830 $ beta1( j ) )result( 7 ) = ulpinv
831 160 CONTINUE
832*
833 DO 180 j = 1, n
834 DO 170 jc = 1, n
835 IF( z( j, jc ).NE.qe( j, jc ) )
836 $ result( 7 ) = ulpinv
837 170 CONTINUE
838 180 CONTINUE
839*
840* End of Loop -- Check for RESULT(j) > THRESH
841*
842 190 CONTINUE
843*
844 ntestt = ntestt + 7
845*
846* Print out tests which fail.
847*
848 DO 200 jr = 1, 7
849 IF( result( jr ).GE.thresh ) THEN
850*
851* If this is the first test to fail,
852* print a header to the data file.
853*
854 IF( nerrs.EQ.0 ) THEN
855 WRITE( nounit, fmt = 9997 )'ZGV'
856*
857* Matrix types
858*
859 WRITE( nounit, fmt = 9996 )
860 WRITE( nounit, fmt = 9995 )
861 WRITE( nounit, fmt = 9994 )'Orthogonal'
862*
863* Tests performed
864*
865 WRITE( nounit, fmt = 9993 )
866*
867 END IF
868 nerrs = nerrs + 1
869 IF( result( jr ).LT.10000.0d0 ) THEN
870 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
871 $ result( jr )
872 ELSE
873 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
874 $ result( jr )
875 END IF
876 END IF
877 200 CONTINUE
878*
879 210 CONTINUE
880 220 CONTINUE
881*
882* Summary
883*
884 CALL alasvm( 'ZGV', nounit, nerrs, ntestt, 0 )
885*
886 work( 1 ) = maxwrk
887*
888 RETURN
889*
890 9999 FORMAT( ' ZDRGEV: ', a, ' returned INFO=', i6, '.', / 3x, 'N=',
891 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
892*
893 9998 FORMAT( ' ZDRGEV: ', a, ' Eigenvectors from ', a, ' incorrectly ',
894 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 3x,
895 $ 'N=', i4, ', JTYPE=', i3, ', ISEED=(', 3( i4, ',' ), i5,
896 $ ')' )
897*
898 9997 FORMAT( / 1x, a3, ' -- Complex Generalized eigenvalue problem ',
899 $ 'driver' )
900*
901 9996 FORMAT( ' Matrix types (see ZDRGEV for details): ' )
902*
903 9995 FORMAT( ' Special Matrices:', 23x,
904 $ '(J''=transposed Jordan block)',
905 $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
906 $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
907 $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
908 $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
909 $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
910 $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
911 9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
912 $ / ' 16=Transposed Jordan Blocks 19=geometric ',
913 $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
914 $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
915 $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
916 $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
917 $ '23=(small,large) 24=(small,small) 25=(large,large)',
918 $ / ' 26=random O(1) matrices.' )
919*
920 9993 FORMAT( / ' Tests performed: ',
921 $ / ' 1 = max | ( b A - a B )''*l | / const.,',
922 $ / ' 2 = | |VR(i)| - 1 | / ulp,',
923 $ / ' 3 = max | ( b A - a B )*r | / const.',
924 $ / ' 4 = | |VL(i)| - 1 | / ulp,',
925 $ / ' 5 = 0 if W same no matter if r or l computed,',
926 $ / ' 6 = 0 if l same no matter if l computed,',
927 $ / ' 7 = 0 if r same no matter if r computed,', / 1x )
928 9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
929 $ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
930 9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
931 $ 4( i4, ',' ), ' result ', i2, ' is', 1p, d10.3 )
932*
933* End of ZDRGEV
934*
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
subroutine xerbla(srname, info)
Definition cblat2.f:3285
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:215
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:160
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:101
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
subroutine zlarfg(n, alpha, x, incx, tau)
ZLARFG generates an elementary reflector (Householder matrix).
Definition zlarfg.f:104
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:104
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:157
subroutine zget52(left, n, a, lda, b, ldb, e, lde, alpha, beta, work, rwork, result)
ZGET52
Definition zget52.f:162
complex *16 function zlarnd(idist, iseed)
ZLARND
Definition zlarnd.f:75
subroutine zlatm4(itype, n, nz1, nz2, rsign, amagn, rcond, triang, idist, iseed, a, lda)
ZLATM4
Definition zlatm4.f:171
Here is the call graph for this function:
Here is the caller graph for this function: