schkst.f
来自「famous linear algebra library (LAPACK) p」· F 代码 · 共 1,804 行 · 第 1/5 页
F
1,804 行
*
TEMP2 = TWO*( TWO*N-ONE )*ULP*
$ ( ONE+EIGHT*HALF**2 ) / ( ONE-HALF )**4
*
TEMP1 = ZERO
DO 230 J = IL, IU
TEMP1 = MAX( TEMP1, ABS( WR( J-IL+1 )-D4( N-J+
$ 1 ) ) / ( ABSTOL+ABS( WR( J-IL+1 ) ) ) )
230 CONTINUE
*
RESULT( 28 ) = TEMP1 / TEMP2
ELSE
RESULT( 28 ) = ZERO
END IF
ELSE
RESULT( 27 ) = ZERO
RESULT( 28 ) = ZERO
END IF
*
* Call SSTEMR(V,I) to compute D1 and Z, do tests.
*
* Compute D1 and Z
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU )
*
IF( SRANGE ) THEN
NTEST = 29
IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) )
IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) )
IF( IU.LT.IL ) THEN
ITEMP = IU
IU = IL
IL = ITEMP
END IF
CALL SSTEMR( 'V', 'I', N, D5, WORK, VL, VU, IL, IU,
$ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I)', IINFO,
$ N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 29 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Tests 29 and 30
*
CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK,
$ M, RESULT( 29 ) )
*
* Call SSTEMR to compute D2, do tests.
*
* Compute D2
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
*
NTEST = 31
CALL SSTEMR( 'N', 'I', N, D5, WORK, VL, VU, IL, IU,
$ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,I)', IINFO,
$ N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 31 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Test 31
*
TEMP1 = ZERO
TEMP2 = ZERO
*
DO 240 J = 1, IU - IL + 1
TEMP1 = MAX( TEMP1, ABS( D1( J ) ),
$ ABS( D2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) )
240 CONTINUE
*
RESULT( 31 ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
*
* Call SSTEMR(V,V) to compute D1 and Z, do tests.
*
* Compute D1 and Z
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU )
*
NTEST = 32
*
IF( N.GT.0 ) THEN
IF( IL.NE.1 ) THEN
VL = D2( IL ) - MAX( HALF*
$ ( D2( IL )-D2( IL-1 ) ), ULP*ANORM,
$ TWO*RTUNFL )
ELSE
VL = D2( 1 ) - MAX( HALF*( D2( N )-D2( 1 ) ),
$ ULP*ANORM, TWO*RTUNFL )
END IF
IF( IU.NE.N ) THEN
VU = D2( IU ) + MAX( HALF*
$ ( D2( IU+1 )-D2( IU ) ), ULP*ANORM,
$ TWO*RTUNFL )
ELSE
VU = D2( N ) + MAX( HALF*( D2( N )-D2( 1 ) ),
$ ULP*ANORM, TWO*RTUNFL )
END IF
ELSE
VL = ZERO
VU = ONE
END IF
*
CALL SSTEMR( 'V', 'V', N, D5, WORK, VL, VU, IL, IU,
$ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,V)', IINFO,
$ N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 32 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Tests 32 and 33
*
CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK,
$ M, RESULT( 32 ) )
*
* Call SSTEMR to compute D2, do tests.
*
* Compute D2
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
*
NTEST = 34
CALL SSTEMR( 'N', 'V', N, D5, WORK, VL, VU, IL, IU,
$ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,V)', IINFO,
$ N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 34 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Test 34
*
TEMP1 = ZERO
TEMP2 = ZERO
*
DO 250 J = 1, IU - IL + 1
TEMP1 = MAX( TEMP1, ABS( D1( J ) ),
$ ABS( D2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) )
250 CONTINUE
*
RESULT( 34 ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
ELSE
RESULT( 29 ) = ZERO
RESULT( 30 ) = ZERO
RESULT( 31 ) = ZERO
RESULT( 32 ) = ZERO
RESULT( 33 ) = ZERO
RESULT( 34 ) = ZERO
END IF
*
*
* Call SSTEMR(V,A) to compute D1 and Z, do tests.
*
* Compute D1 and Z
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
*
NTEST = 35
*
CALL SSTEMR( 'V', 'A', N, D5, WORK, VL, VU, IL, IU,
$ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A)', IINFO, N,
$ JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 35 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Tests 35 and 36
*
CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, M,
$ RESULT( 35 ) )
*
* Call SSTEMR to compute D2, do tests.
*
* Compute D2
*
CALL SCOPY( N, SD, 1, D5, 1 )
IF( N.GT.0 )
$ CALL SCOPY( N-1, SE, 1, WORK, 1 )
*
NTEST = 37
CALL SSTEMR( 'N', 'A', N, D5, WORK, VL, VU, IL, IU,
$ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC,
$ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ),
$ LIWORK-2*N, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,A)', IINFO, N,
$ JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( 37 ) = ULPINV
GO TO 280
END IF
END IF
*
* Do Test 34
*
TEMP1 = ZERO
TEMP2 = ZERO
*
DO 260 J = 1, N
TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) )
260 CONTINUE
*
RESULT( 37 ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
END IF
270 CONTINUE
280 CONTINUE
NTESTT = NTESTT + NTEST
*
* End of Loop -- Check for RESULT(j) > THRESH
*
*
* Print out tests which fail.
*
DO 290 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 )'SST'
WRITE( NOUNIT, FMT = 9997 )
WRITE( NOUNIT, FMT = 9996 )
WRITE( NOUNIT, FMT = 9995 )'Symmetric'
WRITE( NOUNIT, FMT = 9994 )
*
* Tests performed
*
WRITE( NOUNIT, FMT = 9988 )
END IF
NERRS = NERRS + 1
WRITE( NOUNIT, FMT = 9990 )N, IOLDSD, JTYPE, JR,
$ RESULT( JR )
END IF
290 CONTINUE
300 CONTINUE
310 CONTINUE
*
* Summary
*
CALL SLASUM( 'SST', NOUNIT, NERRS, NTESTT )
RETURN
*
9999 FORMAT( ' SCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
$ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
*
9998 FORMAT( / 1X, A3, ' -- Real Symmetric eigenvalue problem' )
9997 FORMAT( ' Matrix types (see SCHKST 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, ' 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( ' 16=Positive definite, evenly spaced eigenvalues',
$ / ' 17=Positive definite, geometrically spaced eigenvlaues',
$ / ' 18=Positive definite, clustered eigenvalues',
$ / ' 19=Positive definite, small evenly spaced eigenvalues',
$ / ' 20=Positive definite, large evenly spaced eigenvalues',
$ / ' 21=Diagonally dominant tridiagonal, geometrically',
$ ' spaced eigenvalues' )
*
9993 FORMAT( / ' Tests performed: ',
$ '(S is Tridiag, D is diagonal, U and Z are ', A, ',', / 20X,
$ A, ', W is a diagonal matrix of eigenvalues,', / 20X,
$ ' V is U represented by Householder vectors, and', / 20X,
$ ' Y is a matrix of eigenvectors of S.)',
$ / ' SSYTRD, UPLO=''U'':', / ' 1= | A - V S V', A1,
$ ' | / ( |A| n ulp ) ', ' 2= | I - U V', A1,
$ ' | / ( n ulp )', / ' SSYTRD, UPLO=''L'':',
$ / ' 3= | A - V S V', A1, ' | / ( |A| n ulp ) ',
$ ' 4= | I - U V', A1, ' | / ( n ulp )' )
9992 FORMAT( ' SSPTRD, UPLO=''U'':', / ' 5= | A - V S V', A1,
$ ' | / ( |A| n ulp ) ', ' 6= | I - U V', A1,
$ ' | / ( n ulp )', / ' SSPTRD, UPLO=''L'':',
$ / ' 7= | A - V S V', A1, ' | / ( |A| n ulp ) ',
$ ' 8= | I - U V', A1, ' | / ( n ulp )',
$ / ' 9= | S - Z D Z', A1, ' | / ( |S| n ulp ) ',
$ ' 10= | I - Z Z', A1, ' | / ( n ulp )',
$ / ' 11= |D(with Z) - D(w/o Z)| / (|D| ulp) ',
$ ' 12= | D(PWK) - D(QR) | / (|D| ulp)',
$ / ' 13= Sturm sequence test on W ' )
9991 FORMA⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?