zchkhs.f
来自「famous linear algebra library (LAPACK) p」· F 代码 · 共 1,049 行 · 第 1/3 页
F
1,049 行
* 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 DLAFTS).
* 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.)
*
* =====================================================================
*
* .. Parameters ..
DOUBLE PRECISION ZERO, ONE
PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
COMPLEX*16 CZERO, CONE
PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ),
$ CONE = ( 1.0D+0, 0.0D+0 ) )
INTEGER MAXTYP
PARAMETER ( MAXTYP = 21 )
* ..
* .. Local Scalars ..
LOGICAL BADNN, MATCH
INTEGER I, IHI, IINFO, ILO, IMODE, IN, ITYPE, J, JCOL,
$ JJ, JSIZE, JTYPE, K, MTYPES, N, N1, NERRS,
$ NMATS, NMAX, NTEST, NTESTT
DOUBLE PRECISION ANINV, ANORM, COND, CONDS, OVFL, RTOVFL, RTULP,
$ RTULPI, RTUNFL, TEMP1, TEMP2, ULP, ULPINV, UNFL
* ..
* .. Local Arrays ..
INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
$ KMAGN( MAXTYP ), KMODE( MAXTYP ),
$ KTYPE( MAXTYP )
DOUBLE PRECISION DUMMA( 4 )
COMPLEX*16 CDUMMA( 4 )
* ..
* .. External Functions ..
DOUBLE PRECISION DLAMCH
EXTERNAL DLAMCH
* ..
* .. External Subroutines ..
EXTERNAL DLABAD, DLAFTS, DLASUM, XERBLA, ZCOPY, ZGEHRD,
$ ZGEMM, ZGET10, ZGET22, ZHSEIN, ZHSEQR, ZHST01,
$ ZLACPY, ZLASET, ZLATME, ZLATMR, ZLATMS, ZTREVC,
$ ZUNGHR, ZUNMHR
* ..
* .. Intrinsic Functions ..
INTRINSIC ABS, DBLE, MAX, MIN, SQRT
* ..
* .. Data statements ..
DATA KTYPE / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
DATA KMAGN / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
$ 3, 1, 2, 3 /
DATA KMODE / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
$ 1, 5, 5, 5, 4, 3, 1 /
DATA KCONDS / 3*0, 5*0, 4*1, 6*2, 3*0 /
* ..
* .. Executable Statements ..
*
* Check for errors
*
NTESTT = 0
INFO = 0
*
BADNN = .FALSE.
NMAX = 0
DO 10 J = 1, NSIZES
NMAX = MAX( NMAX, NN( J ) )
IF( NN( J ).LT.0 )
$ BADNN = .TRUE.
10 CONTINUE
*
* Check for errors
*
IF( NSIZES.LT.0 ) THEN
INFO = -1
ELSE IF( BADNN ) THEN
INFO = -2
ELSE IF( NTYPES.LT.0 ) THEN
INFO = -3
ELSE IF( THRESH.LT.ZERO ) THEN
INFO = -6
ELSE IF( LDA.LE.1 .OR. LDA.LT.NMAX ) THEN
INFO = -9
ELSE IF( LDU.LE.1 .OR. LDU.LT.NMAX ) THEN
INFO = -14
ELSE IF( 4*NMAX*NMAX+2.GT.NWORK ) THEN
INFO = -26
END IF
*
IF( INFO.NE.0 ) THEN
CALL XERBLA( 'ZCHKHS', -INFO )
RETURN
END IF
*
* Quick return if possible
*
IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
$ RETURN
*
* More important constants
*
UNFL = DLAMCH( 'Safe minimum' )
OVFL = DLAMCH( 'Overflow' )
CALL DLABAD( UNFL, OVFL )
ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' )
ULPINV = ONE / ULP
RTUNFL = SQRT( UNFL )
RTOVFL = SQRT( OVFL )
RTULP = SQRT( ULP )
RTULPI = ONE / RTULP
*
* Loop over sizes, types
*
NERRS = 0
NMATS = 0
*
DO 260 JSIZE = 1, NSIZES
N = NN( JSIZE )
N1 = MAX( 1, N )
ANINV = ONE / DBLE( N1 )
*
IF( NSIZES.NE.1 ) THEN
MTYPES = MIN( MAXTYP, NTYPES )
ELSE
MTYPES = MIN( MAXTYP+1, NTYPES )
END IF
*
DO 250 JTYPE = 1, MTYPES
IF( .NOT.DOTYPE( JTYPE ) )
$ GO TO 250
NMATS = NMATS + 1
NTEST = 0
*
* Save ISEED in case of an error.
*
DO 20 J = 1, 4
IOLDSD( J ) = ISEED( J )
20 CONTINUE
*
* Initialize RESULT
*
DO 30 J = 1, 14
RESULT( J ) = ZERO
30 CONTINUE
*
* Compute "A"
*
* Control parameters:
*
* KMAGN KCONDS KMODE KTYPE
* =1 O(1) 1 clustered 1 zero
* =2 large large clustered 2 identity
* =3 small exponential Jordan
* =4 arithmetic diagonal, (w/ eigenvalues)
* =5 random log hermitian, w/ eigenvalues
* =6 random general, w/ eigenvalues
* =7 random diagonal
* =8 random hermitian
* =9 random general
* =10 random triangular
*
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 ZLASET( 'Full', LDA, N, CZERO, CZERO, A, LDA )
IINFO = 0
COND = ULPINV
*
* Special Matrices
*
IF( ITYPE.EQ.1 ) THEN
*
* Zero
*
IINFO = 0
ELSE IF( ITYPE.EQ.2 ) THEN
*
* Identity
*
DO 80 JCOL = 1, N
A( JCOL, JCOL ) = ANORM
80 CONTINUE
*
ELSE IF( ITYPE.EQ.3 ) THEN
*
* Jordan Block
*
DO 90 JCOL = 1, N
A( JCOL, JCOL ) = ANORM
IF( JCOL.GT.1 )
$ A( JCOL, JCOL-1 ) = ONE
90 CONTINUE
*
ELSE IF( ITYPE.EQ.4 ) THEN
*
* Diagonal Matrix, [Eigen]values Specified
*
CALL ZLATMR( N, N, 'D', ISEED, 'N', WORK, IMODE, COND,
$ CONE, 'T', 'N', WORK( N+1 ), 1, ONE,
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
$ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
*
ELSE IF( ITYPE.EQ.5 ) THEN
*
* Hermitian, eigenvalues specified
*
CALL ZLATMS( N, N, 'D', ISEED, 'H', RWORK, IMODE, COND,
$ ANORM, N, N, 'N', A, LDA, WORK, IINFO )
*
ELSE IF( ITYPE.EQ.6 ) THEN
*
* General, eigenvalues specified
*
IF( KCONDS( JTYPE ).EQ.1 ) THEN
CONDS = ONE
ELSE IF( KCONDS( JTYPE ).EQ.2 ) THEN
CONDS = RTULPI
ELSE
CONDS = ZERO
END IF
*
CALL ZLATME( N, 'D', ISEED, WORK, IMODE, COND, CONE, ' ',
$ 'T', 'T', 'T', RWORK, 4, CONDS, N, N, ANORM,
$ A, LDA, WORK( N+1 ), IINFO )
*
ELSE IF( ITYPE.EQ.7 ) THEN
*
* Diagonal, random eigenvalues
*
CALL ZLATMR( N, N, 'D', ISEED, 'N', WORK, 6, ONE, CONE,
$ 'T', 'N', WORK( N+1 ), 1, ONE,
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
$ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
*
ELSE IF( ITYPE.EQ.8 ) THEN
*
* Hermitian, random eigenvalues
*
CALL ZLATMR( N, N, 'D', ISEED, 'H', WORK, 6, ONE, CONE,
$ 'T', 'N', WORK( N+1 ), 1, ONE,
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
$ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
*
ELSE IF( ITYPE.EQ.9 ) THEN
*
* General, random eigenvalues
*
CALL ZLATMR( N, N, 'D', ISEED, 'N', WORK, 6, ONE, CONE,
$ 'T', 'N', WORK( N+1 ), 1, ONE,
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
$ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
*
ELSE IF( ITYPE.EQ.10 ) THEN
*
* Triangular, random eigenvalues
*
CALL ZLATMR( N, N, 'D', ISEED, 'N', WORK, 6, ONE, CONE,
$ 'T', 'N', WORK( N+1 ), 1, ONE,
$ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, 0,
$ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
*
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 ZGEHRD to compute H and U, do tests.
*
CALL ZLACPY( ' ', N, N, A, LDA, H, LDA )
NTEST = 1
*
ILO = 1
IHI = N
*
CALL ZGEHRD( N, ILO, IHI, H, LDA, WORK, WORK( N+1 ),
$ NWORK-N, IINFO )
*
IF( IINFO.NE.0 ) THEN
RESULT( 1 ) = ULPINV
WRITE( NOUNIT, FMT = 9999 )'ZGEHRD', IINFO, N, JTYPE,
$ IOLDSD
INFO = ABS( IINFO )
GO TO 240
END IF
*
DO 120 J = 1, N - 1
UU( J+1, J ) = CZERO
DO 110 I = J + 2, N
U( I, J ) = H( I, J )
UU( I, J ) = H( I, J )
H( I, J ) = CZERO
110 CONTINUE
120 CONTINUE
CALL ZCOPY( N-1, WORK, 1, TAU, 1 )
CALL ZUNGHR( N, ILO, IHI, U, LDU, WORK, WORK( N+1 ),
$ NWORK-N, IINFO )
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?