sdrvst.f
来自「famous linear algebra library (LAPACK) p」· F 代码 · 共 1,699 行 · 第 1/5 页
F
1,699 行
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
RESULT( NTEST+1 ) = ULPINV
RESULT( NTEST+2 ) = ULPINV
GO TO 660
END IF
END IF
*
* Do tests 25 and 26 (or +54)
*
CALL SSYT21( 1, UPLO, N, 0, V, LDU, D1, D2, A, LDU, Z,
$ LDU, TAU, WORK, RESULT( NTEST ) )
*
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
NTEST = NTEST + 2
SRNAMT = 'SSYEV'
CALL SSYEV( 'N', UPLO, N, A, LDU, D3, WORK, LWORK,
$ IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEV(N,' // UPLO // ')',
$ IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 660
END IF
END IF
*
* Do test 27 (or +54)
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 650 J = 1, N
TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D3( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D3( J ) ) )
650 CONTINUE
RESULT( NTEST ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
660 CONTINUE
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
NTEST = NTEST + 1
*
IF( N.GT.0 ) THEN
TEMP3 = MAX( ABS( D1( 1 ) ), ABS( D1( N ) ) )
IF( IL.NE.1 ) THEN
VL = D1( IL ) - MAX( HALF*( D1( IL )-D1( IL-1 ) ),
$ TEN*ULP*TEMP3, TEN*RTUNFL )
ELSE IF( N.GT.0 ) THEN
VL = D1( 1 ) - MAX( HALF*( D1( N )-D1( 1 ) ),
$ TEN*ULP*TEMP3, TEN*RTUNFL )
END IF
IF( IU.NE.N ) THEN
VU = D1( IU ) + MAX( HALF*( D1( IU+1 )-D1( IU ) ),
$ TEN*ULP*TEMP3, TEN*RTUNFL )
ELSE IF( N.GT.0 ) THEN
VU = D1( N ) + MAX( HALF*( D1( N )-D1( 1 ) ),
$ TEN*ULP*TEMP3, TEN*RTUNFL )
END IF
ELSE
TEMP3 = ZERO
VL = ZERO
VU = ONE
END IF
*
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'V', 'A', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M, WA1, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(V,A,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
RESULT( NTEST+1 ) = ULPINV
RESULT( NTEST+2 ) = ULPINV
GO TO 680
END IF
END IF
*
* Do tests 28 and 29 (or +54)
*
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
CALL SSYT21( 1, UPLO, N, 0, A, LDU, D1, D2, Z, LDU, V,
$ LDU, TAU, WORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'N', 'A', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(N,A,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 680
END IF
END IF
*
* Do test 30 (or +54)
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 670 J = 1, N
TEMP1 = MAX( TEMP1, ABS( WA1( J ) ), ABS( WA2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( WA1( J )-WA2( J ) ) )
670 CONTINUE
RESULT( NTEST ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
680 CONTINUE
*
NTEST = NTEST + 1
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'V', 'I', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(V,I,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
RESULT( NTEST+1 ) = ULPINV
RESULT( NTEST+2 ) = ULPINV
GO TO 690
END IF
END IF
*
* Do tests 31 and 32 (or +54)
*
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
CALL SSYT22( 1, UPLO, N, M2, 0, A, LDU, WA2, D2, Z, LDU,
$ V, LDU, TAU, WORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'N', 'I', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M3, WA3, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(N,I,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 690
END IF
END IF
*
* Do test 33 (or +54)
*
TEMP1 = SSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL )
TEMP2 = SSXT1( 1, WA3, M3, WA2, M2, ABSTOL, ULP, UNFL )
RESULT( NTEST ) = ( TEMP1+TEMP2 ) /
$ MAX( UNFL, ULP*TEMP3 )
690 CONTINUE
*
NTEST = NTEST + 1
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'V', 'V', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(V,V,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
RESULT( NTEST+1 ) = ULPINV
RESULT( NTEST+2 ) = ULPINV
GO TO 700
END IF
END IF
*
* Do tests 34 and 35 (or +54)
*
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
CALL SSYT22( 1, UPLO, N, M2, 0, A, LDU, WA2, D2, Z, LDU,
$ V, LDU, TAU, WORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
SRNAMT = 'SSYEVX'
CALL SSYEVX( 'N', 'V', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M3, WA3, Z, LDU, WORK, LWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSYEVX(N,V,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 700
END IF
END IF
*
IF( M3.EQ.0 .AND. N.GT.0 ) THEN
RESULT( NTEST ) = ULPINV
GO TO 700
END IF
*
* Do test 36 (or +54)
*
TEMP1 = SSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL )
TEMP2 = SSXT1( 1, WA3, M3, WA2, M2, ABSTOL, ULP, UNFL )
IF( N.GT.0 ) THEN
TEMP3 = MAX( ABS( WA1( 1 ) ), ABS( WA1( N ) ) )
ELSE
TEMP3 = ZERO
END IF
RESULT( NTEST ) = ( TEMP1+TEMP2 ) /
$ MAX( UNFL, TEMP3*ULP )
*
700 CONTINUE
*
* 5) Call SSPEV and SSPEVX.
*
CALL SLACPY( ' ', N, N, V, LDU, A, LDA )
*
* Load array WORK with the upper or lower triangular
* part of the matrix in packed form.
*
IF( IUPLO.EQ.1 ) THEN
INDX = 1
DO 720 J = 1, N
DO 710 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
710 CONTINUE
720 CONTINUE
ELSE
INDX = 1
DO 740 J = 1, N
DO 730 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
730 CONTINUE
740 CONTINUE
END IF
*
NTEST = NTEST + 1
SRNAMT = 'SSPEV'
CALL SSPEV( 'V', UPLO, N, WORK, D1, Z, LDU, V, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSPEV(V,' // UPLO // ')',
$ IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
RESULT( NTEST+1 ) = ULPINV
RESULT( NTEST+2 ) = ULPINV
GO TO 800
END IF
END IF
*
* Do tests 37 and 38 (or +54)
*
CALL SSYT21( 1, UPLO, N, 0, A, LDA, D1, D2, Z, LDU, V,
$ LDU, TAU, WORK, RESULT( NTEST ) )
*
IF( IUPLO.EQ.1 ) THEN
INDX = 1
DO 760 J = 1, N
DO 750 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
750 CONTINUE
760 CONTINUE
ELSE
INDX = 1
DO 780 J = 1, N
DO 770 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
770 CONTINUE
780 CONTINUE
END IF
*
NTEST = NTEST + 2
SRNAMT = 'SSPEV'
CALL SSPEV( 'N', UPLO, N, WORK, D3, Z, LDU, V, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'SSPEV(N,' // UPLO // ')',
$ IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 800
END IF
END IF
*
* Do test 39 (or +54)
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 790 J = 1, N
TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D3( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D3( J ) ) )
790 CONTINUE
RESULT( NTEST ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
* Load array WORK with the upper or lower triangular part
* of the matrix in packed form.
*
800 CONTINUE
IF(⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?