LAPACK  3.6.1 LAPACK: Linear Algebra PACKage
 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( * ) WR2, real, dimension( * ) WI2, 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.

WR2    - REAL array, dimension (max(NN))
WI2    - REAL array, dimension (max(NN))
The real and imaginary parts of the eigenvalues of A,
as computed when T is computed but not Z.
On exit, WR2 + WI2*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.)```
Date
November 2015

Definition at line 414 of file schkhs.f.

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

Here is the call graph for this function:

Here is the caller graph for this function: