LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
schkhs.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine schkhs (NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, NOUNIT, A, LDA, H, T1, T2, U, LDU, Z, UZ, WR1, WI1, WR3, WI3, EVECTL, EVECTR, EVECTY, EVECTX, UU, TAU, WORK, NWORK, IWORK, SELECT, RESULT, INFO)
 SCHKHS More...
 

Function/Subroutine Documentation

subroutine schkhs ( integer  NSIZES,
integer, dimension( * )  NN,
integer  NTYPES,
logical, dimension( * )  DOTYPE,
integer, dimension( 4 )  ISEED,
real  THRESH,
integer  NOUNIT,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( lda, * )  H,
real, dimension( lda, * )  T1,
real, dimension( lda, * )  T2,
real, dimension( ldu, * )  U,
integer  LDU,
real, dimension( ldu, * )  Z,
real, dimension( ldu, * )  UZ,
real, dimension( * )  WR1,
real, dimension( * )  WI1,
real, dimension( * )  WR3,
real, dimension( * )  WI3,
real, dimension( ldu, * )  EVECTL,
real, dimension( ldu, * )  EVECTR,
real, dimension( ldu, * )  EVECTY,
real, dimension( ldu, * )  EVECTX,
real, dimension( ldu, * )  UU,
real, dimension( * )  TAU,
real, dimension( * )  WORK,
integer  NWORK,
integer, dimension( * )  IWORK,
logical, dimension( * )  SELECT,
real, dimension( 14 )  RESULT,
integer  INFO 
)

SCHKHS

Purpose:
    SCHKHS  checks the nonsymmetric eigenvalue problem routines.

            SGEHRD factors A as  U H U' , where ' means transpose,
            H is hessenberg, and U is an orthogonal matrix.

            SORGHR generates the orthogonal matrix U.

            SORMHR multiplies a matrix by the orthogonal matrix U.

            SHSEQR factors H as  Z T Z' , where Z is orthogonal and
            T is "quasi-triangular", and the eigenvalue vector W.

            STREVC computes the left and right eigenvector matrices
            L and R for T.

            SHSEIN computes the left and right eigenvector matrices
            Y and X for H, using inverse iteration.

    When SCHKHS 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, 14
    tests will be performed:

    (1)     | A - U H U**T | / ( |A| n ulp )

    (2)     | I - UU**T | / ( n ulp )

    (3)     | H - Z T Z**T | / ( |H| n ulp )

    (4)     | I - ZZ**T | / ( n ulp )

    (5)     | A - UZ H (UZ)**T | / ( |A| n ulp )

    (6)     | I - UZ (UZ)**T | / ( n ulp )

    (7)     | T(Z computed) - T(Z not computed) | / ( |T| ulp )

    (8)     | W(Z computed) - W(Z not computed) | / ( |W| ulp )

    (9)     | TR - RW | / ( |T| |R| ulp )

    (10)    | L**H T - W**H L | / ( |T| |L| ulp )

    (11)    | HX - XW | / ( |H| |X| ulp )

    (12)    | Y**H H - W**H Y | / ( |H| |Y| ulp )

    (13)    | AX - XW | / ( |A| |X| ulp )

    (14)    | Y**H A - W**H Y | / ( |A| |Y| ulp )

    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 SQRT( overflow threshold )
    (8)  Same as (4), but multiplied by SQRT( 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 SQRT( overflow threshold )
    (18) Same as (16), but multiplied by SQRT( underflow threshold )

    (19) Nonsymmetric matrix with random entries chosen from (-1,1).
    (20) Same as (19), but multiplied by SQRT( overflow threshold )
    (21) Same as (19), but multiplied by SQRT( underflow threshold )
  NSIZES - INTEGER
           The number of sizes of matrices to use.  If it is zero,
           SCHKHS does nothing.  It must be at least zero.
           Not modified.

  NN     - 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.
           Not modified.

  NTYPES - INTEGER
           The number of elements in DOTYPE.   If it is zero, SCHKHS
           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. .
           Not modified.

  DOTYPE - 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.
           Not modified.

  ISEED  - 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 SCHKHS to continue the same random number
           sequence.
           Modified.

  THRESH - REAL
           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.
           Not modified.

  NOUNIT - INTEGER
           The FORTRAN unit number for printing out error messages
           (e.g., if a routine returns IINFO not equal to 0.)
           Not modified.

  A      - REAL 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.
           Modified.

  LDA    - INTEGER
           The leading dimension of A, H, T1 and T2.  It must be at
           least 1 and at least max( NN ).
           Not modified.

  H      - REAL array, dimension (LDA,max(NN))
           The upper hessenberg matrix computed by SGEHRD.  On exit,
           H contains the Hessenberg form of the matrix in A.
           Modified.

  T1     - REAL array, dimension (LDA,max(NN))
           The Schur (="quasi-triangular") matrix computed by SHSEQR
           if Z is computed.  On exit, T1 contains the Schur form of
           the matrix in A.
           Modified.

  T2     - REAL array, dimension (LDA,max(NN))
           The Schur matrix computed by SHSEQR when Z is not computed.
           This should be identical to T1.
           Modified.

  LDU    - INTEGER
           The leading dimension of U, Z, UZ and UU.  It must be at
           least 1 and at least max( NN ).
           Not modified.

  U      - REAL array, dimension (LDU,max(NN))
           The orthogonal matrix computed by SGEHRD.
           Modified.

  Z      - REAL array, dimension (LDU,max(NN))
           The orthogonal matrix computed by SHSEQR.
           Modified.

  UZ     - REAL array, dimension (LDU,max(NN))
           The product of U times Z.
           Modified.

  WR1    - REAL array, dimension (max(NN))
  WI1    - REAL array, dimension (max(NN))
           The real and imaginary parts of the eigenvalues of A,
           as computed when Z is computed.
           On exit, WR1 + WI1*i are the eigenvalues of the matrix in A.
           Modified.

  WR3    - REAL array, dimension (max(NN))
  WI3    - REAL array, dimension (max(NN))
           Like WR1, WI1, these arrays contain the eigenvalues of A,
           but those computed when SHSEQR only computes the
           eigenvalues, i.e., not the Schur vectors and no more of the
           Schur form than is necessary for computing the
           eigenvalues.
           Modified.

  EVECTL - REAL array, dimension (LDU,max(NN))
           The (upper triangular) left eigenvector matrix for the
           matrix in T1.  For complex conjugate pairs, the real part
           is stored in one row and the imaginary part in the next.
           Modified.

  EVECTR - REAL array, dimension (LDU,max(NN))
           The (upper triangular) right eigenvector matrix for the
           matrix in T1.  For complex conjugate pairs, the real part
           is stored in one column and the imaginary part in the next.
           Modified.

  EVECTY - REAL array, dimension (LDU,max(NN))
           The left eigenvector matrix for the
           matrix in H.  For complex conjugate pairs, the real part
           is stored in one row and the imaginary part in the next.
           Modified.

  EVECTX - REAL array, dimension (LDU,max(NN))
           The right eigenvector matrix for the
           matrix in H.  For complex conjugate pairs, the real part
           is stored in one column and the imaginary part in the next.
           Modified.

  UU     - REAL array, dimension (LDU,max(NN))
           Details of the orthogonal matrix computed by SGEHRD.
           Modified.

  TAU    - REAL array, dimension(max(NN))
           Further details of the orthogonal matrix computed by SGEHRD.
           Modified.

  WORK   - REAL array, dimension (NWORK)
           Workspace.
           Modified.

  NWORK  - INTEGER
           The number of entries in WORK.  NWORK >= 4*NN(j)*NN(j) + 2.

  IWORK  - INTEGER array, dimension (max(NN))
           Workspace.
           Modified.

  SELECT - LOGICAL array, dimension (max(NN))
           Workspace.
           Modified.

  RESULT - REAL array, dimension (14)
           The values computed by the fourteen tests described above.
           The values are currently limited to 1/ulp, to avoid
           overflow.
           Modified.

  INFO   - 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) ).
           -14: LDU < 1 or LDU < NMAX.
           -28: NWORK too small.
           If  SLATMR, SLATMS, or SLATME returns an error code, the
               absolute value of it is returned.
           If 1, then SHSEQR could not find all the shifts.
           If 2, then the EISPACK code (for small blocks) failed.
           If >2, then 30*N iterations were not enough to find an
               eigenvalue or to decompose the problem.
           Modified.

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

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

     ZERO, ONE       Real 0 and 1.
     MAXTYP          The number of types defined.
     MTEST           The number of tests defined: care must be taken
                     that (1) the size of RESULT, (2) the number of
                     tests actually performed, and (3) MTEST agree.
     NTEST           The number of tests performed on this matrix
                     so far.  This should be less than MTEST, and
                     equal to it by the last test.  It will be less
                     if any of the routines being tested indicates
                     that it could not compute the matrices that
                     would be tested.
     NMAX            Largest value in NN.
     NMATS           The number of matrices generated so far.
     NERRS           The number of tests which have exceeded THRESH
                     so far (computed by SLAFTS).
     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.
     RTOVFL, RTUNFL,
     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)       Selects 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.
Date
November 2011

Definition at line 407 of file schkhs.f.

407 *
408 * -- LAPACK test routine (version 3.4.0) --
409 * -- LAPACK is a software package provided by Univ. of Tennessee, --
410 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
411 * November 2011
412 *
413 * .. Scalar Arguments ..
414  INTEGER info, lda, ldu, nounit, nsizes, ntypes, nwork
415  REAL thresh
416 * ..
417 * .. Array Arguments ..
418  LOGICAL dotype( * ), select( * )
419  INTEGER iseed( 4 ), iwork( * ), nn( * )
420  REAL a( lda, * ), evectl( ldu, * ),
421  $ evectr( ldu, * ), evectx( ldu, * ),
422  $ evecty( ldu, * ), h( lda, * ), result( 14 ),
423  $ t1( lda, * ), t2( lda, * ), tau( * ),
424  $ u( ldu, * ), uu( ldu, * ), uz( ldu, * ),
425  $ wi1( * ), wi3( * ), work( * ), wr1( * ),
426  $ wr3( * ), z( ldu, * )
427 * ..
428 *
429 * =====================================================================
430 *
431 * .. Parameters ..
432  REAL zero, one
433  parameter( zero = 0.0, one = 1.0 )
434  INTEGER maxtyp
435  parameter( maxtyp = 21 )
436 * ..
437 * .. Local Scalars ..
438  LOGICAL badnn, match
439  INTEGER i, ihi, iinfo, ilo, imode, in, itype, j, jcol,
440  $ jj, jsize, jtype, k, mtypes, n, n1, nerrs,
441  $ nmats, nmax, nselc, nselr, ntest, ntestt
442  REAL aninv, anorm, cond, conds, ovfl, rtovfl, rtulp,
443  $ rtulpi, rtunfl, temp1, temp2, ulp, ulpinv, unfl
444 * ..
445 * .. Local Arrays ..
446  CHARACTER adumma( 1 )
447  INTEGER idumma( 1 ), ioldsd( 4 ), kconds( maxtyp ),
448  $ kmagn( maxtyp ), kmode( maxtyp ),
449  $ ktype( maxtyp )
450  REAL dumma( 6 )
451 * ..
452 * .. External Functions ..
453  REAL slamch
454  EXTERNAL slamch
455 * ..
456 * .. External Subroutines ..
457  EXTERNAL scopy, sgehrd, sgemm, sget10, sget22, shsein,
460  $ strevc, xerbla
461 * ..
462 * .. Intrinsic Functions ..
463  INTRINSIC abs, max, min, REAL, sqrt
464 * ..
465 * .. Data statements ..
466  DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
467  DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
468  $ 3, 1, 2, 3 /
469  DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
470  $ 1, 5, 5, 5, 4, 3, 1 /
471  DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
472 * ..
473 * .. Executable Statements ..
474 *
475 * Check for errors
476 *
477  ntestt = 0
478  info = 0
479 *
480  badnn = .false.
481  nmax = 0
482  DO 10 j = 1, nsizes
483  nmax = max( nmax, nn( j ) )
484  IF( nn( j ).LT.0 )
485  $ badnn = .true.
486  10 CONTINUE
487 *
488 * Check for errors
489 *
490  IF( nsizes.LT.0 ) THEN
491  info = -1
492  ELSE IF( badnn ) THEN
493  info = -2
494  ELSE IF( ntypes.LT.0 ) THEN
495  info = -3
496  ELSE IF( thresh.LT.zero ) THEN
497  info = -6
498  ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
499  info = -9
500  ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax ) THEN
501  info = -14
502  ELSE IF( 4*nmax*nmax+2.GT.nwork ) THEN
503  info = -28
504  END IF
505 *
506  IF( info.NE.0 ) THEN
507  CALL xerbla( 'SCHKHS', -info )
508  RETURN
509  END IF
510 *
511 * Quick return if possible
512 *
513  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
514  $ RETURN
515 *
516 * More important constants
517 *
518  unfl = slamch( 'Safe minimum' )
519  ovfl = slamch( 'Overflow' )
520  CALL slabad( unfl, ovfl )
521  ulp = slamch( 'Epsilon' )*slamch( 'Base' )
522  ulpinv = one / ulp
523  rtunfl = sqrt( unfl )
524  rtovfl = sqrt( ovfl )
525  rtulp = sqrt( ulp )
526  rtulpi = one / rtulp
527 *
528 * Loop over sizes, types
529 *
530  nerrs = 0
531  nmats = 0
532 *
533  DO 270 jsize = 1, nsizes
534  n = nn( jsize )
535  IF( n.EQ.0 )
536  $ GO TO 270
537  n1 = max( 1, n )
538  aninv = one / REAL( n1 )
539 *
540  IF( nsizes.NE.1 ) THEN
541  mtypes = min( maxtyp, ntypes )
542  ELSE
543  mtypes = min( maxtyp+1, ntypes )
544  END IF
545 *
546  DO 260 jtype = 1, mtypes
547  IF( .NOT.dotype( jtype ) )
548  $ GO TO 260
549  nmats = nmats + 1
550  ntest = 0
551 *
552 * Save ISEED in case of an error.
553 *
554  DO 20 j = 1, 4
555  ioldsd( j ) = iseed( j )
556  20 CONTINUE
557 *
558 * Initialize RESULT
559 *
560  DO 30 j = 1, 14
561  result( j ) = zero
562  30 CONTINUE
563 *
564 * Compute "A"
565 *
566 * Control parameters:
567 *
568 * KMAGN KCONDS KMODE KTYPE
569 * =1 O(1) 1 clustered 1 zero
570 * =2 large large clustered 2 identity
571 * =3 small exponential Jordan
572 * =4 arithmetic diagonal, (w/ eigenvalues)
573 * =5 random log symmetric, w/ eigenvalues
574 * =6 random general, w/ eigenvalues
575 * =7 random diagonal
576 * =8 random symmetric
577 * =9 random general
578 * =10 random triangular
579 *
580  IF( mtypes.GT.maxtyp )
581  $ GO TO 100
582 *
583  itype = ktype( jtype )
584  imode = kmode( jtype )
585 *
586 * Compute norm
587 *
588  GO TO ( 40, 50, 60 )kmagn( jtype )
589 *
590  40 CONTINUE
591  anorm = one
592  GO TO 70
593 *
594  50 CONTINUE
595  anorm = ( rtovfl*ulp )*aninv
596  GO TO 70
597 *
598  60 CONTINUE
599  anorm = rtunfl*n*ulpinv
600  GO TO 70
601 *
602  70 CONTINUE
603 *
604  CALL slaset( 'Full', lda, n, zero, zero, a, lda )
605  iinfo = 0
606  cond = ulpinv
607 *
608 * Special Matrices
609 *
610  IF( itype.EQ.1 ) THEN
611 *
612 * Zero
613 *
614  iinfo = 0
615 *
616  ELSE IF( itype.EQ.2 ) THEN
617 *
618 * Identity
619 *
620  DO 80 jcol = 1, n
621  a( jcol, jcol ) = anorm
622  80 CONTINUE
623 *
624  ELSE IF( itype.EQ.3 ) THEN
625 *
626 * Jordan Block
627 *
628  DO 90 jcol = 1, n
629  a( jcol, jcol ) = anorm
630  IF( jcol.GT.1 )
631  $ a( jcol, jcol-1 ) = one
632  90 CONTINUE
633 *
634  ELSE IF( itype.EQ.4 ) THEN
635 *
636 * Diagonal Matrix, [Eigen]values Specified
637 *
638  CALL slatms( n, n, 'S', iseed, 'S', work, imode, cond,
639  $ anorm, 0, 0, 'N', a, lda, work( n+1 ),
640  $ iinfo )
641 *
642  ELSE IF( itype.EQ.5 ) THEN
643 *
644 * Symmetric, eigenvalues specified
645 *
646  CALL slatms( n, n, 'S', iseed, 'S', work, imode, cond,
647  $ anorm, n, n, 'N', a, lda, work( n+1 ),
648  $ iinfo )
649 *
650  ELSE IF( itype.EQ.6 ) THEN
651 *
652 * General, eigenvalues specified
653 *
654  IF( kconds( jtype ).EQ.1 ) THEN
655  conds = one
656  ELSE IF( kconds( jtype ).EQ.2 ) THEN
657  conds = rtulpi
658  ELSE
659  conds = zero
660  END IF
661 *
662  adumma( 1 ) = ' '
663  CALL slatme( n, 'S', iseed, work, imode, cond, one,
664  $ adumma, 'T', 'T', 'T', work( n+1 ), 4,
665  $ conds, n, n, anorm, a, lda, work( 2*n+1 ),
666  $ iinfo )
667 *
668  ELSE IF( itype.EQ.7 ) THEN
669 *
670 * Diagonal, random eigenvalues
671 *
672  CALL slatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
673  $ 'T', 'N', work( n+1 ), 1, one,
674  $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
675  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
676 *
677  ELSE IF( itype.EQ.8 ) THEN
678 *
679 * Symmetric, random eigenvalues
680 *
681  CALL slatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
682  $ 'T', 'N', work( n+1 ), 1, one,
683  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
684  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
685 *
686  ELSE IF( itype.EQ.9 ) THEN
687 *
688 * General, random eigenvalues
689 *
690  CALL slatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
691  $ 'T', 'N', work( n+1 ), 1, one,
692  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
693  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
694 *
695  ELSE IF( itype.EQ.10 ) THEN
696 *
697 * Triangular, random eigenvalues
698 *
699  CALL slatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
700  $ 'T', 'N', work( n+1 ), 1, one,
701  $ work( 2*n+1 ), 1, one, 'N', idumma, n, 0,
702  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
703 *
704  ELSE
705 *
706  iinfo = 1
707  END IF
708 *
709  IF( iinfo.NE.0 ) THEN
710  WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
711  $ ioldsd
712  info = abs( iinfo )
713  RETURN
714  END IF
715 *
716  100 CONTINUE
717 *
718 * Call SGEHRD to compute H and U, do tests.
719 *
720  CALL slacpy( ' ', n, n, a, lda, h, lda )
721 *
722  ntest = 1
723 *
724  ilo = 1
725  ihi = n
726 *
727  CALL sgehrd( n, ilo, ihi, h, lda, work, work( n+1 ),
728  $ nwork-n, iinfo )
729 *
730  IF( iinfo.NE.0 ) THEN
731  result( 1 ) = ulpinv
732  WRITE( nounit, fmt = 9999 )'SGEHRD', iinfo, n, jtype,
733  $ ioldsd
734  info = abs( iinfo )
735  GO TO 250
736  END IF
737 *
738  DO 120 j = 1, n - 1
739  uu( j+1, j ) = zero
740  DO 110 i = j + 2, n
741  u( i, j ) = h( i, j )
742  uu( i, j ) = h( i, j )
743  h( i, j ) = zero
744  110 CONTINUE
745  120 CONTINUE
746  CALL scopy( n-1, work, 1, tau, 1 )
747  CALL sorghr( n, ilo, ihi, u, ldu, work, work( n+1 ),
748  $ nwork-n, iinfo )
749  ntest = 2
750 *
751  CALL shst01( n, ilo, ihi, a, lda, h, lda, u, ldu, work,
752  $ nwork, result( 1 ) )
753 *
754 * Call SHSEQR to compute T1, T2 and Z, do tests.
755 *
756 * Eigenvalues only (WR3,WI3)
757 *
758  CALL slacpy( ' ', n, n, h, lda, t2, lda )
759  ntest = 3
760  result( 3 ) = ulpinv
761 *
762  CALL shseqr( 'E', 'N', n, ilo, ihi, t2, lda, wr3, wi3, uz,
763  $ ldu, work, nwork, iinfo )
764  IF( iinfo.NE.0 ) THEN
765  WRITE( nounit, fmt = 9999 )'SHSEQR(E)', iinfo, n, jtype,
766  $ ioldsd
767  IF( iinfo.LE.n+2 ) THEN
768  info = abs( iinfo )
769  GO TO 250
770  END IF
771  END IF
772 *
773 * Eigenvalues (WR1,WI1) and Full Schur Form (T2)
774 *
775  CALL slacpy( ' ', n, n, h, lda, t2, lda )
776 *
777  CALL shseqr( 'S', 'N', n, ilo, ihi, t2, lda, wr1, wi1, uz,
778  $ ldu, work, nwork, iinfo )
779  IF( iinfo.NE.0 .AND. iinfo.LE.n+2 ) THEN
780  WRITE( nounit, fmt = 9999 )'SHSEQR(S)', iinfo, n, jtype,
781  $ ioldsd
782  info = abs( iinfo )
783  GO TO 250
784  END IF
785 *
786 * Eigenvalues (WR1,WI1), Schur Form (T1), and Schur vectors
787 * (UZ)
788 *
789  CALL slacpy( ' ', n, n, h, lda, t1, lda )
790  CALL slacpy( ' ', n, n, u, ldu, uz, lda )
791 *
792  CALL shseqr( 'S', 'V', n, ilo, ihi, t1, lda, wr1, wi1, uz,
793  $ ldu, work, nwork, iinfo )
794  IF( iinfo.NE.0 .AND. iinfo.LE.n+2 ) THEN
795  WRITE( nounit, fmt = 9999 )'SHSEQR(V)', iinfo, n, jtype,
796  $ ioldsd
797  info = abs( iinfo )
798  GO TO 250
799  END IF
800 *
801 * Compute Z = U' UZ
802 *
803  CALL sgemm( 'T', 'N', n, n, n, one, u, ldu, uz, ldu, zero,
804  $ z, ldu )
805  ntest = 8
806 *
807 * Do Tests 3: | H - Z T Z' | / ( |H| n ulp )
808 * and 4: | I - Z Z' | / ( n ulp )
809 *
810  CALL shst01( n, ilo, ihi, h, lda, t1, lda, z, ldu, work,
811  $ nwork, result( 3 ) )
812 *
813 * Do Tests 5: | A - UZ T (UZ)' | / ( |A| n ulp )
814 * and 6: | I - UZ (UZ)' | / ( n ulp )
815 *
816  CALL shst01( n, ilo, ihi, a, lda, t1, lda, uz, ldu, work,
817  $ nwork, result( 5 ) )
818 *
819 * Do Test 7: | T2 - T1 | / ( |T| n ulp )
820 *
821  CALL sget10( n, n, t2, lda, t1, lda, work, result( 7 ) )
822 *
823 * Do Test 8: | W3 - W1 | / ( max(|W1|,|W3|) ulp )
824 *
825  temp1 = zero
826  temp2 = zero
827  DO 130 j = 1, n
828  temp1 = max( temp1, abs( wr1( j ) )+abs( wi1( j ) ),
829  $ abs( wr3( j ) )+abs( wi3( j ) ) )
830  temp2 = max( temp2, abs( wr1( j )-wr3( j ) )+
831  $ abs( wr1( j )-wr3( j ) ) )
832  130 CONTINUE
833 *
834  result( 8 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
835 *
836 * Compute the Left and Right Eigenvectors of T
837 *
838 * Compute the Right eigenvector Matrix:
839 *
840  ntest = 9
841  result( 9 ) = ulpinv
842 *
843 * Select last max(N/4,1) real, max(N/4,1) complex eigenvectors
844 *
845  nselc = 0
846  nselr = 0
847  j = n
848  140 CONTINUE
849  IF( wi1( j ).EQ.zero ) THEN
850  IF( nselr.LT.max( n / 4, 1 ) ) THEN
851  nselr = nselr + 1
852  SELECT( j ) = .true.
853  ELSE
854  SELECT( j ) = .false.
855  END IF
856  j = j - 1
857  ELSE
858  IF( nselc.LT.max( n / 4, 1 ) ) THEN
859  nselc = nselc + 1
860  SELECT( j ) = .true.
861  SELECT( j-1 ) = .false.
862  ELSE
863  SELECT( j ) = .false.
864  SELECT( j-1 ) = .false.
865  END IF
866  j = j - 2
867  END IF
868  IF( j.GT.0 )
869  $ GO TO 140
870 *
871  CALL strevc( 'Right', 'All', SELECT, n, t1, lda, dumma, ldu,
872  $ evectr, ldu, n, in, work, iinfo )
873  IF( iinfo.NE.0 ) THEN
874  WRITE( nounit, fmt = 9999 )'STREVC(R,A)', iinfo, n,
875  $ jtype, ioldsd
876  info = abs( iinfo )
877  GO TO 250
878  END IF
879 *
880 * Test 9: | TR - RW | / ( |T| |R| ulp )
881 *
882  CALL sget22( 'N', 'N', 'N', n, t1, lda, evectr, ldu, wr1,
883  $ wi1, work, dumma( 1 ) )
884  result( 9 ) = dumma( 1 )
885  IF( dumma( 2 ).GT.thresh ) THEN
886  WRITE( nounit, fmt = 9998 )'Right', 'STREVC',
887  $ dumma( 2 ), n, jtype, ioldsd
888  END IF
889 *
890 * Compute selected right eigenvectors and confirm that
891 * they agree with previous right eigenvectors
892 *
893  CALL strevc( 'Right', 'Some', SELECT, n, t1, lda, dumma,
894  $ ldu, evectl, ldu, n, in, work, iinfo )
895  IF( iinfo.NE.0 ) THEN
896  WRITE( nounit, fmt = 9999 )'STREVC(R,S)', iinfo, n,
897  $ jtype, ioldsd
898  info = abs( iinfo )
899  GO TO 250
900  END IF
901 *
902  k = 1
903  match = .true.
904  DO 170 j = 1, n
905  IF( SELECT( j ) .AND. wi1( j ).EQ.zero ) THEN
906  DO 150 jj = 1, n
907  IF( evectr( jj, j ).NE.evectl( jj, k ) ) THEN
908  match = .false.
909  GO TO 180
910  END IF
911  150 CONTINUE
912  k = k + 1
913  ELSE IF( SELECT( j ) .AND. wi1( j ).NE.zero ) THEN
914  DO 160 jj = 1, n
915  IF( evectr( jj, j ).NE.evectl( jj, k ) .OR.
916  $ evectr( jj, j+1 ).NE.evectl( jj, k+1 ) ) THEN
917  match = .false.
918  GO TO 180
919  END IF
920  160 CONTINUE
921  k = k + 2
922  END IF
923  170 CONTINUE
924  180 CONTINUE
925  IF( .NOT.match )
926  $ WRITE( nounit, fmt = 9997 )'Right', 'STREVC', n, jtype,
927  $ ioldsd
928 *
929 * Compute the Left eigenvector Matrix:
930 *
931  ntest = 10
932  result( 10 ) = ulpinv
933  CALL strevc( 'Left', 'All', SELECT, n, t1, lda, evectl, ldu,
934  $ dumma, ldu, n, in, work, iinfo )
935  IF( iinfo.NE.0 ) THEN
936  WRITE( nounit, fmt = 9999 )'STREVC(L,A)', iinfo, n,
937  $ jtype, ioldsd
938  info = abs( iinfo )
939  GO TO 250
940  END IF
941 *
942 * Test 10: | LT - WL | / ( |T| |L| ulp )
943 *
944  CALL sget22( 'Trans', 'N', 'Conj', n, t1, lda, evectl, ldu,
945  $ wr1, wi1, work, dumma( 3 ) )
946  result( 10 ) = dumma( 3 )
947  IF( dumma( 4 ).GT.thresh ) THEN
948  WRITE( nounit, fmt = 9998 )'Left', 'STREVC', dumma( 4 ),
949  $ n, jtype, ioldsd
950  END IF
951 *
952 * Compute selected left eigenvectors and confirm that
953 * they agree with previous left eigenvectors
954 *
955  CALL strevc( 'Left', 'Some', SELECT, n, t1, lda, evectr,
956  $ ldu, dumma, ldu, n, in, work, iinfo )
957  IF( iinfo.NE.0 ) THEN
958  WRITE( nounit, fmt = 9999 )'STREVC(L,S)', iinfo, n,
959  $ jtype, ioldsd
960  info = abs( iinfo )
961  GO TO 250
962  END IF
963 *
964  k = 1
965  match = .true.
966  DO 210 j = 1, n
967  IF( SELECT( j ) .AND. wi1( j ).EQ.zero ) THEN
968  DO 190 jj = 1, n
969  IF( evectl( jj, j ).NE.evectr( jj, k ) ) THEN
970  match = .false.
971  GO TO 220
972  END IF
973  190 CONTINUE
974  k = k + 1
975  ELSE IF( SELECT( j ) .AND. wi1( j ).NE.zero ) THEN
976  DO 200 jj = 1, n
977  IF( evectl( jj, j ).NE.evectr( jj, k ) .OR.
978  $ evectl( jj, j+1 ).NE.evectr( jj, k+1 ) ) THEN
979  match = .false.
980  GO TO 220
981  END IF
982  200 CONTINUE
983  k = k + 2
984  END IF
985  210 CONTINUE
986  220 CONTINUE
987  IF( .NOT.match )
988  $ WRITE( nounit, fmt = 9997 )'Left', 'STREVC', n, jtype,
989  $ ioldsd
990 *
991 * Call SHSEIN for Right eigenvectors of H, do test 11
992 *
993  ntest = 11
994  result( 11 ) = ulpinv
995  DO 230 j = 1, n
996  SELECT( j ) = .true.
997  230 CONTINUE
998 *
999  CALL shsein( 'Right', 'Qr', 'Ninitv', SELECT, n, h, lda,
1000  $ wr3, wi3, dumma, ldu, evectx, ldu, n1, in,
1001  $ work, iwork, iwork, iinfo )
1002  IF( iinfo.NE.0 ) THEN
1003  WRITE( nounit, fmt = 9999 )'SHSEIN(R)', iinfo, n, jtype,
1004  $ ioldsd
1005  info = abs( iinfo )
1006  IF( iinfo.LT.0 )
1007  $ GO TO 250
1008  ELSE
1009 *
1010 * Test 11: | HX - XW | / ( |H| |X| ulp )
1011 *
1012 * (from inverse iteration)
1013 *
1014  CALL sget22( 'N', 'N', 'N', n, h, lda, evectx, ldu, wr3,
1015  $ wi3, work, dumma( 1 ) )
1016  IF( dumma( 1 ).LT.ulpinv )
1017  $ result( 11 ) = dumma( 1 )*aninv
1018  IF( dumma( 2 ).GT.thresh ) THEN
1019  WRITE( nounit, fmt = 9998 )'Right', 'SHSEIN',
1020  $ dumma( 2 ), n, jtype, ioldsd
1021  END IF
1022  END IF
1023 *
1024 * Call SHSEIN for Left eigenvectors of H, do test 12
1025 *
1026  ntest = 12
1027  result( 12 ) = ulpinv
1028  DO 240 j = 1, n
1029  SELECT( j ) = .true.
1030  240 CONTINUE
1031 *
1032  CALL shsein( 'Left', 'Qr', 'Ninitv', SELECT, n, h, lda, wr3,
1033  $ wi3, evecty, ldu, dumma, ldu, n1, in, work,
1034  $ iwork, iwork, iinfo )
1035  IF( iinfo.NE.0 ) THEN
1036  WRITE( nounit, fmt = 9999 )'SHSEIN(L)', iinfo, n, jtype,
1037  $ ioldsd
1038  info = abs( iinfo )
1039  IF( iinfo.LT.0 )
1040  $ GO TO 250
1041  ELSE
1042 *
1043 * Test 12: | YH - WY | / ( |H| |Y| ulp )
1044 *
1045 * (from inverse iteration)
1046 *
1047  CALL sget22( 'C', 'N', 'C', n, h, lda, evecty, ldu, wr3,
1048  $ wi3, work, dumma( 3 ) )
1049  IF( dumma( 3 ).LT.ulpinv )
1050  $ result( 12 ) = dumma( 3 )*aninv
1051  IF( dumma( 4 ).GT.thresh ) THEN
1052  WRITE( nounit, fmt = 9998 )'Left', 'SHSEIN',
1053  $ dumma( 4 ), n, jtype, ioldsd
1054  END IF
1055  END IF
1056 *
1057 * Call SORMHR for Right eigenvectors of A, do test 13
1058 *
1059  ntest = 13
1060  result( 13 ) = ulpinv
1061 *
1062  CALL sormhr( 'Left', 'No transpose', n, n, ilo, ihi, uu,
1063  $ ldu, tau, evectx, ldu, work, nwork, iinfo )
1064  IF( iinfo.NE.0 ) THEN
1065  WRITE( nounit, fmt = 9999 )'SORMHR(R)', iinfo, n, jtype,
1066  $ ioldsd
1067  info = abs( iinfo )
1068  IF( iinfo.LT.0 )
1069  $ GO TO 250
1070  ELSE
1071 *
1072 * Test 13: | AX - XW | / ( |A| |X| ulp )
1073 *
1074 * (from inverse iteration)
1075 *
1076  CALL sget22( 'N', 'N', 'N', n, a, lda, evectx, ldu, wr3,
1077  $ wi3, work, dumma( 1 ) )
1078  IF( dumma( 1 ).LT.ulpinv )
1079  $ result( 13 ) = dumma( 1 )*aninv
1080  END IF
1081 *
1082 * Call SORMHR for Left eigenvectors of A, do test 14
1083 *
1084  ntest = 14
1085  result( 14 ) = ulpinv
1086 *
1087  CALL sormhr( 'Left', 'No transpose', n, n, ilo, ihi, uu,
1088  $ ldu, tau, evecty, ldu, work, nwork, iinfo )
1089  IF( iinfo.NE.0 ) THEN
1090  WRITE( nounit, fmt = 9999 )'SORMHR(L)', iinfo, n, jtype,
1091  $ ioldsd
1092  info = abs( iinfo )
1093  IF( iinfo.LT.0 )
1094  $ GO TO 250
1095  ELSE
1096 *
1097 * Test 14: | YA - WY | / ( |A| |Y| ulp )
1098 *
1099 * (from inverse iteration)
1100 *
1101  CALL sget22( 'C', 'N', 'C', n, a, lda, evecty, ldu, wr3,
1102  $ wi3, work, dumma( 3 ) )
1103  IF( dumma( 3 ).LT.ulpinv )
1104  $ result( 14 ) = dumma( 3 )*aninv
1105  END IF
1106 *
1107 * End of Loop -- Check for RESULT(j) > THRESH
1108 *
1109  250 CONTINUE
1110 *
1111  ntestt = ntestt + ntest
1112  CALL slafts( 'SHS', n, n, jtype, ntest, result, ioldsd,
1113  $ thresh, nounit, nerrs )
1114 *
1115  260 CONTINUE
1116  270 CONTINUE
1117 *
1118 * Summary
1119 *
1120  CALL slasum( 'SHS', nounit, nerrs, ntestt )
1121 *
1122  RETURN
1123 *
1124  9999 FORMAT( ' SCHKHS: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1125  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1126  9998 FORMAT( ' SCHKHS: ', a, ' Eigenvectors from ', a, ' incorrectly ',
1127  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
1128  $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
1129  $ ')' )
1130  9997 FORMAT( ' SCHKHS: Selected ', a, ' Eigenvectors from ', a,
1131  $ ' do not match other eigenvectors ', 9x, 'N=', i6,
1132  $ ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1133 *
1134 * End of SCHKHS
1135 *
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: slaset.f:112
subroutine slatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
SLATMS
Definition: slatms.f:323
subroutine slatmr(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)
SLATMR
Definition: slatmr.f:473
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine shst01(N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK, LWORK, RESULT)
SHST01
Definition: shst01.f:136
subroutine shseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, WR, WI, Z, LDZ, WORK, LWORK, INFO)
SHSEQR
Definition: shseqr.f:318
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53
subroutine sget22(TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, WR, WI, WORK, RESULT)
SGET22
Definition: sget22.f:169
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
subroutine sormhr(SIDE, TRANS, M, N, ILO, IHI, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMHR
Definition: sormhr.f:181
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
Definition: sgemm.f:189
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:76
subroutine slasum(TYPE, IOUNIT, IE, NRUN)
SLASUM
Definition: slasum.f:42
subroutine sorghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
SORGHR
Definition: sorghr.f:128
subroutine slatme(N, DIST, ISEED, D, MODE, COND, DMAX, EI, RSIGN, UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM, A, LDA, WORK, INFO)
SLATME
Definition: slatme.f:334
subroutine strevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
STREVC
Definition: strevc.f:224
subroutine sgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
SGEHRD
Definition: sgehrd.f:170
subroutine shsein(SIDE, EIGSRC, INITV, SELECT, N, H, LDH, WR, WI, VL, LDVL, VR, LDVR, MM, M, WORK, IFAILL, IFAILR, INFO)
SHSEIN
Definition: shsein.f:265
subroutine sget10(M, N, A, LDA, B, LDB, WORK, RESULT)
SGET10
Definition: sget10.f:95
subroutine slafts(TYPE, M, N, IMAT, NTESTS, RESULT, ISEED, THRESH, IOUNIT, IE)
SLAFTS
Definition: slafts.f:101

Here is the call graph for this function:

Here is the caller graph for this function: