626 lines
22 KiB
FortranFixed
626 lines
22 KiB
FortranFixed
SUBROUTINE SCHKSB( NSIZES, NN, NWDTHS, KK, NTYPES, DOTYPE, ISEED,
|
|
$ THRESH, NOUNIT, A, LDA, SD, SE, U, LDU, WORK,
|
|
$ LWORK, RESULT, INFO )
|
|
*
|
|
* -- LAPACK test routine (version 3.1) --
|
|
* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
|
|
* November 2006
|
|
*
|
|
* .. Scalar Arguments ..
|
|
INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES,
|
|
$ NWDTHS
|
|
REAL THRESH
|
|
* ..
|
|
* .. Array Arguments ..
|
|
LOGICAL DOTYPE( * )
|
|
INTEGER ISEED( 4 ), KK( * ), NN( * )
|
|
REAL A( LDA, * ), RESULT( * ), SD( * ), SE( * ),
|
|
$ U( LDU, * ), WORK( * )
|
|
* ..
|
|
*
|
|
* Purpose
|
|
* =======
|
|
*
|
|
* SCHKSB tests the reduction of a symmetric band matrix to tridiagonal
|
|
* form, used with the symmetric eigenvalue problem.
|
|
*
|
|
* SSBTRD factors a symmetric band matrix A as U S U' , where ' means
|
|
* transpose, S is symmetric tridiagonal, and U is orthogonal.
|
|
* SSBTRD can use either just the lower or just the upper triangle
|
|
* of A; SCHKSB checks both cases.
|
|
*
|
|
* When SCHKSB is called, a number of matrix "sizes" ("n's"), a number
|
|
* of bandwidths ("k's"), and a number of matrix "types" are
|
|
* specified. For each size ("n"), each bandwidth ("k") less than or
|
|
* equal to "n", and each type of matrix, one matrix will be generated
|
|
* and used to test the symmetric banded reduction routine. For each
|
|
* matrix, a number of tests will be performed:
|
|
*
|
|
* (1) | A - V S V' | / ( |A| n ulp ) computed by SSBTRD with
|
|
* UPLO='U'
|
|
*
|
|
* (2) | I - UU' | / ( n ulp )
|
|
*
|
|
* (3) | A - V S V' | / ( |A| n ulp ) computed by SSBTRD with
|
|
* UPLO='L'
|
|
*
|
|
* (4) | I - UU' | / ( n 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 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 orthogonal 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 orthogonal 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 orthogonal 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) Symmetric 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 )
|
|
*
|
|
* Arguments
|
|
* =========
|
|
*
|
|
* NSIZES (input) INTEGER
|
|
* The number of sizes of matrices to use. If it is zero,
|
|
* SCHKSB does nothing. It must be at least zero.
|
|
*
|
|
* NN (input) 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.
|
|
*
|
|
* NWDTHS (input) INTEGER
|
|
* The number of bandwidths to use. If it is zero,
|
|
* SCHKSB does nothing. It must be at least zero.
|
|
*
|
|
* KK (input) INTEGER array, dimension (NWDTHS)
|
|
* An array containing the bandwidths to be used for the band
|
|
* matrices. The values must be at least zero.
|
|
*
|
|
* NTYPES (input) INTEGER
|
|
* The number of elements in DOTYPE. If it is zero, SCHKSB
|
|
* 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. .
|
|
*
|
|
* DOTYPE (input) 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.
|
|
*
|
|
* ISEED (input/output) 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 SCHKSB to continue the same random number
|
|
* sequence.
|
|
*
|
|
* THRESH (input) 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.
|
|
*
|
|
* NOUNIT (input) INTEGER
|
|
* The FORTRAN unit number for printing out error messages
|
|
* (e.g., if a routine returns IINFO not equal to 0.)
|
|
*
|
|
* A (input/workspace) REAL array, dimension
|
|
* (LDA, max(NN))
|
|
* Used to hold the matrix whose eigenvalues are to be
|
|
* computed.
|
|
*
|
|
* LDA (input) INTEGER
|
|
* The leading dimension of A. It must be at least 2 (not 1!)
|
|
* and at least max( KK )+1.
|
|
*
|
|
* SD (workspace) REAL array, dimension (max(NN))
|
|
* Used to hold the diagonal of the tridiagonal matrix computed
|
|
* by SSBTRD.
|
|
*
|
|
* SE (workspace) REAL array, dimension (max(NN))
|
|
* Used to hold the off-diagonal of the tridiagonal matrix
|
|
* computed by SSBTRD.
|
|
*
|
|
* U (workspace) REAL array, dimension (LDU, max(NN))
|
|
* Used to hold the orthogonal matrix computed by SSBTRD.
|
|
*
|
|
* LDU (input) INTEGER
|
|
* The leading dimension of U. It must be at least 1
|
|
* and at least max( NN ).
|
|
*
|
|
* WORK (workspace) REAL array, dimension (LWORK)
|
|
*
|
|
* LWORK (input) INTEGER
|
|
* The number of entries in WORK. This must be at least
|
|
* max( LDA+1, max(NN)+1 )*max(NN).
|
|
*
|
|
* RESULT (output) REAL array, dimension (4)
|
|
* The values computed by the tests described above.
|
|
* The values are currently limited to 1/ulp, to avoid
|
|
* overflow.
|
|
*
|
|
* INFO (output) INTEGER
|
|
* If 0, then everything ran OK.
|
|
*
|
|
*-----------------------------------------------------------------------
|
|
*
|
|
* 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.
|
|
* 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) )
|
|
*
|
|
* =====================================================================
|
|
*
|
|
* .. Parameters ..
|
|
REAL ZERO, ONE, TWO, TEN
|
|
PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0,
|
|
$ TEN = 10.0E0 )
|
|
REAL HALF
|
|
PARAMETER ( HALF = ONE / TWO )
|
|
INTEGER MAXTYP
|
|
PARAMETER ( MAXTYP = 15 )
|
|
* ..
|
|
* .. Local Scalars ..
|
|
LOGICAL BADNN, BADNNB
|
|
INTEGER I, IINFO, IMODE, ITYPE, J, JC, JCOL, JR, JSIZE,
|
|
$ JTYPE, JWIDTH, K, KMAX, MTYPES, N, NERRS,
|
|
$ NMATS, NMAX, NTEST, NTESTT
|
|
REAL ANINV, ANORM, COND, OVFL, RTOVFL, RTUNFL,
|
|
$ TEMP1, ULP, ULPINV, UNFL
|
|
* ..
|
|
* .. Local Arrays ..
|
|
INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KMAGN( MAXTYP ),
|
|
$ KMODE( MAXTYP ), KTYPE( MAXTYP )
|
|
* ..
|
|
* .. External Functions ..
|
|
REAL SLAMCH
|
|
EXTERNAL SLAMCH
|
|
* ..
|
|
* .. External Subroutines ..
|
|
EXTERNAL SLACPY, SLASUM, SLATMR, SLATMS, SLASET, SSBT21,
|
|
$ SSBTRD, XERBLA
|
|
* ..
|
|
* .. Intrinsic Functions ..
|
|
INTRINSIC ABS, MAX, MIN, REAL, SQRT
|
|
* ..
|
|
* .. Data statements ..
|
|
DATA KTYPE / 1, 2, 5*4, 5*5, 3*8 /
|
|
DATA KMAGN / 2*1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
|
|
$ 2, 3 /
|
|
DATA KMODE / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
|
|
$ 0, 0 /
|
|
* ..
|
|
* .. Executable Statements ..
|
|
*
|
|
* Check for errors
|
|
*
|
|
NTESTT = 0
|
|
INFO = 0
|
|
*
|
|
* Important constants
|
|
*
|
|
BADNN = .FALSE.
|
|
NMAX = 1
|
|
DO 10 J = 1, NSIZES
|
|
NMAX = MAX( NMAX, NN( J ) )
|
|
IF( NN( J ).LT.0 )
|
|
$ BADNN = .TRUE.
|
|
10 CONTINUE
|
|
*
|
|
BADNNB = .FALSE.
|
|
KMAX = 0
|
|
DO 20 J = 1, NSIZES
|
|
KMAX = MAX( KMAX, KK( J ) )
|
|
IF( KK( J ).LT.0 )
|
|
$ BADNNB = .TRUE.
|
|
20 CONTINUE
|
|
KMAX = MIN( NMAX-1, KMAX )
|
|
*
|
|
* Check for errors
|
|
*
|
|
IF( NSIZES.LT.0 ) THEN
|
|
INFO = -1
|
|
ELSE IF( BADNN ) THEN
|
|
INFO = -2
|
|
ELSE IF( NWDTHS.LT.0 ) THEN
|
|
INFO = -3
|
|
ELSE IF( BADNNB ) THEN
|
|
INFO = -4
|
|
ELSE IF( NTYPES.LT.0 ) THEN
|
|
INFO = -5
|
|
ELSE IF( LDA.LT.KMAX+1 ) THEN
|
|
INFO = -11
|
|
ELSE IF( LDU.LT.NMAX ) THEN
|
|
INFO = -15
|
|
ELSE IF( ( MAX( LDA, NMAX )+1 )*NMAX.GT.LWORK ) THEN
|
|
INFO = -17
|
|
END IF
|
|
*
|
|
IF( INFO.NE.0 ) THEN
|
|
CALL XERBLA( 'SCHKSB', -INFO )
|
|
RETURN
|
|
END IF
|
|
*
|
|
* Quick return if possible
|
|
*
|
|
IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 .OR. NWDTHS.EQ.0 )
|
|
$ RETURN
|
|
*
|
|
* More Important constants
|
|
*
|
|
UNFL = SLAMCH( 'Safe minimum' )
|
|
OVFL = ONE / UNFL
|
|
ULP = SLAMCH( 'Epsilon' )*SLAMCH( 'Base' )
|
|
ULPINV = ONE / ULP
|
|
RTUNFL = SQRT( UNFL )
|
|
RTOVFL = SQRT( OVFL )
|
|
*
|
|
* Loop over sizes, types
|
|
*
|
|
NERRS = 0
|
|
NMATS = 0
|
|
*
|
|
DO 190 JSIZE = 1, NSIZES
|
|
N = NN( JSIZE )
|
|
ANINV = ONE / REAL( MAX( 1, N ) )
|
|
*
|
|
DO 180 JWIDTH = 1, NWDTHS
|
|
K = KK( JWIDTH )
|
|
IF( K.GT.N )
|
|
$ GO TO 180
|
|
K = MAX( 0, MIN( N-1, K ) )
|
|
*
|
|
IF( NSIZES.NE.1 ) THEN
|
|
MTYPES = MIN( MAXTYP, NTYPES )
|
|
ELSE
|
|
MTYPES = MIN( MAXTYP+1, NTYPES )
|
|
END IF
|
|
*
|
|
DO 170 JTYPE = 1, MTYPES
|
|
IF( .NOT.DOTYPE( JTYPE ) )
|
|
$ GO TO 170
|
|
NMATS = NMATS + 1
|
|
NTEST = 0
|
|
*
|
|
DO 30 J = 1, 4
|
|
IOLDSD( J ) = ISEED( J )
|
|
30 CONTINUE
|
|
*
|
|
* Compute "A".
|
|
* Store as "Upper"; later, we will copy to other format.
|
|
*
|
|
* Control parameters:
|
|
*
|
|
* KMAGN KMODE KTYPE
|
|
* =1 O(1) clustered 1 zero
|
|
* =2 large clustered 2 identity
|
|
* =3 small exponential (none)
|
|
* =4 arithmetic diagonal, (w/ eigenvalues)
|
|
* =5 random log symmetric, w/ eigenvalues
|
|
* =6 random (none)
|
|
* =7 random diagonal
|
|
* =8 random symmetric
|
|
* =9 positive definite
|
|
* =10 diagonally dominant tridiagonal
|
|
*
|
|
IF( MTYPES.GT.MAXTYP )
|
|
$ GO TO 100
|
|
*
|
|
ITYPE = KTYPE( JTYPE )
|
|
IMODE = KMODE( JTYPE )
|
|
*
|
|
* Compute norm
|
|
*
|
|
GO TO ( 40, 50, 60 )KMAGN( JTYPE )
|
|
*
|
|
40 CONTINUE
|
|
ANORM = ONE
|
|
GO TO 70
|
|
*
|
|
50 CONTINUE
|
|
ANORM = ( RTOVFL*ULP )*ANINV
|
|
GO TO 70
|
|
*
|
|
60 CONTINUE
|
|
ANORM = RTUNFL*N*ULPINV
|
|
GO TO 70
|
|
*
|
|
70 CONTINUE
|
|
*
|
|
CALL SLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA )
|
|
IINFO = 0
|
|
IF( JTYPE.LE.15 ) THEN
|
|
COND = ULPINV
|
|
ELSE
|
|
COND = ULPINV*ANINV / TEN
|
|
END IF
|
|
*
|
|
* Special Matrices -- Identity & Jordan block
|
|
*
|
|
* Zero
|
|
*
|
|
IF( ITYPE.EQ.1 ) THEN
|
|
IINFO = 0
|
|
*
|
|
ELSE IF( ITYPE.EQ.2 ) THEN
|
|
*
|
|
* Identity
|
|
*
|
|
DO 80 JCOL = 1, N
|
|
A( K+1, JCOL ) = ANORM
|
|
80 CONTINUE
|
|
*
|
|
ELSE IF( ITYPE.EQ.4 ) THEN
|
|
*
|
|
* Diagonal Matrix, [Eigen]values Specified
|
|
*
|
|
CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
|
|
$ ANORM, 0, 0, 'Q', A( K+1, 1 ), LDA,
|
|
$ WORK( N+1 ), IINFO )
|
|
*
|
|
ELSE IF( ITYPE.EQ.5 ) THEN
|
|
*
|
|
* Symmetric, eigenvalues specified
|
|
*
|
|
CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
|
|
$ ANORM, K, K, 'Q', A, LDA, WORK( N+1 ),
|
|
$ IINFO )
|
|
*
|
|
ELSE IF( ITYPE.EQ.7 ) THEN
|
|
*
|
|
* Diagonal, random eigenvalues
|
|
*
|
|
CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
|
|
$ 'T', 'N', WORK( N+1 ), 1, ONE,
|
|
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
|
|
$ ZERO, ANORM, 'Q', A( K+1, 1 ), LDA,
|
|
$ IDUMMA, IINFO )
|
|
*
|
|
ELSE IF( ITYPE.EQ.8 ) THEN
|
|
*
|
|
* Symmetric, random eigenvalues
|
|
*
|
|
CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
|
|
$ 'T', 'N', WORK( N+1 ), 1, ONE,
|
|
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, K, K,
|
|
$ ZERO, ANORM, 'Q', A, LDA, IDUMMA, IINFO )
|
|
*
|
|
ELSE IF( ITYPE.EQ.9 ) THEN
|
|
*
|
|
* Positive definite, eigenvalues specified.
|
|
*
|
|
CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND,
|
|
$ ANORM, K, K, 'Q', A, LDA, WORK( N+1 ),
|
|
$ IINFO )
|
|
*
|
|
ELSE IF( ITYPE.EQ.10 ) THEN
|
|
*
|
|
* Positive definite tridiagonal, eigenvalues specified.
|
|
*
|
|
IF( N.GT.1 )
|
|
$ K = MAX( 1, K )
|
|
CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND,
|
|
$ ANORM, 1, 1, 'Q', A( K, 1 ), LDA,
|
|
$ WORK( N+1 ), IINFO )
|
|
DO 90 I = 2, N
|
|
TEMP1 = ABS( A( K, I ) ) /
|
|
$ SQRT( ABS( A( K+1, I-1 )*A( K+1, I ) ) )
|
|
IF( TEMP1.GT.HALF ) THEN
|
|
A( K, I ) = HALF*SQRT( ABS( A( K+1,
|
|
$ I-1 )*A( K+1, I ) ) )
|
|
END IF
|
|
90 CONTINUE
|
|
*
|
|
ELSE
|
|
*
|
|
IINFO = 1
|
|
END IF
|
|
*
|
|
IF( IINFO.NE.0 ) THEN
|
|
WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N,
|
|
$ JTYPE, IOLDSD
|
|
INFO = ABS( IINFO )
|
|
RETURN
|
|
END IF
|
|
*
|
|
100 CONTINUE
|
|
*
|
|
* Call SSBTRD to compute S and U from upper triangle.
|
|
*
|
|
CALL SLACPY( ' ', K+1, N, A, LDA, WORK, LDA )
|
|
*
|
|
NTEST = 1
|
|
CALL SSBTRD( 'V', 'U', N, K, WORK, LDA, SD, SE, U, LDU,
|
|
$ WORK( LDA*N+1 ), IINFO )
|
|
*
|
|
IF( IINFO.NE.0 ) THEN
|
|
WRITE( NOUNIT, FMT = 9999 )'SSBTRD(U)', IINFO, N,
|
|
$ JTYPE, IOLDSD
|
|
INFO = ABS( IINFO )
|
|
IF( IINFO.LT.0 ) THEN
|
|
RETURN
|
|
ELSE
|
|
RESULT( 1 ) = ULPINV
|
|
GO TO 150
|
|
END IF
|
|
END IF
|
|
*
|
|
* Do tests 1 and 2
|
|
*
|
|
CALL SSBT21( 'Upper', N, K, 1, A, LDA, SD, SE, U, LDU,
|
|
$ WORK, RESULT( 1 ) )
|
|
*
|
|
* Convert A from Upper-Triangle-Only storage to
|
|
* Lower-Triangle-Only storage.
|
|
*
|
|
DO 120 JC = 1, N
|
|
DO 110 JR = 0, MIN( K, N-JC )
|
|
A( JR+1, JC ) = A( K+1-JR, JC+JR )
|
|
110 CONTINUE
|
|
120 CONTINUE
|
|
DO 140 JC = N + 1 - K, N
|
|
DO 130 JR = MIN( K, N-JC ) + 1, K
|
|
A( JR+1, JC ) = ZERO
|
|
130 CONTINUE
|
|
140 CONTINUE
|
|
*
|
|
* Call SSBTRD to compute S and U from lower triangle
|
|
*
|
|
CALL SLACPY( ' ', K+1, N, A, LDA, WORK, LDA )
|
|
*
|
|
NTEST = 3
|
|
CALL SSBTRD( 'V', 'L', N, K, WORK, LDA, SD, SE, U, LDU,
|
|
$ WORK( LDA*N+1 ), IINFO )
|
|
*
|
|
IF( IINFO.NE.0 ) THEN
|
|
WRITE( NOUNIT, FMT = 9999 )'SSBTRD(L)', IINFO, N,
|
|
$ JTYPE, IOLDSD
|
|
INFO = ABS( IINFO )
|
|
IF( IINFO.LT.0 ) THEN
|
|
RETURN
|
|
ELSE
|
|
RESULT( 3 ) = ULPINV
|
|
GO TO 150
|
|
END IF
|
|
END IF
|
|
NTEST = 4
|
|
*
|
|
* Do tests 3 and 4
|
|
*
|
|
CALL SSBT21( 'Lower', N, K, 1, A, LDA, SD, SE, U, LDU,
|
|
$ WORK, RESULT( 3 ) )
|
|
*
|
|
* End of Loop -- Check for RESULT(j) > THRESH
|
|
*
|
|
150 CONTINUE
|
|
NTESTT = NTESTT + NTEST
|
|
*
|
|
* Print out tests which fail.
|
|
*
|
|
DO 160 JR = 1, NTEST
|
|
IF( RESULT( JR ).GE.THRESH ) THEN
|
|
*
|
|
* If this is the first test to fail,
|
|
* print a header to the data file.
|
|
*
|
|
IF( NERRS.EQ.0 ) THEN
|
|
WRITE( NOUNIT, FMT = 9998 )'SSB'
|
|
WRITE( NOUNIT, FMT = 9997 )
|
|
WRITE( NOUNIT, FMT = 9996 )
|
|
WRITE( NOUNIT, FMT = 9995 )'Symmetric'
|
|
WRITE( NOUNIT, FMT = 9994 )'orthogonal', '''',
|
|
$ 'transpose', ( '''', J = 1, 4 )
|
|
END IF
|
|
NERRS = NERRS + 1
|
|
WRITE( NOUNIT, FMT = 9993 )N, K, IOLDSD, JTYPE,
|
|
$ JR, RESULT( JR )
|
|
END IF
|
|
160 CONTINUE
|
|
*
|
|
170 CONTINUE
|
|
180 CONTINUE
|
|
190 CONTINUE
|
|
*
|
|
* Summary
|
|
*
|
|
CALL SLASUM( 'SSB', NOUNIT, NERRS, NTESTT )
|
|
RETURN
|
|
*
|
|
9999 FORMAT( ' SCHKSB: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
|
|
$ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
|
|
*
|
|
9998 FORMAT( / 1X, A3,
|
|
$ ' -- Real Symmetric Banded Tridiagonal Reduction Routines' )
|
|
9997 FORMAT( ' Matrix types (see SCHKSB for details): ' )
|
|
*
|
|
9996 FORMAT( / ' Special Matrices:',
|
|
$ / ' 1=Zero matrix. ',
|
|
$ ' 5=Diagonal: clustered entries.',
|
|
$ / ' 2=Identity matrix. ',
|
|
$ ' 6=Diagonal: large, evenly spaced.',
|
|
$ / ' 3=Diagonal: evenly spaced entries. ',
|
|
$ ' 7=Diagonal: small, evenly spaced.',
|
|
$ / ' 4=Diagonal: geometr. spaced entries.' )
|
|
9995 FORMAT( ' Dense ', A, ' Banded Matrices:',
|
|
$ / ' 8=Evenly spaced eigenvals. ',
|
|
$ ' 12=Small, evenly spaced eigenvals.',
|
|
$ / ' 9=Geometrically spaced eigenvals. ',
|
|
$ ' 13=Matrix with random O(1) entries.',
|
|
$ / ' 10=Clustered eigenvalues. ',
|
|
$ ' 14=Matrix with large random entries.',
|
|
$ / ' 11=Large, evenly spaced eigenvals. ',
|
|
$ ' 15=Matrix with small random entries.' )
|
|
*
|
|
9994 FORMAT( / ' Tests performed: (S is Tridiag, U is ', A, ',',
|
|
$ / 20X, A, ' means ', A, '.', / ' UPLO=''U'':',
|
|
$ / ' 1= | A - U S U', A1, ' | / ( |A| n ulp ) ',
|
|
$ ' 2= | I - U U', A1, ' | / ( n ulp )', / ' UPLO=''L'':',
|
|
$ / ' 3= | A - U S U', A1, ' | / ( |A| n ulp ) ',
|
|
$ ' 4= | I - U U', A1, ' | / ( n ulp )' )
|
|
9993 FORMAT( ' N=', I5, ', K=', I4, ', seed=', 4( I4, ',' ), ' type ',
|
|
$ I2, ', test(', I2, ')=', G10.3 )
|
|
*
|
|
* End of SCHKSB
|
|
*
|
|
END
|