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

◆ ddrvev()

subroutine ddrvev ( integer  nsizes,
integer, dimension( * )  nn,
integer  ntypes,
logical, dimension( * )  dotype,
integer, dimension( 4 )  iseed,
double precision  thresh,
integer  nounit,
double precision, dimension( lda, * )  a,
integer  lda,
double precision, dimension( lda, * )  h,
double precision, dimension( * )  wr,
double precision, dimension( * )  wi,
double precision, dimension( * )  wr1,
double precision, dimension( * )  wi1,
double precision, dimension( ldvl, * )  vl,
integer  ldvl,
double precision, dimension( ldvr, * )  vr,
integer  ldvr,
double precision, dimension( ldlre, * )  lre,
integer  ldlre,
double precision, dimension( 7 )  result,
double precision, dimension( * )  work,
integer  nwork,
integer, dimension( * )  iwork,
integer  info 
)

DDRVEV

Purpose:
    DDRVEV  checks the nonsymmetric eigenvalue problem driver DGEEV.

    When DDRVEV 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, one matrix will be generated and used
    to test the nonsymmetric eigenroutines.  For each matrix, 7
    tests will be performed:

    (1)     | A * VR - VR * W | / ( n |A| ulp )

      Here VR is the matrix of unit right eigenvectors.
      W is a block diagonal matrix, with a 1x1 block for each
      real eigenvalue and a 2x2 block for each complex conjugate
      pair.  If eigenvalues j and j+1 are a complex conjugate pair,
      so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the
      2 x 2 block corresponding to the pair will be:

              (  wr  wi  )
              ( -wi  wr  )

      Such a block multiplying an n x 2 matrix  ( ur ui ) on the
      right will be the same as multiplying  ur + i*ui  by  wr + i*wi.

    (2)     | A**H * VL - VL * W**H | / ( n |A| ulp )

      Here VL is the matrix of unit left eigenvectors, A**H is the
      conjugate transpose of A, and W is as above.

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

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

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

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

    (5)     W(full) = W(partial)

      W(full) denotes the eigenvalues computed when both VR and VL
      are also computed, and W(partial) denotes the eigenvalues
      computed when only W, only W and VR, or only W and VL are
      computed.

    (6)     VR(full) = VR(partial)

      VR(full) denotes the right eigenvectors computed when both VR
      and VL are computed, and VR(partial) denotes the result
      when only VR is computed.

     (7)     VL(full) = VL(partial)

      VL(full) denotes the left eigenvectors computed when both VR
      and VL are also computed, and VL(partial) denotes the result
      when only VL is computed.

    The "sizes" 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)  The zero matrix.
    (2)  The identity matrix.
    (3)  A (transposed) Jordan block, with 1's on the diagonal.

    (4)  A diagonal matrix with evenly spaced entries
         1, ..., ULP  and random signs.
         (ULP = (first number larger than 1) - 1 )
    (5)  A diagonal matrix with geometrically spaced entries
         1, ..., ULP  and random signs.
    (6)  A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
         and random signs.

    (7)  Same as (4), but multiplied by a constant near
         the overflow threshold
    (8)  Same as (4), but multiplied by a constant near
         the underflow threshold

    (9)  A matrix of the form  U' T U, where U is orthogonal and
         T has evenly spaced entries 1, ..., ULP with random signs
         on the diagonal and random O(1) entries in the upper
         triangle.

    (10) A matrix of the form  U' T U, where U is orthogonal and
         T has geometrically spaced entries 1, ..., ULP with random
         signs on the diagonal and random O(1) entries in the upper
         triangle.

    (11) A matrix of the form  U' T U, where U is orthogonal and
         T has "clustered" entries 1, ULP,..., ULP with random
         signs on the diagonal and random O(1) entries in the upper
         triangle.

    (12) A matrix of the form  U' T U, where U is orthogonal and
         T has real or complex conjugate paired eigenvalues randomly
         chosen from ( ULP, 1 ) and random O(1) entries in the upper
         triangle.

    (13) A matrix of the form  X' T X, where X has condition
         SQRT( ULP ) and T has evenly spaced entries 1, ..., ULP
         with random signs on the diagonal and random O(1) entries
         in the upper triangle.

    (14) A matrix of the form  X' T X, where X has condition
         SQRT( ULP ) and T has geometrically spaced entries
         1, ..., ULP with random signs on the diagonal and random
         O(1) entries in the upper triangle.

    (15) A matrix of the form  X' T X, where X has condition
         SQRT( ULP ) and T has "clustered" entries 1, ULP,..., ULP
         with random signs on the diagonal and random O(1) entries
         in the upper triangle.

    (16) A matrix of the form  X' T X, where X has condition
         SQRT( ULP ) and T has real or complex conjugate paired
         eigenvalues randomly chosen from ( ULP, 1 ) and random
         O(1) entries in the upper triangle.

    (17) Same as (16), but multiplied by a constant
         near the overflow threshold
    (18) Same as (16), but multiplied by a constant
         near the underflow threshold

    (19) Nonsymmetric matrix with random entries chosen from (-1,1).
         If N is at least 4, all entries in first two rows and last
         row, and first column and last two columns are zero.
    (20) Same as (19), but multiplied by a constant
         near the overflow threshold
    (21) Same as (19), but multiplied by a constant
         near the underflow threshold
Parameters
[in]NSIZES
          NSIZES is INTEGER
          The number of sizes of matrices to use.  If it is zero,
          DDRVEV does nothing.  It must be at least zero.
[in]NN
          NN is INTEGER array, dimension (NSIZES)
          An array containing the sizes to be used for the matrices.
          Zero values will be skipped.  The values must be at least
          zero.
[in]NTYPES
          NTYPES is INTEGER
          The number of elements in DOTYPE.   If it is zero, DDRVEV
          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 DDRVEV 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 INFO not equal to 0.)
[out]A
          A is DOUBLE PRECISION array, dimension (LDA, max(NN))
          Used to hold the matrix whose eigenvalues are to be
          computed.  On exit, A contains the last matrix actually used.
[in]LDA
          LDA is INTEGER
          The leading dimension of A, and H. LDA must be at
          least 1 and at least max(NN).
[out]H
          H is DOUBLE PRECISION array, dimension (LDA, max(NN))
          Another copy of the test matrix A, modified by DGEEV.
[out]WR
          WR is DOUBLE PRECISION array, dimension (max(NN))
[out]WI
          WI is DOUBLE PRECISION array, dimension (max(NN))

          The real and imaginary parts of the eigenvalues of A.
          On exit, WR + WI*i are the eigenvalues of the matrix in A.
[out]WR1
          WR1 is DOUBLE PRECISION array, dimension (max(NN))
[out]WI1
          WI1 is DOUBLE PRECISION array, dimension (max(NN))

          Like WR, WI, these arrays contain the eigenvalues of A,
          but those computed when DGEEV only computes a partial
          eigendecomposition, i.e. not the eigenvalues and left
          and right eigenvectors.
[out]VL
          VL is DOUBLE PRECISION array, dimension (LDVL, max(NN))
          VL holds the computed left eigenvectors.
[in]LDVL
          LDVL is INTEGER
          Leading dimension of VL. Must be at least max(1,max(NN)).
[out]VR
          VR is DOUBLE PRECISION array, dimension (LDVR, max(NN))
          VR holds the computed right eigenvectors.
[in]LDVR
          LDVR is INTEGER
          Leading dimension of VR. Must be at least max(1,max(NN)).
[out]LRE
          LRE is DOUBLE PRECISION array, dimension (LDLRE,max(NN))
          LRE holds the computed right or left eigenvectors.
[in]LDLRE
          LDLRE is INTEGER
          Leading dimension of LRE. Must be at least max(1,max(NN)).
[out]RESULT
          RESULT is DOUBLE PRECISION array, dimension (7)
          The values computed by the seven tests described above.
          The values are currently limited to 1/ulp, to avoid overflow.
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (NWORK)
[in]NWORK
          NWORK is INTEGER
          The number of entries in WORK.  This must be at least
          5*NN(j)+2*NN(j)**2 for all j.
[out]IWORK
          IWORK is INTEGER array, dimension (max(NN))
[out]INFO
          INFO is INTEGER
          If 0, then everything ran OK.
           -1: NSIZES < 0
           -2: Some NN(j) < 0
           -3: NTYPES < 0
           -6: THRESH < 0
           -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
          -16: LDVL < 1 or LDVL < NMAX, where NMAX is max( NN(j) ).
          -18: LDVR < 1 or LDVR < NMAX, where NMAX is max( NN(j) ).
          -20: LDLRE < 1 or LDLRE < NMAX, where NMAX is max( NN(j) ).
          -23: NWORK too small.
          If  DLATMR, SLATMS, SLATME or DGEEV returns an error code,
              the absolute value of it is returned.

-----------------------------------------------------------------------

     Some Local Variables and Parameters:
     ---- ----- --------- --- ----------

     ZERO, ONE       Real 0 and 1.
     MAXTYP          The number of types defined.
     NMAX            Largest value in NN.
     NERRS           The number of tests which have exceeded THRESH
     COND, CONDS,
     IMODE           Values to be passed to the matrix generators.
     ANORM           Norm of A; passed to matrix generators.

     OVFL, UNFL      Overflow and underflow thresholds.
     ULP, ULPINV     Finest relative precision and its inverse.
     RTULP, RTULPI   Square roots of the previous 4 values.

             The following four arrays decode JTYPE:
     KTYPE(j)        The general type (1-10) for type "j".
     KMODE(j)        The MODE value to be passed to the matrix
                     generator for type "j".
     KMAGN(j)        The order of magnitude ( O(1),
                     O(overflow^(1/2) ), O(underflow^(1/2) )
     KCONDS(j)       Selectw whether CONDS is to be 1 or
                     1/sqrt(ulp).  (0 means irrelevant.)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 402 of file ddrvev.f.

406*
407* -- LAPACK test routine --
408* -- LAPACK is a software package provided by Univ. of Tennessee, --
409* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
410*
411* .. Scalar Arguments ..
412 INTEGER INFO, LDA, LDLRE, LDVL, LDVR, NOUNIT, NSIZES,
413 $ NTYPES, NWORK
414 DOUBLE PRECISION THRESH
415* ..
416* .. Array Arguments ..
417 LOGICAL DOTYPE( * )
418 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
419 DOUBLE PRECISION A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
420 $ RESULT( 7 ), VL( LDVL, * ), VR( LDVR, * ),
421 $ WI( * ), WI1( * ), WORK( * ), WR( * ), WR1( * )
422* ..
423*
424* =====================================================================
425*
426* .. Parameters ..
427 DOUBLE PRECISION ZERO, ONE
428 parameter( zero = 0.0d0, one = 1.0d0 )
429 DOUBLE PRECISION TWO
430 parameter( two = 2.0d0 )
431 INTEGER MAXTYP
432 parameter( maxtyp = 21 )
433* ..
434* .. Local Scalars ..
435 LOGICAL BADNN
436 CHARACTER*3 PATH
437 INTEGER IINFO, IMODE, ITYPE, IWK, J, JCOL, JJ, JSIZE,
438 $ JTYPE, MTYPES, N, NERRS, NFAIL, NMAX, NNWORK,
439 $ NTEST, NTESTF, NTESTT
440 DOUBLE PRECISION ANORM, COND, CONDS, OVFL, RTULP, RTULPI, TNRM,
441 $ ULP, ULPINV, UNFL, VMX, VRMX, VTST
442* ..
443* .. Local Arrays ..
444 CHARACTER ADUMMA( 1 )
445 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
446 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
447 $ KTYPE( MAXTYP )
448 DOUBLE PRECISION DUM( 1 ), RES( 2 )
449* ..
450* .. External Functions ..
451 DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2
452 EXTERNAL dlamch, dlapy2, dnrm2
453* ..
454* .. External Subroutines ..
455 EXTERNAL dgeev, dget22, dlacpy, dlaset, dlasum,
457* ..
458* .. Intrinsic Functions ..
459 INTRINSIC abs, max, min, sqrt
460* ..
461* .. Data statements ..
462 DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
463 DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
464 $ 3, 1, 2, 3 /
465 DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
466 $ 1, 5, 5, 5, 4, 3, 1 /
467 DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
468* ..
469* .. Executable Statements ..
470*
471 path( 1: 1 ) = 'Double precision'
472 path( 2: 3 ) = 'EV'
473*
474* Check for errors
475*
476 ntestt = 0
477 ntestf = 0
478 info = 0
479*
480* Important constants
481*
482 badnn = .false.
483 nmax = 0
484 DO 10 j = 1, nsizes
485 nmax = max( nmax, nn( j ) )
486 IF( nn( j ).LT.0 )
487 $ badnn = .true.
488 10 CONTINUE
489*
490* Check for errors
491*
492 IF( nsizes.LT.0 ) THEN
493 info = -1
494 ELSE IF( badnn ) THEN
495 info = -2
496 ELSE IF( ntypes.LT.0 ) THEN
497 info = -3
498 ELSE IF( thresh.LT.zero ) THEN
499 info = -6
500 ELSE IF( nounit.LE.0 ) THEN
501 info = -7
502 ELSE IF( lda.LT.1 .OR. lda.LT.nmax ) THEN
503 info = -9
504 ELSE IF( ldvl.LT.1 .OR. ldvl.LT.nmax ) THEN
505 info = -16
506 ELSE IF( ldvr.LT.1 .OR. ldvr.LT.nmax ) THEN
507 info = -18
508 ELSE IF( ldlre.LT.1 .OR. ldlre.LT.nmax ) THEN
509 info = -20
510 ELSE IF( 5*nmax+2*nmax**2.GT.nwork ) THEN
511 info = -23
512 END IF
513*
514 IF( info.NE.0 ) THEN
515 CALL xerbla( 'DDRVEV', -info )
516 RETURN
517 END IF
518*
519* Quick return if nothing to do
520*
521 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
522 $ RETURN
523*
524* More Important constants
525*
526 unfl = dlamch( 'Safe minimum' )
527 ovfl = one / unfl
528 ulp = dlamch( 'Precision' )
529 ulpinv = one / ulp
530 rtulp = sqrt( ulp )
531 rtulpi = one / rtulp
532*
533* Loop over sizes, types
534*
535 nerrs = 0
536*
537 DO 270 jsize = 1, nsizes
538 n = nn( jsize )
539 IF( nsizes.NE.1 ) THEN
540 mtypes = min( maxtyp, ntypes )
541 ELSE
542 mtypes = min( maxtyp+1, ntypes )
543 END IF
544*
545 DO 260 jtype = 1, mtypes
546 IF( .NOT.dotype( jtype ) )
547 $ GO TO 260
548*
549* Save ISEED in case of an error.
550*
551 DO 20 j = 1, 4
552 ioldsd( j ) = iseed( j )
553 20 CONTINUE
554*
555* Compute "A"
556*
557* Control parameters:
558*
559* KMAGN KCONDS KMODE KTYPE
560* =1 O(1) 1 clustered 1 zero
561* =2 large large clustered 2 identity
562* =3 small exponential Jordan
563* =4 arithmetic diagonal, (w/ eigenvalues)
564* =5 random log symmetric, w/ eigenvalues
565* =6 random general, w/ eigenvalues
566* =7 random diagonal
567* =8 random symmetric
568* =9 random general
569* =10 random triangular
570*
571 IF( mtypes.GT.maxtyp )
572 $ GO TO 90
573*
574 itype = ktype( jtype )
575 imode = kmode( jtype )
576*
577* Compute norm
578*
579 GO TO ( 30, 40, 50 )kmagn( jtype )
580*
581 30 CONTINUE
582 anorm = one
583 GO TO 60
584*
585 40 CONTINUE
586 anorm = ovfl*ulp
587 GO TO 60
588*
589 50 CONTINUE
590 anorm = unfl*ulpinv
591 GO TO 60
592*
593 60 CONTINUE
594*
595 CALL dlaset( 'Full', lda, n, zero, zero, a, lda )
596 iinfo = 0
597 cond = ulpinv
598*
599* Special Matrices -- Identity & Jordan block
600*
601* Zero
602*
603 IF( itype.EQ.1 ) THEN
604 iinfo = 0
605*
606 ELSE IF( itype.EQ.2 ) THEN
607*
608* Identity
609*
610 DO 70 jcol = 1, n
611 a( jcol, jcol ) = anorm
612 70 CONTINUE
613*
614 ELSE IF( itype.EQ.3 ) THEN
615*
616* Jordan Block
617*
618 DO 80 jcol = 1, n
619 a( jcol, jcol ) = anorm
620 IF( jcol.GT.1 )
621 $ a( jcol, jcol-1 ) = one
622 80 CONTINUE
623*
624 ELSE IF( itype.EQ.4 ) THEN
625*
626* Diagonal Matrix, [Eigen]values Specified
627*
628 CALL dlatms( n, n, 'S', iseed, 'S', work, imode, cond,
629 $ anorm, 0, 0, 'N', a, lda, work( n+1 ),
630 $ iinfo )
631*
632 ELSE IF( itype.EQ.5 ) THEN
633*
634* Symmetric, eigenvalues specified
635*
636 CALL dlatms( n, n, 'S', iseed, 'S', work, imode, cond,
637 $ anorm, n, n, 'N', a, lda, work( n+1 ),
638 $ iinfo )
639*
640 ELSE IF( itype.EQ.6 ) THEN
641*
642* General, eigenvalues specified
643*
644 IF( kconds( jtype ).EQ.1 ) THEN
645 conds = one
646 ELSE IF( kconds( jtype ).EQ.2 ) THEN
647 conds = rtulpi
648 ELSE
649 conds = zero
650 END IF
651*
652 adumma( 1 ) = ' '
653 CALL dlatme( n, 'S', iseed, work, imode, cond, one,
654 $ adumma, 'T', 'T', 'T', work( n+1 ), 4,
655 $ conds, n, n, anorm, a, lda, work( 2*n+1 ),
656 $ iinfo )
657*
658 ELSE IF( itype.EQ.7 ) THEN
659*
660* Diagonal, random eigenvalues
661*
662 CALL dlatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
663 $ 'T', 'N', work( n+1 ), 1, one,
664 $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
665 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
666*
667 ELSE IF( itype.EQ.8 ) THEN
668*
669* Symmetric, random eigenvalues
670*
671 CALL dlatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
672 $ 'T', 'N', work( n+1 ), 1, one,
673 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
674 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
675*
676 ELSE IF( itype.EQ.9 ) THEN
677*
678* General, random eigenvalues
679*
680 CALL dlatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
681 $ 'T', 'N', work( n+1 ), 1, one,
682 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
683 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
684 IF( n.GE.4 ) THEN
685 CALL dlaset( 'Full', 2, n, zero, zero, a, lda )
686 CALL dlaset( 'Full', n-3, 1, zero, zero, a( 3, 1 ),
687 $ lda )
688 CALL dlaset( 'Full', n-3, 2, zero, zero, a( 3, n-1 ),
689 $ lda )
690 CALL dlaset( 'Full', 1, n, zero, zero, a( n, 1 ),
691 $ lda )
692 END IF
693*
694 ELSE IF( itype.EQ.10 ) THEN
695*
696* Triangular, random eigenvalues
697*
698 CALL dlatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
699 $ 'T', 'N', work( n+1 ), 1, one,
700 $ work( 2*n+1 ), 1, one, 'N', idumma, n, 0,
701 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
702*
703 ELSE
704*
705 iinfo = 1
706 END IF
707*
708 IF( iinfo.NE.0 ) THEN
709 WRITE( nounit, fmt = 9993 )'Generator', iinfo, n, jtype,
710 $ ioldsd
711 info = abs( iinfo )
712 RETURN
713 END IF
714*
715 90 CONTINUE
716*
717* Test for minimal and generous workspace
718*
719 DO 250 iwk = 1, 2
720 IF( iwk.EQ.1 ) THEN
721 nnwork = 4*n
722 ELSE
723 nnwork = 5*n + 2*n**2
724 END IF
725 nnwork = max( nnwork, 1 )
726*
727* Initialize RESULT
728*
729 DO 100 j = 1, 7
730 result( j ) = -one
731 100 CONTINUE
732*
733* Compute eigenvalues and eigenvectors, and test them
734*
735 CALL dlacpy( 'F', n, n, a, lda, h, lda )
736 CALL dgeev( 'V', 'V', n, h, lda, wr, wi, vl, ldvl, vr,
737 $ ldvr, work, nnwork, iinfo )
738 IF( iinfo.NE.0 ) THEN
739 result( 1 ) = ulpinv
740 WRITE( nounit, fmt = 9993 )'DGEEV1', iinfo, n, jtype,
741 $ ioldsd
742 info = abs( iinfo )
743 GO TO 220
744 END IF
745*
746* Do Test (1)
747*
748 CALL dget22( 'N', 'N', 'N', n, a, lda, vr, ldvr, wr, wi,
749 $ work, res )
750 result( 1 ) = res( 1 )
751*
752* Do Test (2)
753*
754 CALL dget22( 'T', 'N', 'T', n, a, lda, vl, ldvl, wr, wi,
755 $ work, res )
756 result( 2 ) = res( 1 )
757*
758* Do Test (3)
759*
760 DO 120 j = 1, n
761 tnrm = one
762 IF( wi( j ).EQ.zero ) THEN
763 tnrm = dnrm2( n, vr( 1, j ), 1 )
764 ELSE IF( wi( j ).GT.zero ) THEN
765 tnrm = dlapy2( dnrm2( n, vr( 1, j ), 1 ),
766 $ dnrm2( n, vr( 1, j+1 ), 1 ) )
767 END IF
768 result( 3 ) = max( result( 3 ),
769 $ min( ulpinv, abs( tnrm-one ) / ulp ) )
770 IF( wi( j ).GT.zero ) THEN
771 vmx = zero
772 vrmx = zero
773 DO 110 jj = 1, n
774 vtst = dlapy2( vr( jj, j ), vr( jj, j+1 ) )
775 IF( vtst.GT.vmx )
776 $ vmx = vtst
777 IF( vr( jj, j+1 ).EQ.zero .AND.
778 $ abs( vr( jj, j ) ).GT.vrmx )
779 $ vrmx = abs( vr( jj, j ) )
780 110 CONTINUE
781 IF( vrmx / vmx.LT.one-two*ulp )
782 $ result( 3 ) = ulpinv
783 END IF
784 120 CONTINUE
785*
786* Do Test (4)
787*
788 DO 140 j = 1, n
789 tnrm = one
790 IF( wi( j ).EQ.zero ) THEN
791 tnrm = dnrm2( n, vl( 1, j ), 1 )
792 ELSE IF( wi( j ).GT.zero ) THEN
793 tnrm = dlapy2( dnrm2( n, vl( 1, j ), 1 ),
794 $ dnrm2( n, vl( 1, j+1 ), 1 ) )
795 END IF
796 result( 4 ) = max( result( 4 ),
797 $ min( ulpinv, abs( tnrm-one ) / ulp ) )
798 IF( wi( j ).GT.zero ) THEN
799 vmx = zero
800 vrmx = zero
801 DO 130 jj = 1, n
802 vtst = dlapy2( vl( jj, j ), vl( jj, j+1 ) )
803 IF( vtst.GT.vmx )
804 $ vmx = vtst
805 IF( vl( jj, j+1 ).EQ.zero .AND.
806 $ abs( vl( jj, j ) ).GT.vrmx )
807 $ vrmx = abs( vl( jj, j ) )
808 130 CONTINUE
809 IF( vrmx / vmx.LT.one-two*ulp )
810 $ result( 4 ) = ulpinv
811 END IF
812 140 CONTINUE
813*
814* Compute eigenvalues only, and test them
815*
816 CALL dlacpy( 'F', n, n, a, lda, h, lda )
817 CALL dgeev( 'N', 'N', n, h, lda, wr1, wi1, dum, 1, dum,
818 $ 1, work, nnwork, iinfo )
819 IF( iinfo.NE.0 ) THEN
820 result( 1 ) = ulpinv
821 WRITE( nounit, fmt = 9993 )'DGEEV2', iinfo, n, jtype,
822 $ ioldsd
823 info = abs( iinfo )
824 GO TO 220
825 END IF
826*
827* Do Test (5)
828*
829 DO 150 j = 1, n
830 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
831 $ result( 5 ) = ulpinv
832 150 CONTINUE
833*
834* Compute eigenvalues and right eigenvectors, and test them
835*
836 CALL dlacpy( 'F', n, n, a, lda, h, lda )
837 CALL dgeev( 'N', 'V', n, h, lda, wr1, wi1, dum, 1, lre,
838 $ ldlre, work, nnwork, iinfo )
839 IF( iinfo.NE.0 ) THEN
840 result( 1 ) = ulpinv
841 WRITE( nounit, fmt = 9993 )'DGEEV3', iinfo, n, jtype,
842 $ ioldsd
843 info = abs( iinfo )
844 GO TO 220
845 END IF
846*
847* Do Test (5) again
848*
849 DO 160 j = 1, n
850 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
851 $ result( 5 ) = ulpinv
852 160 CONTINUE
853*
854* Do Test (6)
855*
856 DO 180 j = 1, n
857 DO 170 jj = 1, n
858 IF( vr( j, jj ).NE.lre( j, jj ) )
859 $ result( 6 ) = ulpinv
860 170 CONTINUE
861 180 CONTINUE
862*
863* Compute eigenvalues and left eigenvectors, and test them
864*
865 CALL dlacpy( 'F', n, n, a, lda, h, lda )
866 CALL dgeev( 'V', 'N', n, h, lda, wr1, wi1, lre, ldlre,
867 $ dum, 1, work, nnwork, iinfo )
868 IF( iinfo.NE.0 ) THEN
869 result( 1 ) = ulpinv
870 WRITE( nounit, fmt = 9993 )'DGEEV4', iinfo, n, jtype,
871 $ ioldsd
872 info = abs( iinfo )
873 GO TO 220
874 END IF
875*
876* Do Test (5) again
877*
878 DO 190 j = 1, n
879 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
880 $ result( 5 ) = ulpinv
881 190 CONTINUE
882*
883* Do Test (7)
884*
885 DO 210 j = 1, n
886 DO 200 jj = 1, n
887 IF( vl( j, jj ).NE.lre( j, jj ) )
888 $ result( 7 ) = ulpinv
889 200 CONTINUE
890 210 CONTINUE
891*
892* End of Loop -- Check for RESULT(j) > THRESH
893*
894 220 CONTINUE
895*
896 ntest = 0
897 nfail = 0
898 DO 230 j = 1, 7
899 IF( result( j ).GE.zero )
900 $ ntest = ntest + 1
901 IF( result( j ).GE.thresh )
902 $ nfail = nfail + 1
903 230 CONTINUE
904*
905 IF( nfail.GT.0 )
906 $ ntestf = ntestf + 1
907 IF( ntestf.EQ.1 ) THEN
908 WRITE( nounit, fmt = 9999 )path
909 WRITE( nounit, fmt = 9998 )
910 WRITE( nounit, fmt = 9997 )
911 WRITE( nounit, fmt = 9996 )
912 WRITE( nounit, fmt = 9995 )thresh
913 ntestf = 2
914 END IF
915*
916 DO 240 j = 1, 7
917 IF( result( j ).GE.thresh ) THEN
918 WRITE( nounit, fmt = 9994 )n, iwk, ioldsd, jtype,
919 $ j, result( j )
920 END IF
921 240 CONTINUE
922*
923 nerrs = nerrs + nfail
924 ntestt = ntestt + ntest
925*
926 250 CONTINUE
927 260 CONTINUE
928 270 CONTINUE
929*
930* Summary
931*
932 CALL dlasum( path, nounit, nerrs, ntestt )
933*
934 9999 FORMAT( / 1x, a3, ' -- Real Eigenvalue-Eigenvector Decomposition',
935 $ ' Driver', / ' Matrix types (see DDRVEV for details): ' )
936*
937 9998 FORMAT( / ' Special Matrices:', / ' 1=Zero matrix. ',
938 $ ' ', ' 5=Diagonal: geometr. spaced entries.',
939 $ / ' 2=Identity matrix. ', ' 6=Diagona',
940 $ 'l: clustered entries.', / ' 3=Transposed Jordan block. ',
941 $ ' ', ' 7=Diagonal: large, evenly spaced.', / ' ',
942 $ '4=Diagonal: evenly spaced entries. ', ' 8=Diagonal: s',
943 $ 'mall, evenly spaced.' )
944 9997 FORMAT( ' Dense, Non-Symmetric Matrices:', / ' 9=Well-cond., ev',
945 $ 'enly spaced eigenvals.', ' 14=Ill-cond., geomet. spaced e',
946 $ 'igenals.', / ' 10=Well-cond., geom. spaced eigenvals. ',
947 $ ' 15=Ill-conditioned, clustered e.vals.', / ' 11=Well-cond',
948 $ 'itioned, clustered e.vals. ', ' 16=Ill-cond., random comp',
949 $ 'lex ', / ' 12=Well-cond., random complex ', 6x, ' ',
950 $ ' 17=Ill-cond., large rand. complx ', / ' 13=Ill-condi',
951 $ 'tioned, evenly spaced. ', ' 18=Ill-cond., small rand.',
952 $ ' complx ' )
953 9996 FORMAT( ' 19=Matrix with random O(1) entries. ', ' 21=Matrix ',
954 $ 'with small random entries.', / ' 20=Matrix with large ran',
955 $ 'dom entries. ', / )
956 9995 FORMAT( ' Tests performed with test threshold =', f8.2,
957 $ / / ' 1 = | A VR - VR W | / ( n |A| ulp ) ',
958 $ / ' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
959 $ / ' 3 = | |VR(i)| - 1 | / ulp ',
960 $ / ' 4 = | |VL(i)| - 1 | / ulp ',
961 $ / ' 5 = 0 if W same no matter if VR or VL computed,',
962 $ ' 1/ulp otherwise', /
963 $ ' 6 = 0 if VR same no matter if VL computed,',
964 $ ' 1/ulp otherwise', /
965 $ ' 7 = 0 if VL same no matter if VR computed,',
966 $ ' 1/ulp otherwise', / )
967 9994 FORMAT( ' N=', i5, ', IWK=', i2, ', seed=', 4( i4, ',' ),
968 $ ' type ', i2, ', test(', i2, ')=', g10.3 )
969 9993 FORMAT( ' DDRVEV: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
970 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
971*
972 RETURN
973*
974* End of DDRVEV
975*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dget22(transa, transe, transw, n, a, lda, e, lde, wr, wi, work, result)
DGET22
Definition dget22.f:168
subroutine dlasum(type, iounit, ie, nrun)
DLASUM
Definition dlasum.f:43
subroutine dlatme(n, dist, iseed, d, mode, cond, dmax, ei, rsign, upper, sim, ds, modes, conds, kl, ku, anorm, a, lda, work, info)
DLATME
Definition dlatme.f:332
subroutine dlatmr(m, n, dist, iseed, sym, d, mode, cond, dmax, rsign, grade, dl, model, condl, dr, moder, condr, pivtng, ipivot, kl, ku, sparse, anorm, pack, a, lda, iwork, info)
DLATMR
Definition dlatmr.f:471
subroutine dlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
DLATMS
Definition dlatms.f:321
subroutine dgeev(jobvl, jobvr, n, a, lda, wr, wi, vl, ldvl, vr, ldvr, work, lwork, info)
DGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition dgeev.f:192
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
double precision function dlapy2(x, y)
DLAPY2 returns sqrt(x2+y2).
Definition dlapy2.f:63
subroutine dlaset(uplo, m, n, alpha, beta, a, lda)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition dlaset.f:110
real(wp) function dnrm2(n, x, incx)
DNRM2
Definition dnrm2.f90:89
Here is the call graph for this function:
Here is the caller graph for this function: