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

◆ cchkst2stg()

subroutine cchkst2stg ( integer  NSIZES,
integer, dimension( * )  NN,
integer  NTYPES,
logical, dimension( * )  DOTYPE,
integer, dimension( 4 )  ISEED,
real  THRESH,
integer  NOUNIT,
complex, dimension( lda, * )  A,
integer  LDA,
complex, dimension( * )  AP,
real, dimension( * )  SD,
real, dimension( * )  SE,
real, dimension( * )  D1,
real, dimension( * )  D2,
real, dimension( * )  D3,
real, dimension( * )  D4,
real, dimension( * )  D5,
real, dimension( * )  WA1,
real, dimension( * )  WA2,
real, dimension( * )  WA3,
real, dimension( * )  WR,
complex, dimension( ldu, * )  U,
integer  LDU,
complex, dimension( ldu, * )  V,
complex, dimension( * )  VP,
complex, dimension( * )  TAU,
complex, dimension( ldu, * )  Z,
complex, dimension( * )  WORK,
integer  LWORK,
real, dimension( * )  RWORK,
integer  LRWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
real, dimension( * )  RESULT,
integer  INFO 
)

CCHKST2STG

Purpose:
 CCHKST2STG  checks the Hermitian eigenvalue problem routines
 using the 2-stage reduction techniques. Since the generation
 of Q or the vectors is not available in this release, we only 
 compare the eigenvalue resulting when using the 2-stage to the 
 one considered as reference using the standard 1-stage reduction
 CHETRD. For that, we call the standard CHETRD and compute D1 using 
 DSTEQR, then we call the 2-stage CHETRD_2STAGE with Upper and Lower
 and we compute D2 and D3 using DSTEQR and then we replaced tests
 3 and 4 by tests 11 and 12. test 1 and 2 remain to verify that 
 the 1-stage results are OK and can be trusted.
 This testing routine will converge to the CCHKST in the next 
 release when vectors and generation of Q will be implemented.

    CHETRD factors A as  U S U* , where * means conjugate transpose,
    S is real symmetric tridiagonal, and U is unitary.
    CHETRD can use either just the lower or just the upper triangle
    of A; CCHKST2STG checks both cases.
    U is represented as a product of Householder
    transformations, whose vectors are stored in the first
    n-1 columns of V, and whose scale factors are in TAU.

    CHPTRD does the same as CHETRD, except that A and V are stored
    in "packed" format.

    CUNGTR constructs the matrix U from the contents of V and TAU.

    CUPGTR constructs the matrix U from the contents of VP and TAU.

    CSTEQR factors S as  Z D1 Z* , where Z is the unitary
    matrix of eigenvectors and D1 is a diagonal matrix with
    the eigenvalues on the diagonal.  D2 is the matrix of
    eigenvalues computed when Z is not computed.

    SSTERF computes D3, the matrix of eigenvalues, by the
    PWK method, which does not yield eigenvectors.

    CPTEQR factors S as  Z4 D4 Z4* , for a
    Hermitian positive definite tridiagonal matrix.
    D5 is the matrix of eigenvalues computed when Z is not
    computed.

    SSTEBZ computes selected eigenvalues.  WA1, WA2, and
    WA3 will denote eigenvalues computed to high
    absolute accuracy, with different range options.
    WR will denote eigenvalues computed to high relative
    accuracy.

    CSTEIN computes Y, the eigenvectors of S, given the
    eigenvalues.

    CSTEDC factors S as Z D1 Z* , where Z is the unitary
    matrix of eigenvectors and D1 is a diagonal matrix with
    the eigenvalues on the diagonal ('I' option). It may also
    update an input unitary matrix, usually the output
    from CHETRD/CUNGTR or CHPTRD/CUPGTR ('V' option). It may
    also just compute eigenvalues ('N' option).

    CSTEMR factors S as Z D1 Z* , where Z is the unitary
    matrix of eigenvectors and D1 is a diagonal matrix with
    the eigenvalues on the diagonal ('I' option).  CSTEMR
    uses the Relatively Robust Representation whenever possible.

 When CCHKST2STG 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 Hermitian eigenroutines.  For each matrix, a number
 of tests will be performed:

 (1)     | A - V S V* | / ( |A| n ulp ) CHETRD( UPLO='U', ... )

 (2)     | I - UV* | / ( n ulp )        CUNGTR( UPLO='U', ... )

 (3)     | A - V S V* | / ( |A| n ulp ) CHETRD( UPLO='L', ... )
         replaced by | D1 - D2 | / ( |D1| ulp ) where D1 is the 
         eigenvalue matrix computed using S and D2 is the 
         eigenvalue matrix computed using S_2stage the output of
         CHETRD_2STAGE("N", "U",....). D1 and D2 are computed 
         via DSTEQR('N',...) 

 (4)     | I - UV* | / ( n ulp )        CUNGTR( UPLO='L', ... )
         replaced by | D1 - D3 | / ( |D1| ulp ) where D1 is the 
         eigenvalue matrix computed using S and D3 is the 
         eigenvalue matrix computed using S_2stage the output of
         CHETRD_2STAGE("N", "L",....). D1 and D3 are computed 
         via DSTEQR('N',...)  

 (5-8)   Same as 1-4, but for CHPTRD and CUPGTR.

 (9)     | S - Z D Z* | / ( |S| n ulp ) CSTEQR('V',...)

 (10)    | I - ZZ* | / ( n ulp )        CSTEQR('V',...)

 (11)    | D1 - D2 | / ( |D1| ulp )        CSTEQR('N',...)

 (12)    | D1 - D3 | / ( |D1| ulp )        SSTERF

 (13)    0 if the true eigenvalues (computed by sturm count)
         of S are within THRESH of
         those in D1.  2*THRESH if they are not.  (Tested using
         SSTECH)

 For S positive definite,

 (14)    | S - Z4 D4 Z4* | / ( |S| n ulp ) CPTEQR('V',...)

 (15)    | I - Z4 Z4* | / ( n ulp )        CPTEQR('V',...)

 (16)    | D4 - D5 | / ( 100 |D4| ulp )       CPTEQR('N',...)

 When S is also diagonally dominant by the factor gamma < 1,

 (17)    max | D4(i) - WR(i) | / ( |D4(i)| omega ) ,
          i
         omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4
                                              SSTEBZ( 'A', 'E', ...)

 (18)    | WA1 - D3 | / ( |D3| ulp )          SSTEBZ( 'A', 'E', ...)

 (19)    ( max { min | WA2(i)-WA3(j) | } +
            i     j
           max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp )
            i     j
                                              SSTEBZ( 'I', 'E', ...)

 (20)    | S - Y WA1 Y* | / ( |S| n ulp )  SSTEBZ, CSTEIN

 (21)    | I - Y Y* | / ( n ulp )          SSTEBZ, CSTEIN

 (22)    | S - Z D Z* | / ( |S| n ulp )    CSTEDC('I')

 (23)    | I - ZZ* | / ( n ulp )           CSTEDC('I')

 (24)    | S - Z D Z* | / ( |S| n ulp )    CSTEDC('V')

 (25)    | I - ZZ* | / ( n ulp )           CSTEDC('V')

 (26)    | D1 - D2 | / ( |D1| ulp )           CSTEDC('V') and
                                              CSTEDC('N')

 Test 27 is disabled at the moment because CSTEMR does not
 guarantee high relatvie accuracy.

 (27)    max | D6(i) - WR(i) | / ( |D6(i)| omega ) ,
          i
         omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4
                                              CSTEMR('V', 'A')

 (28)    max | D6(i) - WR(i) | / ( |D6(i)| omega ) ,
          i
         omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4
                                              CSTEMR('V', 'I')

 Tests 29 through 34 are disable at present because CSTEMR
 does not handle partial spectrum requests.

 (29)    | S - Z D Z* | / ( |S| n ulp )    CSTEMR('V', 'I')

 (30)    | I - ZZ* | / ( n ulp )           CSTEMR('V', 'I')

 (31)    ( max { min | WA2(i)-WA3(j) | } +
            i     j
           max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp )
            i     j
         CSTEMR('N', 'I') vs. CSTEMR('V', 'I')

 (32)    | S - Z D Z* | / ( |S| n ulp )    CSTEMR('V', 'V')

 (33)    | I - ZZ* | / ( n ulp )           CSTEMR('V', 'V')

 (34)    ( max { min | WA2(i)-WA3(j) | } +
            i     j
           max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp )
            i     j
         CSTEMR('N', 'V') vs. CSTEMR('V', 'V')

 (35)    | S - Z D Z* | / ( |S| n ulp )    CSTEMR('V', 'A')

 (36)    | I - ZZ* | / ( n ulp )           CSTEMR('V', 'A')

 (37)    ( max { min | WA2(i)-WA3(j) | } +
            i     j
           max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp )
            i     j
         CSTEMR('N', 'A') vs. CSTEMR('V', 'A')

 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 diagonal matrix with evenly spaced entries
      1, ..., ULP  and random signs.
      (ULP = (first number larger than 1) - 1 )
 (4)  A diagonal matrix with geometrically spaced entries
      1, ..., ULP  and random signs.
 (5)  A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
      and random signs.

 (6)  Same as (4), but multiplied by SQRT( overflow threshold )
 (7)  Same as (4), but multiplied by SQRT( underflow threshold )

 (8)  A matrix of the form  U* D U, where U is unitary and
      D has evenly spaced entries 1, ..., ULP with random signs
      on the diagonal.

 (9)  A matrix of the form  U* D U, where U is unitary and
      D has geometrically spaced entries 1, ..., ULP with random
      signs on the diagonal.

 (10) A matrix of the form  U* D U, where U is unitary and
      D has "clustered" entries 1, ULP,..., ULP with random
      signs on the diagonal.

 (11) Same as (8), but multiplied by SQRT( overflow threshold )
 (12) Same as (8), but multiplied by SQRT( underflow threshold )

 (13) Hermitian matrix with random entries chosen from (-1,1).
 (14) Same as (13), but multiplied by SQRT( overflow threshold )
 (15) Same as (13), but multiplied by SQRT( underflow threshold )
 (16) Same as (8), but diagonal elements are all positive.
 (17) Same as (9), but diagonal elements are all positive.
 (18) Same as (10), but diagonal elements are all positive.
 (19) Same as (16), but multiplied by SQRT( overflow threshold )
 (20) Same as (16), but multiplied by SQRT( underflow threshold )
 (21) A diagonally dominant tridiagonal matrix with geometrically
      spaced diagonal entries 1, ..., ULP.
Parameters
[in]NSIZES
          NSIZES is INTEGER
          The number of sizes of matrices to use.  If it is zero,
          CCHKST2STG 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, CCHKST2STG
          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 CCHKST2STG to continue the same random number
          sequence.
[in]THRESH
          THRESH is 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.
[in]NOUNIT
          NOUNIT is INTEGER
          The FORTRAN unit number for printing out error messages
          (e.g., if a routine returns IINFO not equal to 0.)
[in,out]A
          A is COMPLEX array of
                                  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.  It must be at
          least 1 and at least max( NN ).
[out]AP
          AP is COMPLEX array of
                      dimension( max(NN)*max(NN+1)/2 )
          The matrix A stored in packed format.
[out]SD
          SD is REAL array of
                             dimension( max(NN) )
          The diagonal of the tridiagonal matrix computed by CHETRD.
          On exit, SD and SE contain the tridiagonal form of the
          matrix in A.
[out]SE
          SE is REAL array of
                             dimension( max(NN) )
          The off-diagonal of the tridiagonal matrix computed by
          CHETRD.  On exit, SD and SE contain the tridiagonal form of
          the matrix in A.
[out]D1
          D1 is REAL array of
                             dimension( max(NN) )
          The eigenvalues of A, as computed by CSTEQR simlutaneously
          with Z.  On exit, the eigenvalues in D1 correspond with the
          matrix in A.
[out]D2
          D2 is REAL array of
                             dimension( max(NN) )
          The eigenvalues of A, as computed by CSTEQR if Z is not
          computed.  On exit, the eigenvalues in D2 correspond with
          the matrix in A.
[out]D3
          D3 is REAL array of
                             dimension( max(NN) )
          The eigenvalues of A, as computed by SSTERF.  On exit, the
          eigenvalues in D3 correspond with the matrix in A.
[out]D4
          D4 is REAL array of
                             dimension( max(NN) )
          The eigenvalues of A, as computed by CPTEQR(V).
          CPTEQR factors S as  Z4 D4 Z4*
          On exit, the eigenvalues in D4 correspond with the matrix in A.
[out]D5
          D5 is REAL array of
                             dimension( max(NN) )
          The eigenvalues of A, as computed by CPTEQR(N)
          when Z is not computed. On exit, the
          eigenvalues in D4 correspond with the matrix in A.
[out]WA1
          WA1 is REAL array of
                             dimension( max(NN) )
          All eigenvalues of A, computed to high
          absolute accuracy, with different range options.
          as computed by SSTEBZ.
[out]WA2
          WA2 is REAL array of
                             dimension( max(NN) )
          Selected eigenvalues of A, computed to high
          absolute accuracy, with different range options.
          as computed by SSTEBZ.
          Choose random values for IL and IU, and ask for the
          IL-th through IU-th eigenvalues.
[out]WA3
          WA3 is REAL array of
                             dimension( max(NN) )
          Selected eigenvalues of A, computed to high
          absolute accuracy, with different range options.
          as computed by SSTEBZ.
          Determine the values VL and VU of the IL-th and IU-th
          eigenvalues and ask for all eigenvalues in this range.
[out]WR
          WR is REAL array of
                             dimension( max(NN) )
          All eigenvalues of A, computed to high
          absolute accuracy, with different options.
          as computed by SSTEBZ.
[out]U
          U is COMPLEX array of
                             dimension( LDU, max(NN) ).
          The unitary matrix computed by CHETRD + CUNGTR.
[in]LDU
          LDU is INTEGER
          The leading dimension of U, Z, and V.  It must be at least 1
          and at least max( NN ).
[out]V
          V is COMPLEX array of
                             dimension( LDU, max(NN) ).
          The Housholder vectors computed by CHETRD in reducing A to
          tridiagonal form.  The vectors computed with UPLO='U' are
          in the upper triangle, and the vectors computed with UPLO='L'
          are in the lower triangle.  (As described in CHETRD, the
          sub- and superdiagonal are not set to 1, although the
          true Householder vector has a 1 in that position.  The
          routines that use V, such as CUNGTR, set those entries to
          1 before using them, and then restore them later.)
[out]VP
          VP is COMPLEX array of
                      dimension( max(NN)*max(NN+1)/2 )
          The matrix V stored in packed format.
[out]TAU
          TAU is COMPLEX array of
                             dimension( max(NN) )
          The Householder factors computed by CHETRD in reducing A
          to tridiagonal form.
[out]Z
          Z is COMPLEX array of
                             dimension( LDU, max(NN) ).
          The unitary matrix of eigenvectors computed by CSTEQR,
          CPTEQR, and CSTEIN.
[out]WORK
          WORK is COMPLEX array of
                      dimension( LWORK )
[in]LWORK
          LWORK is INTEGER
          The number of entries in WORK.  This must be at least
          1 + 4 * Nmax + 2 * Nmax * lg Nmax + 3 * Nmax**2
          where Nmax = max( NN(j), 2 ) and lg = log base 2.
[out]IWORK
          IWORK is INTEGER array,
          Workspace.
[out]LIWORK
          LIWORK is INTEGER
          The number of entries in IWORK.  This must be at least
                  6 + 6*Nmax + 5 * Nmax * lg Nmax
          where Nmax = max( NN(j), 2 ) and lg = log base 2.
[out]RWORK
          RWORK is REAL array
[in]LRWORK
          LRWORK is INTEGER
          The number of entries in LRWORK (dimension( ??? )
[out]RESULT
          RESULT is REAL array, dimension (26)
          The values computed by the tests described above.
          The values are currently limited to 1/ulp, to avoid
          overflow.
[out]INFO
          INFO is INTEGER
          If 0, then everything ran OK.
           -1: NSIZES < 0
           -2: Some NN(j) < 0
           -3: NTYPES < 0
           -5: THRESH < 0
           -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
          -23: LDU < 1 or LDU < NMAX.
          -29: LWORK too small.
          If  CLATMR, CLATMS, CHETRD, CUNGTR, CSTEQR, SSTERF,
              or CUNMC2 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.
       NTEST           The number of tests performed, or which can
                       be performed so far, for the current matrix.
       NTESTT          The total number of tests performed so far.
       NBLOCK          Blocksize as returned by ENVIR.
       NMAX            Largest value in NN.
       NMATS           The number of matrices generated so far.
       NERRS           The number of tests which have exceeded THRESH
                       so far.
       COND, 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  Square roots of the previous 2 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) )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 620 of file cchkst2stg.f.

625*
626* -- LAPACK test routine --
627* -- LAPACK is a software package provided by Univ. of Tennessee, --
628* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
629*
630* .. Scalar Arguments ..
631 INTEGER INFO, LDA, LDU, LIWORK, LRWORK, LWORK, NOUNIT,
632 $ NSIZES, NTYPES
633 REAL THRESH
634* ..
635* .. Array Arguments ..
636 LOGICAL DOTYPE( * )
637 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
638 REAL D1( * ), D2( * ), D3( * ), D4( * ), D5( * ),
639 $ RESULT( * ), RWORK( * ), SD( * ), SE( * ),
640 $ WA1( * ), WA2( * ), WA3( * ), WR( * )
641 COMPLEX A( LDA, * ), AP( * ), TAU( * ), U( LDU, * ),
642 $ V( LDU, * ), VP( * ), WORK( * ), Z( LDU, * )
643* ..
644*
645* =====================================================================
646*
647* .. Parameters ..
648 REAL ZERO, ONE, TWO, EIGHT, TEN, HUN
649 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0,
650 $ eight = 8.0e0, ten = 10.0e0, hun = 100.0e0 )
651 COMPLEX CZERO, CONE
652 parameter( czero = ( 0.0e+0, 0.0e+0 ),
653 $ cone = ( 1.0e+0, 0.0e+0 ) )
654 REAL HALF
655 parameter( half = one / two )
656 INTEGER MAXTYP
657 parameter( maxtyp = 21 )
658 LOGICAL CRANGE
659 parameter( crange = .false. )
660 LOGICAL CREL
661 parameter( crel = .false. )
662* ..
663* .. Local Scalars ..
664 LOGICAL BADNN, TRYRAC
665 INTEGER I, IINFO, IL, IMODE, INDE, INDRWK, ITEMP,
666 $ ITYPE, IU, J, JC, JR, JSIZE, JTYPE, LGN,
667 $ LIWEDC, LOG2UI, LRWEDC, LWEDC, M, M2, M3,
668 $ MTYPES, N, NAP, NBLOCK, NERRS, NMATS, NMAX,
669 $ NSPLIT, NTEST, NTESTT, LH, LW
670 REAL ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL,
671 $ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP,
672 $ ULPINV, UNFL, VL, VU
673* ..
674* .. Local Arrays ..
675 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
676 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
677 $ KTYPE( MAXTYP )
678 REAL DUMMA( 1 )
679* ..
680* .. External Functions ..
681 INTEGER ILAENV
682 REAL SLAMCH, SLARND, SSXT1
683 EXTERNAL ilaenv, slamch, slarnd, ssxt1
684* ..
685* .. External Subroutines ..
686 EXTERNAL scopy, slabad, slasum, sstebz, sstech, ssterf,
691* ..
692* .. Intrinsic Functions ..
693 INTRINSIC abs, real, conjg, int, log, max, min, sqrt
694* ..
695* .. Data statements ..
696 DATA ktype / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8,
697 $ 8, 8, 9, 9, 9, 9, 9, 10 /
698 DATA kmagn / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
699 $ 2, 3, 1, 1, 1, 2, 3, 1 /
700 DATA kmode / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
701 $ 0, 0, 4, 3, 1, 4, 4, 3 /
702* ..
703* .. Executable Statements ..
704*
705* Keep ftnchek happy
706 idumma( 1 ) = 1
707*
708* Check for errors
709*
710 ntestt = 0
711 info = 0
712*
713* Important constants
714*
715 badnn = .false.
716 tryrac = .true.
717 nmax = 1
718 DO 10 j = 1, nsizes
719 nmax = max( nmax, nn( j ) )
720 IF( nn( j ).LT.0 )
721 $ badnn = .true.
722 10 CONTINUE
723*
724 nblock = ilaenv( 1, 'CHETRD', 'L', nmax, -1, -1, -1 )
725 nblock = min( nmax, max( 1, nblock ) )
726*
727* Check for errors
728*
729 IF( nsizes.LT.0 ) THEN
730 info = -1
731 ELSE IF( badnn ) THEN
732 info = -2
733 ELSE IF( ntypes.LT.0 ) THEN
734 info = -3
735 ELSE IF( lda.LT.nmax ) THEN
736 info = -9
737 ELSE IF( ldu.LT.nmax ) THEN
738 info = -23
739 ELSE IF( 2*max( 2, nmax )**2.GT.lwork ) THEN
740 info = -29
741 END IF
742*
743 IF( info.NE.0 ) THEN
744 CALL xerbla( 'CCHKST2STG', -info )
745 RETURN
746 END IF
747*
748* Quick return if possible
749*
750 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
751 $ RETURN
752*
753* More Important constants
754*
755 unfl = slamch( 'Safe minimum' )
756 ovfl = one / unfl
757 CALL slabad( unfl, ovfl )
758 ulp = slamch( 'Epsilon' )*slamch( 'Base' )
759 ulpinv = one / ulp
760 log2ui = int( log( ulpinv ) / log( two ) )
761 rtunfl = sqrt( unfl )
762 rtovfl = sqrt( ovfl )
763*
764* Loop over sizes, types
765*
766 DO 20 i = 1, 4
767 iseed2( i ) = iseed( i )
768 20 CONTINUE
769 nerrs = 0
770 nmats = 0
771*
772 DO 310 jsize = 1, nsizes
773 n = nn( jsize )
774 IF( n.GT.0 ) THEN
775 lgn = int( log( real( n ) ) / log( two ) )
776 IF( 2**lgn.LT.n )
777 $ lgn = lgn + 1
778 IF( 2**lgn.LT.n )
779 $ lgn = lgn + 1
780 lwedc = 1 + 4*n + 2*n*lgn + 4*n**2
781 lrwedc = 1 + 3*n + 2*n*lgn + 4*n**2
782 liwedc = 6 + 6*n + 5*n*lgn
783 ELSE
784 lwedc = 8
785 lrwedc = 7
786 liwedc = 12
787 END IF
788 nap = ( n*( n+1 ) ) / 2
789 aninv = one / real( max( 1, n ) )
790*
791 IF( nsizes.NE.1 ) THEN
792 mtypes = min( maxtyp, ntypes )
793 ELSE
794 mtypes = min( maxtyp+1, ntypes )
795 END IF
796*
797 DO 300 jtype = 1, mtypes
798 IF( .NOT.dotype( jtype ) )
799 $ GO TO 300
800 nmats = nmats + 1
801 ntest = 0
802*
803 DO 30 j = 1, 4
804 ioldsd( j ) = iseed( j )
805 30 CONTINUE
806*
807* Compute "A"
808*
809* Control parameters:
810*
811* KMAGN KMODE KTYPE
812* =1 O(1) clustered 1 zero
813* =2 large clustered 2 identity
814* =3 small exponential (none)
815* =4 arithmetic diagonal, (w/ eigenvalues)
816* =5 random log Hermitian, w/ eigenvalues
817* =6 random (none)
818* =7 random diagonal
819* =8 random Hermitian
820* =9 positive definite
821* =10 diagonally dominant tridiagonal
822*
823 IF( mtypes.GT.maxtyp )
824 $ GO TO 100
825*
826 itype = ktype( jtype )
827 imode = kmode( jtype )
828*
829* Compute norm
830*
831 GO TO ( 40, 50, 60 )kmagn( jtype )
832*
833 40 CONTINUE
834 anorm = one
835 GO TO 70
836*
837 50 CONTINUE
838 anorm = ( rtovfl*ulp )*aninv
839 GO TO 70
840*
841 60 CONTINUE
842 anorm = rtunfl*n*ulpinv
843 GO TO 70
844*
845 70 CONTINUE
846*
847 CALL claset( 'Full', lda, n, czero, czero, a, lda )
848 iinfo = 0
849 IF( jtype.LE.15 ) THEN
850 cond = ulpinv
851 ELSE
852 cond = ulpinv*aninv / ten
853 END IF
854*
855* Special Matrices -- Identity & Jordan block
856*
857* Zero
858*
859 IF( itype.EQ.1 ) THEN
860 iinfo = 0
861*
862 ELSE IF( itype.EQ.2 ) THEN
863*
864* Identity
865*
866 DO 80 jc = 1, n
867 a( jc, jc ) = anorm
868 80 CONTINUE
869*
870 ELSE IF( itype.EQ.4 ) THEN
871*
872* Diagonal Matrix, [Eigen]values Specified
873*
874 CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
875 $ anorm, 0, 0, 'N', a, lda, work, iinfo )
876*
877*
878 ELSE IF( itype.EQ.5 ) THEN
879*
880* Hermitian, eigenvalues specified
881*
882 CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
883 $ anorm, n, n, 'N', a, lda, work, iinfo )
884*
885 ELSE IF( itype.EQ.7 ) THEN
886*
887* Diagonal, random eigenvalues
888*
889 CALL clatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
890 $ 'T', 'N', work( n+1 ), 1, one,
891 $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
892 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
893*
894 ELSE IF( itype.EQ.8 ) THEN
895*
896* Hermitian, random eigenvalues
897*
898 CALL clatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
899 $ 'T', 'N', work( n+1 ), 1, one,
900 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
901 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
902*
903 ELSE IF( itype.EQ.9 ) THEN
904*
905* Positive definite, eigenvalues specified.
906*
907 CALL clatms( n, n, 'S', iseed, 'P', rwork, imode, cond,
908 $ anorm, n, n, 'N', a, lda, work, iinfo )
909*
910 ELSE IF( itype.EQ.10 ) THEN
911*
912* Positive definite tridiagonal, eigenvalues specified.
913*
914 CALL clatms( n, n, 'S', iseed, 'P', rwork, imode, cond,
915 $ anorm, 1, 1, 'N', a, lda, work, iinfo )
916 DO 90 i = 2, n
917 temp1 = abs( a( i-1, i ) )
918 temp2 = sqrt( abs( a( i-1, i-1 )*a( i, i ) ) )
919 IF( temp1.GT.half*temp2 ) THEN
920 a( i-1, i ) = a( i-1, i )*
921 $ ( half*temp2 / ( unfl+temp1 ) )
922 a( i, i-1 ) = conjg( a( i-1, i ) )
923 END IF
924 90 CONTINUE
925*
926 ELSE
927*
928 iinfo = 1
929 END IF
930*
931 IF( iinfo.NE.0 ) THEN
932 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
933 $ ioldsd
934 info = abs( iinfo )
935 RETURN
936 END IF
937*
938 100 CONTINUE
939*
940* Call CHETRD and CUNGTR to compute S and U from
941* upper triangle.
942*
943 CALL clacpy( 'U', n, n, a, lda, v, ldu )
944*
945 ntest = 1
946 CALL chetrd( 'U', n, v, ldu, sd, se, tau, work, lwork,
947 $ iinfo )
948*
949 IF( iinfo.NE.0 ) THEN
950 WRITE( nounit, fmt = 9999 )'CHETRD(U)', iinfo, n, jtype,
951 $ ioldsd
952 info = abs( iinfo )
953 IF( iinfo.LT.0 ) THEN
954 RETURN
955 ELSE
956 result( 1 ) = ulpinv
957 GO TO 280
958 END IF
959 END IF
960*
961 CALL clacpy( 'U', n, n, v, ldu, u, ldu )
962*
963 ntest = 2
964 CALL cungtr( 'U', n, u, ldu, tau, work, lwork, iinfo )
965 IF( iinfo.NE.0 ) THEN
966 WRITE( nounit, fmt = 9999 )'CUNGTR(U)', iinfo, n, jtype,
967 $ ioldsd
968 info = abs( iinfo )
969 IF( iinfo.LT.0 ) THEN
970 RETURN
971 ELSE
972 result( 2 ) = ulpinv
973 GO TO 280
974 END IF
975 END IF
976*
977* Do tests 1 and 2
978*
979 CALL chet21( 2, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
980 $ ldu, tau, work, rwork, result( 1 ) )
981 CALL chet21( 3, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
982 $ ldu, tau, work, rwork, result( 2 ) )
983*
984* Compute D1 the eigenvalues resulting from the tridiagonal
985* form using the standard 1-stage algorithm and use it as a
986* reference to compare with the 2-stage technique
987*
988* Compute D1 from the 1-stage and used as reference for the
989* 2-stage
990*
991 CALL scopy( n, sd, 1, d1, 1 )
992 IF( n.GT.0 )
993 $ CALL scopy( n-1, se, 1, rwork, 1 )
994*
995 CALL csteqr( 'N', n, d1, rwork, work, ldu, rwork( n+1 ),
996 $ iinfo )
997 IF( iinfo.NE.0 ) THEN
998 WRITE( nounit, fmt = 9999 )'CSTEQR(N)', iinfo, n, jtype,
999 $ ioldsd
1000 info = abs( iinfo )
1001 IF( iinfo.LT.0 ) THEN
1002 RETURN
1003 ELSE
1004 result( 3 ) = ulpinv
1005 GO TO 280
1006 END IF
1007 END IF
1008*
1009* 2-STAGE TRD Upper case is used to compute D2.
1010* Note to set SD and SE to zero to be sure not reusing
1011* the one from above. Compare it with D1 computed
1012* using the 1-stage.
1013*
1014 CALL slaset( 'Full', n, 1, zero, zero, sd, n )
1015 CALL slaset( 'Full', n, 1, zero, zero, se, n )
1016 CALL clacpy( 'U', n, n, a, lda, v, ldu )
1017 lh = max(1, 4*n)
1018 lw = lwork - lh
1019 CALL chetrd_2stage( 'N', "U", n, v, ldu, sd, se, tau,
1020 $ work, lh, work( lh+1 ), lw, iinfo )
1021*
1022* Compute D2 from the 2-stage Upper case
1023*
1024 CALL scopy( n, sd, 1, d2, 1 )
1025 IF( n.GT.0 )
1026 $ CALL scopy( n-1, se, 1, rwork, 1 )
1027*
1028 ntest = 3
1029 CALL csteqr( 'N', n, d2, rwork, work, ldu, rwork( n+1 ),
1030 $ iinfo )
1031 IF( iinfo.NE.0 ) THEN
1032 WRITE( nounit, fmt = 9999 )'CSTEQR(N)', iinfo, n, jtype,
1033 $ ioldsd
1034 info = abs( iinfo )
1035 IF( iinfo.LT.0 ) THEN
1036 RETURN
1037 ELSE
1038 result( 3 ) = ulpinv
1039 GO TO 280
1040 END IF
1041 END IF
1042*
1043* 2-STAGE TRD Lower case is used to compute D3.
1044* Note to set SD and SE to zero to be sure not reusing
1045* the one from above. Compare it with D1 computed
1046* using the 1-stage.
1047*
1048 CALL slaset( 'Full', n, 1, zero, zero, sd, n )
1049 CALL slaset( 'Full', n, 1, zero, zero, se, n )
1050 CALL clacpy( 'L', n, n, a, lda, v, ldu )
1051 CALL chetrd_2stage( 'N', "L", n, v, ldu, sd, se, tau,
1052 $ work, lh, work( lh+1 ), lw, iinfo )
1053*
1054* Compute D3 from the 2-stage Upper case
1055*
1056 CALL scopy( n, sd, 1, d3, 1 )
1057 IF( n.GT.0 )
1058 $ CALL scopy( n-1, se, 1, rwork, 1 )
1059*
1060 ntest = 4
1061 CALL csteqr( 'N', n, d3, rwork, work, ldu, rwork( n+1 ),
1062 $ iinfo )
1063 IF( iinfo.NE.0 ) THEN
1064 WRITE( nounit, fmt = 9999 )'CSTEQR(N)', iinfo, n, jtype,
1065 $ ioldsd
1066 info = abs( iinfo )
1067 IF( iinfo.LT.0 ) THEN
1068 RETURN
1069 ELSE
1070 result( 4 ) = ulpinv
1071 GO TO 280
1072 END IF
1073 END IF
1074*
1075* Do Tests 3 and 4 which are similar to 11 and 12 but with the
1076* D1 computed using the standard 1-stage reduction as reference
1077*
1078 ntest = 4
1079 temp1 = zero
1080 temp2 = zero
1081 temp3 = zero
1082 temp4 = zero
1083*
1084 DO 151 j = 1, n
1085 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1086 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1087 temp3 = max( temp3, abs( d1( j ) ), abs( d3( j ) ) )
1088 temp4 = max( temp4, abs( d1( j )-d3( j ) ) )
1089 151 CONTINUE
1090*
1091 result( 3 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1092 result( 4 ) = temp4 / max( unfl, ulp*max( temp3, temp4 ) )
1093*
1094* Store the upper triangle of A in AP
1095*
1096 i = 0
1097 DO 120 jc = 1, n
1098 DO 110 jr = 1, jc
1099 i = i + 1
1100 ap( i ) = a( jr, jc )
1101 110 CONTINUE
1102 120 CONTINUE
1103*
1104* Call CHPTRD and CUPGTR to compute S and U from AP
1105*
1106 CALL ccopy( nap, ap, 1, vp, 1 )
1107*
1108 ntest = 5
1109 CALL chptrd( 'U', n, vp, sd, se, tau, iinfo )
1110*
1111 IF( iinfo.NE.0 ) THEN
1112 WRITE( nounit, fmt = 9999 )'CHPTRD(U)', iinfo, n, jtype,
1113 $ ioldsd
1114 info = abs( iinfo )
1115 IF( iinfo.LT.0 ) THEN
1116 RETURN
1117 ELSE
1118 result( 5 ) = ulpinv
1119 GO TO 280
1120 END IF
1121 END IF
1122*
1123 ntest = 6
1124 CALL cupgtr( 'U', n, vp, tau, u, ldu, work, iinfo )
1125 IF( iinfo.NE.0 ) THEN
1126 WRITE( nounit, fmt = 9999 )'CUPGTR(U)', iinfo, n, jtype,
1127 $ ioldsd
1128 info = abs( iinfo )
1129 IF( iinfo.LT.0 ) THEN
1130 RETURN
1131 ELSE
1132 result( 6 ) = ulpinv
1133 GO TO 280
1134 END IF
1135 END IF
1136*
1137* Do tests 5 and 6
1138*
1139 CALL chpt21( 2, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1140 $ work, rwork, result( 5 ) )
1141 CALL chpt21( 3, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1142 $ work, rwork, result( 6 ) )
1143*
1144* Store the lower triangle of A in AP
1145*
1146 i = 0
1147 DO 140 jc = 1, n
1148 DO 130 jr = jc, n
1149 i = i + 1
1150 ap( i ) = a( jr, jc )
1151 130 CONTINUE
1152 140 CONTINUE
1153*
1154* Call CHPTRD and CUPGTR to compute S and U from AP
1155*
1156 CALL ccopy( nap, ap, 1, vp, 1 )
1157*
1158 ntest = 7
1159 CALL chptrd( 'L', n, vp, sd, se, tau, iinfo )
1160*
1161 IF( iinfo.NE.0 ) THEN
1162 WRITE( nounit, fmt = 9999 )'CHPTRD(L)', iinfo, n, jtype,
1163 $ ioldsd
1164 info = abs( iinfo )
1165 IF( iinfo.LT.0 ) THEN
1166 RETURN
1167 ELSE
1168 result( 7 ) = ulpinv
1169 GO TO 280
1170 END IF
1171 END IF
1172*
1173 ntest = 8
1174 CALL cupgtr( 'L', n, vp, tau, u, ldu, work, iinfo )
1175 IF( iinfo.NE.0 ) THEN
1176 WRITE( nounit, fmt = 9999 )'CUPGTR(L)', iinfo, n, jtype,
1177 $ ioldsd
1178 info = abs( iinfo )
1179 IF( iinfo.LT.0 ) THEN
1180 RETURN
1181 ELSE
1182 result( 8 ) = ulpinv
1183 GO TO 280
1184 END IF
1185 END IF
1186*
1187 CALL chpt21( 2, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1188 $ work, rwork, result( 7 ) )
1189 CALL chpt21( 3, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1190 $ work, rwork, result( 8 ) )
1191*
1192* Call CSTEQR to compute D1, D2, and Z, do tests.
1193*
1194* Compute D1 and Z
1195*
1196 CALL scopy( n, sd, 1, d1, 1 )
1197 IF( n.GT.0 )
1198 $ CALL scopy( n-1, se, 1, rwork, 1 )
1199 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1200*
1201 ntest = 9
1202 CALL csteqr( 'V', n, d1, rwork, z, ldu, rwork( n+1 ),
1203 $ iinfo )
1204 IF( iinfo.NE.0 ) THEN
1205 WRITE( nounit, fmt = 9999 )'CSTEQR(V)', iinfo, n, jtype,
1206 $ ioldsd
1207 info = abs( iinfo )
1208 IF( iinfo.LT.0 ) THEN
1209 RETURN
1210 ELSE
1211 result( 9 ) = ulpinv
1212 GO TO 280
1213 END IF
1214 END IF
1215*
1216* Compute D2
1217*
1218 CALL scopy( n, sd, 1, d2, 1 )
1219 IF( n.GT.0 )
1220 $ CALL scopy( n-1, se, 1, rwork, 1 )
1221*
1222 ntest = 11
1223 CALL csteqr( 'N', n, d2, rwork, work, ldu, rwork( n+1 ),
1224 $ iinfo )
1225 IF( iinfo.NE.0 ) THEN
1226 WRITE( nounit, fmt = 9999 )'CSTEQR(N)', iinfo, n, jtype,
1227 $ ioldsd
1228 info = abs( iinfo )
1229 IF( iinfo.LT.0 ) THEN
1230 RETURN
1231 ELSE
1232 result( 11 ) = ulpinv
1233 GO TO 280
1234 END IF
1235 END IF
1236*
1237* Compute D3 (using PWK method)
1238*
1239 CALL scopy( n, sd, 1, d3, 1 )
1240 IF( n.GT.0 )
1241 $ CALL scopy( n-1, se, 1, rwork, 1 )
1242*
1243 ntest = 12
1244 CALL ssterf( n, d3, rwork, iinfo )
1245 IF( iinfo.NE.0 ) THEN
1246 WRITE( nounit, fmt = 9999 )'SSTERF', iinfo, n, jtype,
1247 $ ioldsd
1248 info = abs( iinfo )
1249 IF( iinfo.LT.0 ) THEN
1250 RETURN
1251 ELSE
1252 result( 12 ) = ulpinv
1253 GO TO 280
1254 END IF
1255 END IF
1256*
1257* Do Tests 9 and 10
1258*
1259 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1260 $ result( 9 ) )
1261*
1262* Do Tests 11 and 12
1263*
1264 temp1 = zero
1265 temp2 = zero
1266 temp3 = zero
1267 temp4 = zero
1268*
1269 DO 150 j = 1, n
1270 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1271 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1272 temp3 = max( temp3, abs( d1( j ) ), abs( d3( j ) ) )
1273 temp4 = max( temp4, abs( d1( j )-d3( j ) ) )
1274 150 CONTINUE
1275*
1276 result( 11 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1277 result( 12 ) = temp4 / max( unfl, ulp*max( temp3, temp4 ) )
1278*
1279* Do Test 13 -- Sturm Sequence Test of Eigenvalues
1280* Go up by factors of two until it succeeds
1281*
1282 ntest = 13
1283 temp1 = thresh*( half-ulp )
1284*
1285 DO 160 j = 0, log2ui
1286 CALL sstech( n, sd, se, d1, temp1, rwork, iinfo )
1287 IF( iinfo.EQ.0 )
1288 $ GO TO 170
1289 temp1 = temp1*two
1290 160 CONTINUE
1291*
1292 170 CONTINUE
1293 result( 13 ) = temp1
1294*
1295* For positive definite matrices ( JTYPE.GT.15 ) call CPTEQR
1296* and do tests 14, 15, and 16 .
1297*
1298 IF( jtype.GT.15 ) THEN
1299*
1300* Compute D4 and Z4
1301*
1302 CALL scopy( n, sd, 1, d4, 1 )
1303 IF( n.GT.0 )
1304 $ CALL scopy( n-1, se, 1, rwork, 1 )
1305 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1306*
1307 ntest = 14
1308 CALL cpteqr( 'V', n, d4, rwork, z, ldu, rwork( n+1 ),
1309 $ iinfo )
1310 IF( iinfo.NE.0 ) THEN
1311 WRITE( nounit, fmt = 9999 )'CPTEQR(V)', iinfo, n,
1312 $ jtype, ioldsd
1313 info = abs( iinfo )
1314 IF( iinfo.LT.0 ) THEN
1315 RETURN
1316 ELSE
1317 result( 14 ) = ulpinv
1318 GO TO 280
1319 END IF
1320 END IF
1321*
1322* Do Tests 14 and 15
1323*
1324 CALL cstt21( n, 0, sd, se, d4, dumma, z, ldu, work,
1325 $ rwork, result( 14 ) )
1326*
1327* Compute D5
1328*
1329 CALL scopy( n, sd, 1, d5, 1 )
1330 IF( n.GT.0 )
1331 $ CALL scopy( n-1, se, 1, rwork, 1 )
1332*
1333 ntest = 16
1334 CALL cpteqr( 'N', n, d5, rwork, z, ldu, rwork( n+1 ),
1335 $ iinfo )
1336 IF( iinfo.NE.0 ) THEN
1337 WRITE( nounit, fmt = 9999 )'CPTEQR(N)', iinfo, n,
1338 $ jtype, ioldsd
1339 info = abs( iinfo )
1340 IF( iinfo.LT.0 ) THEN
1341 RETURN
1342 ELSE
1343 result( 16 ) = ulpinv
1344 GO TO 280
1345 END IF
1346 END IF
1347*
1348* Do Test 16
1349*
1350 temp1 = zero
1351 temp2 = zero
1352 DO 180 j = 1, n
1353 temp1 = max( temp1, abs( d4( j ) ), abs( d5( j ) ) )
1354 temp2 = max( temp2, abs( d4( j )-d5( j ) ) )
1355 180 CONTINUE
1356*
1357 result( 16 ) = temp2 / max( unfl,
1358 $ hun*ulp*max( temp1, temp2 ) )
1359 ELSE
1360 result( 14 ) = zero
1361 result( 15 ) = zero
1362 result( 16 ) = zero
1363 END IF
1364*
1365* Call SSTEBZ with different options and do tests 17-18.
1366*
1367* If S is positive definite and diagonally dominant,
1368* ask for all eigenvalues with high relative accuracy.
1369*
1370 vl = zero
1371 vu = zero
1372 il = 0
1373 iu = 0
1374 IF( jtype.EQ.21 ) THEN
1375 ntest = 17
1376 abstol = unfl + unfl
1377 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se,
1378 $ m, nsplit, wr, iwork( 1 ), iwork( n+1 ),
1379 $ rwork, iwork( 2*n+1 ), iinfo )
1380 IF( iinfo.NE.0 ) THEN
1381 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,rel)', iinfo, n,
1382 $ jtype, ioldsd
1383 info = abs( iinfo )
1384 IF( iinfo.LT.0 ) THEN
1385 RETURN
1386 ELSE
1387 result( 17 ) = ulpinv
1388 GO TO 280
1389 END IF
1390 END IF
1391*
1392* Do test 17
1393*
1394 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1395 $ ( one-half )**4
1396*
1397 temp1 = zero
1398 DO 190 j = 1, n
1399 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1400 $ ( abstol+abs( d4( j ) ) ) )
1401 190 CONTINUE
1402*
1403 result( 17 ) = temp1 / temp2
1404 ELSE
1405 result( 17 ) = zero
1406 END IF
1407*
1408* Now ask for all eigenvalues with high absolute accuracy.
1409*
1410 ntest = 18
1411 abstol = unfl + unfl
1412 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se, m,
1413 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), rwork,
1414 $ iwork( 2*n+1 ), iinfo )
1415 IF( iinfo.NE.0 ) THEN
1416 WRITE( nounit, fmt = 9999 )'SSTEBZ(A)', iinfo, n, jtype,
1417 $ ioldsd
1418 info = abs( iinfo )
1419 IF( iinfo.LT.0 ) THEN
1420 RETURN
1421 ELSE
1422 result( 18 ) = ulpinv
1423 GO TO 280
1424 END IF
1425 END IF
1426*
1427* Do test 18
1428*
1429 temp1 = zero
1430 temp2 = zero
1431 DO 200 j = 1, n
1432 temp1 = max( temp1, abs( d3( j ) ), abs( wa1( j ) ) )
1433 temp2 = max( temp2, abs( d3( j )-wa1( j ) ) )
1434 200 CONTINUE
1435*
1436 result( 18 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1437*
1438* Choose random values for IL and IU, and ask for the
1439* IL-th through IU-th eigenvalues.
1440*
1441 ntest = 19
1442 IF( n.LE.1 ) THEN
1443 il = 1
1444 iu = n
1445 ELSE
1446 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1447 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1448 IF( iu.LT.il ) THEN
1449 itemp = iu
1450 iu = il
1451 il = itemp
1452 END IF
1453 END IF
1454*
1455 CALL sstebz( 'I', 'E', n, vl, vu, il, iu, abstol, sd, se,
1456 $ m2, nsplit, wa2, iwork( 1 ), iwork( n+1 ),
1457 $ rwork, iwork( 2*n+1 ), iinfo )
1458 IF( iinfo.NE.0 ) THEN
1459 WRITE( nounit, fmt = 9999 )'SSTEBZ(I)', iinfo, n, jtype,
1460 $ ioldsd
1461 info = abs( iinfo )
1462 IF( iinfo.LT.0 ) THEN
1463 RETURN
1464 ELSE
1465 result( 19 ) = ulpinv
1466 GO TO 280
1467 END IF
1468 END IF
1469*
1470* Determine the values VL and VU of the IL-th and IU-th
1471* eigenvalues and ask for all eigenvalues in this range.
1472*
1473 IF( n.GT.0 ) THEN
1474 IF( il.NE.1 ) THEN
1475 vl = wa1( il ) - max( half*( wa1( il )-wa1( il-1 ) ),
1476 $ ulp*anorm, two*rtunfl )
1477 ELSE
1478 vl = wa1( 1 ) - max( half*( wa1( n )-wa1( 1 ) ),
1479 $ ulp*anorm, two*rtunfl )
1480 END IF
1481 IF( iu.NE.n ) THEN
1482 vu = wa1( iu ) + max( half*( wa1( iu+1 )-wa1( iu ) ),
1483 $ ulp*anorm, two*rtunfl )
1484 ELSE
1485 vu = wa1( n ) + max( half*( wa1( n )-wa1( 1 ) ),
1486 $ ulp*anorm, two*rtunfl )
1487 END IF
1488 ELSE
1489 vl = zero
1490 vu = one
1491 END IF
1492*
1493 CALL sstebz( 'V', 'E', n, vl, vu, il, iu, abstol, sd, se,
1494 $ m3, nsplit, wa3, iwork( 1 ), iwork( n+1 ),
1495 $ rwork, iwork( 2*n+1 ), iinfo )
1496 IF( iinfo.NE.0 ) THEN
1497 WRITE( nounit, fmt = 9999 )'SSTEBZ(V)', iinfo, n, jtype,
1498 $ ioldsd
1499 info = abs( iinfo )
1500 IF( iinfo.LT.0 ) THEN
1501 RETURN
1502 ELSE
1503 result( 19 ) = ulpinv
1504 GO TO 280
1505 END IF
1506 END IF
1507*
1508 IF( m3.EQ.0 .AND. n.NE.0 ) THEN
1509 result( 19 ) = ulpinv
1510 GO TO 280
1511 END IF
1512*
1513* Do test 19
1514*
1515 temp1 = ssxt1( 1, wa2, m2, wa3, m3, abstol, ulp, unfl )
1516 temp2 = ssxt1( 1, wa3, m3, wa2, m2, abstol, ulp, unfl )
1517 IF( n.GT.0 ) THEN
1518 temp3 = max( abs( wa1( n ) ), abs( wa1( 1 ) ) )
1519 ELSE
1520 temp3 = zero
1521 END IF
1522*
1523 result( 19 ) = ( temp1+temp2 ) / max( unfl, temp3*ulp )
1524*
1525* Call CSTEIN to compute eigenvectors corresponding to
1526* eigenvalues in WA1. (First call SSTEBZ again, to make sure
1527* it returns these eigenvalues in the correct order.)
1528*
1529 ntest = 21
1530 CALL sstebz( 'A', 'B', n, vl, vu, il, iu, abstol, sd, se, m,
1531 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), rwork,
1532 $ iwork( 2*n+1 ), iinfo )
1533 IF( iinfo.NE.0 ) THEN
1534 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,B)', iinfo, n,
1535 $ jtype, ioldsd
1536 info = abs( iinfo )
1537 IF( iinfo.LT.0 ) THEN
1538 RETURN
1539 ELSE
1540 result( 20 ) = ulpinv
1541 result( 21 ) = ulpinv
1542 GO TO 280
1543 END IF
1544 END IF
1545*
1546 CALL cstein( n, sd, se, m, wa1, iwork( 1 ), iwork( n+1 ), z,
1547 $ ldu, rwork, iwork( 2*n+1 ), iwork( 3*n+1 ),
1548 $ iinfo )
1549 IF( iinfo.NE.0 ) THEN
1550 WRITE( nounit, fmt = 9999 )'CSTEIN', iinfo, n, jtype,
1551 $ ioldsd
1552 info = abs( iinfo )
1553 IF( iinfo.LT.0 ) THEN
1554 RETURN
1555 ELSE
1556 result( 20 ) = ulpinv
1557 result( 21 ) = ulpinv
1558 GO TO 280
1559 END IF
1560 END IF
1561*
1562* Do tests 20 and 21
1563*
1564 CALL cstt21( n, 0, sd, se, wa1, dumma, z, ldu, work, rwork,
1565 $ result( 20 ) )
1566*
1567* Call CSTEDC(I) to compute D1 and Z, do tests.
1568*
1569* Compute D1 and Z
1570*
1571 inde = 1
1572 indrwk = inde + n
1573 CALL scopy( n, sd, 1, d1, 1 )
1574 IF( n.GT.0 )
1575 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1576 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1577*
1578 ntest = 22
1579 CALL cstedc( 'I', n, d1, rwork( inde ), z, ldu, work, lwedc,
1580 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1581 IF( iinfo.NE.0 ) THEN
1582 WRITE( nounit, fmt = 9999 )'CSTEDC(I)', iinfo, n, jtype,
1583 $ ioldsd
1584 info = abs( iinfo )
1585 IF( iinfo.LT.0 ) THEN
1586 RETURN
1587 ELSE
1588 result( 22 ) = ulpinv
1589 GO TO 280
1590 END IF
1591 END IF
1592*
1593* Do Tests 22 and 23
1594*
1595 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1596 $ result( 22 ) )
1597*
1598* Call CSTEDC(V) to compute D1 and Z, do tests.
1599*
1600* Compute D1 and Z
1601*
1602 CALL scopy( n, sd, 1, d1, 1 )
1603 IF( n.GT.0 )
1604 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1605 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1606*
1607 ntest = 24
1608 CALL cstedc( 'V', n, d1, rwork( inde ), z, ldu, work, lwedc,
1609 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1610 IF( iinfo.NE.0 ) THEN
1611 WRITE( nounit, fmt = 9999 )'CSTEDC(V)', iinfo, n, jtype,
1612 $ ioldsd
1613 info = abs( iinfo )
1614 IF( iinfo.LT.0 ) THEN
1615 RETURN
1616 ELSE
1617 result( 24 ) = ulpinv
1618 GO TO 280
1619 END IF
1620 END IF
1621*
1622* Do Tests 24 and 25
1623*
1624 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1625 $ result( 24 ) )
1626*
1627* Call CSTEDC(N) to compute D2, do tests.
1628*
1629* Compute D2
1630*
1631 CALL scopy( n, sd, 1, d2, 1 )
1632 IF( n.GT.0 )
1633 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1634 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1635*
1636 ntest = 26
1637 CALL cstedc( 'N', n, d2, rwork( inde ), z, ldu, work, lwedc,
1638 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1639 IF( iinfo.NE.0 ) THEN
1640 WRITE( nounit, fmt = 9999 )'CSTEDC(N)', iinfo, n, jtype,
1641 $ ioldsd
1642 info = abs( iinfo )
1643 IF( iinfo.LT.0 ) THEN
1644 RETURN
1645 ELSE
1646 result( 26 ) = ulpinv
1647 GO TO 280
1648 END IF
1649 END IF
1650*
1651* Do Test 26
1652*
1653 temp1 = zero
1654 temp2 = zero
1655*
1656 DO 210 j = 1, n
1657 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1658 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1659 210 CONTINUE
1660*
1661 result( 26 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1662*
1663* Only test CSTEMR if IEEE compliant
1664*
1665 IF( ilaenv( 10, 'CSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND.
1666 $ ilaenv( 11, 'CSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN
1667*
1668* Call CSTEMR, do test 27 (relative eigenvalue accuracy)
1669*
1670* If S is positive definite and diagonally dominant,
1671* ask for all eigenvalues with high relative accuracy.
1672*
1673 vl = zero
1674 vu = zero
1675 il = 0
1676 iu = 0
1677 IF( jtype.EQ.21 .AND. crel ) THEN
1678 ntest = 27
1679 abstol = unfl + unfl
1680 CALL cstemr( 'V', 'A', n, sd, se, vl, vu, il, iu,
1681 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1682 $ rwork, lrwork, iwork( 2*n+1 ), lwork-2*n,
1683 $ iinfo )
1684 IF( iinfo.NE.0 ) THEN
1685 WRITE( nounit, fmt = 9999 )'CSTEMR(V,A,rel)',
1686 $ iinfo, n, jtype, ioldsd
1687 info = abs( iinfo )
1688 IF( iinfo.LT.0 ) THEN
1689 RETURN
1690 ELSE
1691 result( 27 ) = ulpinv
1692 GO TO 270
1693 END IF
1694 END IF
1695*
1696* Do test 27
1697*
1698 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1699 $ ( one-half )**4
1700*
1701 temp1 = zero
1702 DO 220 j = 1, n
1703 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1704 $ ( abstol+abs( d4( j ) ) ) )
1705 220 CONTINUE
1706*
1707 result( 27 ) = temp1 / temp2
1708*
1709 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1710 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1711 IF( iu.LT.il ) THEN
1712 itemp = iu
1713 iu = il
1714 il = itemp
1715 END IF
1716*
1717 IF( crange ) THEN
1718 ntest = 28
1719 abstol = unfl + unfl
1720 CALL cstemr( 'V', 'I', n, sd, se, vl, vu, il, iu,
1721 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1722 $ rwork, lrwork, iwork( 2*n+1 ),
1723 $ lwork-2*n, iinfo )
1724*
1725 IF( iinfo.NE.0 ) THEN
1726 WRITE( nounit, fmt = 9999 )'CSTEMR(V,I,rel)',
1727 $ iinfo, n, jtype, ioldsd
1728 info = abs( iinfo )
1729 IF( iinfo.LT.0 ) THEN
1730 RETURN
1731 ELSE
1732 result( 28 ) = ulpinv
1733 GO TO 270
1734 END IF
1735 END IF
1736*
1737* Do test 28
1738*
1739 temp2 = two*( two*n-one )*ulp*
1740 $ ( one+eight*half**2 ) / ( one-half )**4
1741*
1742 temp1 = zero
1743 DO 230 j = il, iu
1744 temp1 = max( temp1, abs( wr( j-il+1 )-d4( n-j+
1745 $ 1 ) ) / ( abstol+abs( wr( j-il+1 ) ) ) )
1746 230 CONTINUE
1747*
1748 result( 28 ) = temp1 / temp2
1749 ELSE
1750 result( 28 ) = zero
1751 END IF
1752 ELSE
1753 result( 27 ) = zero
1754 result( 28 ) = zero
1755 END IF
1756*
1757* Call CSTEMR(V,I) to compute D1 and Z, do tests.
1758*
1759* Compute D1 and Z
1760*
1761 CALL scopy( n, sd, 1, d5, 1 )
1762 IF( n.GT.0 )
1763 $ CALL scopy( n-1, se, 1, rwork, 1 )
1764 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1765*
1766 IF( crange ) THEN
1767 ntest = 29
1768 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1769 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1770 IF( iu.LT.il ) THEN
1771 itemp = iu
1772 iu = il
1773 il = itemp
1774 END IF
1775 CALL cstemr( 'V', 'I', n, d5, rwork, vl, vu, il, iu,
1776 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1777 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1778 $ liwork-2*n, iinfo )
1779 IF( iinfo.NE.0 ) THEN
1780 WRITE( nounit, fmt = 9999 )'CSTEMR(V,I)', iinfo,
1781 $ n, jtype, ioldsd
1782 info = abs( iinfo )
1783 IF( iinfo.LT.0 ) THEN
1784 RETURN
1785 ELSE
1786 result( 29 ) = ulpinv
1787 GO TO 280
1788 END IF
1789 END IF
1790*
1791* Do Tests 29 and 30
1792*
1793* Call CSTEMR to compute D2, do tests.
1794*
1795* Compute D2
1796*
1797 CALL scopy( n, sd, 1, d5, 1 )
1798 IF( n.GT.0 )
1799 $ CALL scopy( n-1, se, 1, rwork, 1 )
1800*
1801 ntest = 31
1802 CALL cstemr( 'N', 'I', n, d5, rwork, vl, vu, il, iu,
1803 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1804 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1805 $ liwork-2*n, iinfo )
1806 IF( iinfo.NE.0 ) THEN
1807 WRITE( nounit, fmt = 9999 )'CSTEMR(N,I)', iinfo,
1808 $ n, jtype, ioldsd
1809 info = abs( iinfo )
1810 IF( iinfo.LT.0 ) THEN
1811 RETURN
1812 ELSE
1813 result( 31 ) = ulpinv
1814 GO TO 280
1815 END IF
1816 END IF
1817*
1818* Do Test 31
1819*
1820 temp1 = zero
1821 temp2 = zero
1822*
1823 DO 240 j = 1, iu - il + 1
1824 temp1 = max( temp1, abs( d1( j ) ),
1825 $ abs( d2( j ) ) )
1826 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1827 240 CONTINUE
1828*
1829 result( 31 ) = temp2 / max( unfl,
1830 $ ulp*max( temp1, temp2 ) )
1831*
1832* Call CSTEMR(V,V) to compute D1 and Z, do tests.
1833*
1834* Compute D1 and Z
1835*
1836 CALL scopy( n, sd, 1, d5, 1 )
1837 IF( n.GT.0 )
1838 $ CALL scopy( n-1, se, 1, rwork, 1 )
1839 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1840*
1841 ntest = 32
1842*
1843 IF( n.GT.0 ) THEN
1844 IF( il.NE.1 ) THEN
1845 vl = d2( il ) - max( half*
1846 $ ( d2( il )-d2( il-1 ) ), ulp*anorm,
1847 $ two*rtunfl )
1848 ELSE
1849 vl = d2( 1 ) - max( half*( d2( n )-d2( 1 ) ),
1850 $ ulp*anorm, two*rtunfl )
1851 END IF
1852 IF( iu.NE.n ) THEN
1853 vu = d2( iu ) + max( half*
1854 $ ( d2( iu+1 )-d2( iu ) ), ulp*anorm,
1855 $ two*rtunfl )
1856 ELSE
1857 vu = d2( n ) + max( half*( d2( n )-d2( 1 ) ),
1858 $ ulp*anorm, two*rtunfl )
1859 END IF
1860 ELSE
1861 vl = zero
1862 vu = one
1863 END IF
1864*
1865 CALL cstemr( 'V', 'V', n, d5, rwork, vl, vu, il, iu,
1866 $ m, d1, z, ldu, m, iwork( 1 ), tryrac,
1867 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1868 $ liwork-2*n, iinfo )
1869 IF( iinfo.NE.0 ) THEN
1870 WRITE( nounit, fmt = 9999 )'CSTEMR(V,V)', iinfo,
1871 $ n, jtype, ioldsd
1872 info = abs( iinfo )
1873 IF( iinfo.LT.0 ) THEN
1874 RETURN
1875 ELSE
1876 result( 32 ) = ulpinv
1877 GO TO 280
1878 END IF
1879 END IF
1880*
1881* Do Tests 32 and 33
1882*
1883 CALL cstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work,
1884 $ m, rwork, result( 32 ) )
1885*
1886* Call CSTEMR to compute D2, do tests.
1887*
1888* Compute D2
1889*
1890 CALL scopy( n, sd, 1, d5, 1 )
1891 IF( n.GT.0 )
1892 $ CALL scopy( n-1, se, 1, rwork, 1 )
1893*
1894 ntest = 34
1895 CALL cstemr( 'N', 'V', n, d5, rwork, vl, vu, il, iu,
1896 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1897 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1898 $ liwork-2*n, iinfo )
1899 IF( iinfo.NE.0 ) THEN
1900 WRITE( nounit, fmt = 9999 )'CSTEMR(N,V)', iinfo,
1901 $ n, jtype, ioldsd
1902 info = abs( iinfo )
1903 IF( iinfo.LT.0 ) THEN
1904 RETURN
1905 ELSE
1906 result( 34 ) = ulpinv
1907 GO TO 280
1908 END IF
1909 END IF
1910*
1911* Do Test 34
1912*
1913 temp1 = zero
1914 temp2 = zero
1915*
1916 DO 250 j = 1, iu - il + 1
1917 temp1 = max( temp1, abs( d1( j ) ),
1918 $ abs( d2( j ) ) )
1919 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1920 250 CONTINUE
1921*
1922 result( 34 ) = temp2 / max( unfl,
1923 $ ulp*max( temp1, temp2 ) )
1924 ELSE
1925 result( 29 ) = zero
1926 result( 30 ) = zero
1927 result( 31 ) = zero
1928 result( 32 ) = zero
1929 result( 33 ) = zero
1930 result( 34 ) = zero
1931 END IF
1932*
1933* Call CSTEMR(V,A) to compute D1 and Z, do tests.
1934*
1935* Compute D1 and Z
1936*
1937 CALL scopy( n, sd, 1, d5, 1 )
1938 IF( n.GT.0 )
1939 $ CALL scopy( n-1, se, 1, rwork, 1 )
1940*
1941 ntest = 35
1942*
1943 CALL cstemr( 'V', 'A', n, d5, rwork, vl, vu, il, iu,
1944 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1945 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1946 $ liwork-2*n, iinfo )
1947 IF( iinfo.NE.0 ) THEN
1948 WRITE( nounit, fmt = 9999 )'CSTEMR(V,A)', iinfo, n,
1949 $ jtype, ioldsd
1950 info = abs( iinfo )
1951 IF( iinfo.LT.0 ) THEN
1952 RETURN
1953 ELSE
1954 result( 35 ) = ulpinv
1955 GO TO 280
1956 END IF
1957 END IF
1958*
1959* Do Tests 35 and 36
1960*
1961 CALL cstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work, m,
1962 $ rwork, result( 35 ) )
1963*
1964* Call CSTEMR to compute D2, do tests.
1965*
1966* Compute D2
1967*
1968 CALL scopy( n, sd, 1, d5, 1 )
1969 IF( n.GT.0 )
1970 $ CALL scopy( n-1, se, 1, rwork, 1 )
1971*
1972 ntest = 37
1973 CALL cstemr( 'N', 'A', n, d5, rwork, vl, vu, il, iu,
1974 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1975 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1976 $ liwork-2*n, iinfo )
1977 IF( iinfo.NE.0 ) THEN
1978 WRITE( nounit, fmt = 9999 )'CSTEMR(N,A)', iinfo, n,
1979 $ jtype, ioldsd
1980 info = abs( iinfo )
1981 IF( iinfo.LT.0 ) THEN
1982 RETURN
1983 ELSE
1984 result( 37 ) = ulpinv
1985 GO TO 280
1986 END IF
1987 END IF
1988*
1989* Do Test 37
1990*
1991 temp1 = zero
1992 temp2 = zero
1993*
1994 DO 260 j = 1, n
1995 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1996 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1997 260 CONTINUE
1998*
1999 result( 37 ) = temp2 / max( unfl,
2000 $ ulp*max( temp1, temp2 ) )
2001 END IF
2002 270 CONTINUE
2003 280 CONTINUE
2004 ntestt = ntestt + ntest
2005*
2006* End of Loop -- Check for RESULT(j) > THRESH
2007*
2008* Print out tests which fail.
2009*
2010 DO 290 jr = 1, ntest
2011 IF( result( jr ).GE.thresh ) THEN
2012*
2013* If this is the first test to fail,
2014* print a header to the data file.
2015*
2016 IF( nerrs.EQ.0 ) THEN
2017 WRITE( nounit, fmt = 9998 )'CST'
2018 WRITE( nounit, fmt = 9997 )
2019 WRITE( nounit, fmt = 9996 )
2020 WRITE( nounit, fmt = 9995 )'Hermitian'
2021 WRITE( nounit, fmt = 9994 )
2022*
2023* Tests performed
2024*
2025 WRITE( nounit, fmt = 9987 )
2026 END IF
2027 nerrs = nerrs + 1
2028 IF( result( jr ).LT.10000.0e0 ) THEN
2029 WRITE( nounit, fmt = 9989 )n, jtype, ioldsd, jr,
2030 $ result( jr )
2031 ELSE
2032 WRITE( nounit, fmt = 9988 )n, jtype, ioldsd, jr,
2033 $ result( jr )
2034 END IF
2035 END IF
2036 290 CONTINUE
2037 300 CONTINUE
2038 310 CONTINUE
2039*
2040* Summary
2041*
2042 CALL slasum( 'CST', nounit, nerrs, ntestt )
2043 RETURN
2044*
2045 9999 FORMAT( ' CCHKST2STG: ', a, ' returned INFO=', i6, '.', / 9x,
2046 $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
2047*
2048 9998 FORMAT( / 1x, a3, ' -- Complex Hermitian eigenvalue problem' )
2049 9997 FORMAT( ' Matrix types (see CCHKST2STG for details): ' )
2050*
2051 9996 FORMAT( / ' Special Matrices:',
2052 $ / ' 1=Zero matrix. ',
2053 $ ' 5=Diagonal: clustered entries.',
2054 $ / ' 2=Identity matrix. ',
2055 $ ' 6=Diagonal: large, evenly spaced.',
2056 $ / ' 3=Diagonal: evenly spaced entries. ',
2057 $ ' 7=Diagonal: small, evenly spaced.',
2058 $ / ' 4=Diagonal: geometr. spaced entries.' )
2059 9995 FORMAT( ' Dense ', a, ' Matrices:',
2060 $ / ' 8=Evenly spaced eigenvals. ',
2061 $ ' 12=Small, evenly spaced eigenvals.',
2062 $ / ' 9=Geometrically spaced eigenvals. ',
2063 $ ' 13=Matrix with random O(1) entries.',
2064 $ / ' 10=Clustered eigenvalues. ',
2065 $ ' 14=Matrix with large random entries.',
2066 $ / ' 11=Large, evenly spaced eigenvals. ',
2067 $ ' 15=Matrix with small random entries.' )
2068 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues',
2069 $ / ' 17=Positive definite, geometrically spaced eigenvlaues',
2070 $ / ' 18=Positive definite, clustered eigenvalues',
2071 $ / ' 19=Positive definite, small evenly spaced eigenvalues',
2072 $ / ' 20=Positive definite, large evenly spaced eigenvalues',
2073 $ / ' 21=Diagonally dominant tridiagonal, geometrically',
2074 $ ' spaced eigenvalues' )
2075*
2076 9989 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
2077 $ 4( i4, ',' ), ' result ', i3, ' is', 0p, f8.2 )
2078 9988 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
2079 $ 4( i4, ',' ), ' result ', i3, ' is', 1p, e10.3 )
2080*
2081 9987 FORMAT( / 'Test performed: see CCHKST2STG for details.', / )
2082*
2083* End of CCHKST2STG
2084*
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:74
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:110
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:86
subroutine sstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
SSTEBZ
Definition: sstebz.f:273
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:81
subroutine cstt21(N, KBAND, AD, AE, SD, SE, U, LDU, WORK, RWORK, RESULT)
CSTT21
Definition: cstt21.f:133
subroutine chet21(ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V, LDV, TAU, WORK, RWORK, RESULT)
CHET21
Definition: chet21.f:214
subroutine cstt22(N, M, KBAND, AD, AE, SD, SE, U, LDU, WORK, LDWORK, RWORK, RESULT)
CSTT22
Definition: cstt22.f:145
subroutine chpt21(ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP, TAU, WORK, RWORK, RESULT)
CHPT21
Definition: chpt21.f:228
subroutine clatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
CLATMS
Definition: clatms.f:332
subroutine clatmr(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)
CLATMR
Definition: clatmr.f:490
subroutine chetrd_2stage(VECT, UPLO, N, A, LDA, D, E, TAU, HOUS2, LHOUS2, WORK, LWORK, INFO)
CHETRD_2STAGE
subroutine chetrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
CHETRD
Definition: chetrd.f:192
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: claset.f:106
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:103
subroutine cupgtr(UPLO, N, AP, TAU, Q, LDQ, WORK, INFO)
CUPGTR
Definition: cupgtr.f:114
subroutine chptrd(UPLO, N, AP, D, E, TAU, INFO)
CHPTRD
Definition: chptrd.f:151
subroutine cstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
CSTEIN
Definition: cstein.f:182
subroutine csteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
CSTEQR
Definition: csteqr.f:132
subroutine cstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO)
CSTEDC
Definition: cstedc.f:212
subroutine cstemr(JOBZ, RANGE, N, D, E, VL, VU, IL, IU, M, W, Z, LDZ, NZC, ISUPPZ, TRYRAC, WORK, LWORK, IWORK, LIWORK, INFO)
CSTEMR
Definition: cstemr.f:338
subroutine cungtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
CUNGTR
Definition: cungtr.f:123
subroutine cpteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
CPTEQR
Definition: cpteqr.f:145
real function slarnd(IDIST, ISEED)
SLARND
Definition: slarnd.f:73
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:82
subroutine sstech(N, A, B, EIG, TOL, WORK, INFO)
SSTECH
Definition: sstech.f:101
real function ssxt1(IJOB, D1, N1, D2, N2, ABSTOL, ULP, UNFL)
SSXT1
Definition: ssxt1.f:106
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
subroutine slasum(TYPE, IOUNIT, IE, NRUN)
SLASUM
Definition: slasum.f:41
Here is the call graph for this function:
Here is the caller graph for this function: