zdrvst.f
来自「famous linear algebra library (LAPACK) p」· F 代码 · 共 1,778 行 · 第 1/5 页
F
1,778 行
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, RWORK,
$ IWORK, IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHEEVX(N,A,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 150
END IF
END IF
*
* Do test 6.
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 140 J = 1, N
TEMP1 = MAX( TEMP1, ABS( WA1( J ) ), ABS( WA2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( WA1( J )-WA2( J ) ) )
140 CONTINUE
RESULT( NTEST ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
150 CONTINUE
CALL ZLACPY( ' ', N, N, V, LDU, A, LDA )
*
NTEST = NTEST + 1
*
CALL ZHEEVX( 'V', 'I', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, RWORK,
$ IWORK, IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHEEVX(V,I,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 160
END IF
END IF
*
* Do tests 7 and 8.
*
CALL ZLACPY( ' ', N, N, V, LDU, A, LDA )
*
CALL ZHET22( 1, UPLO, N, M2, 0, A, LDU, WA2, D2, Z, LDU,
$ V, LDU, TAU, WORK, RWORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
*
CALL ZHEEVX( 'N', 'I', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M3, WA3, Z, LDU, WORK, LWORK, RWORK,
$ IWORK, IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHEEVX(N,I,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 160
END IF
END IF
*
* Do test 9.
*
TEMP1 = DSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL )
TEMP2 = DSXT1( 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 )
*
160 CONTINUE
CALL ZLACPY( ' ', N, N, V, LDU, A, LDA )
*
NTEST = NTEST + 1
*
CALL ZHEEVX( 'V', 'V', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, WORK, LWORK, RWORK,
$ IWORK, IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHEEVX(V,V,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 170
END IF
END IF
*
* Do tests 10 and 11.
*
CALL ZLACPY( ' ', N, N, V, LDU, A, LDA )
*
CALL ZHET22( 1, UPLO, N, M2, 0, A, LDU, WA2, D2, Z, LDU,
$ V, LDU, TAU, WORK, RWORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
*
CALL ZHEEVX( 'N', 'V', UPLO, N, A, LDU, VL, VU, IL, IU,
$ ABSTOL, M3, WA3, Z, LDU, WORK, LWORK, RWORK,
$ IWORK, IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHEEVX(N,V,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 170
END IF
END IF
*
IF( M3.EQ.0 .AND. N.GT.0 ) THEN
RESULT( NTEST ) = ULPINV
GO TO 170
END IF
*
* Do test 12.
*
TEMP1 = DSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL )
TEMP2 = DSXT1( 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 )
*
170 CONTINUE
*
* Call ZHPEVD and CHPEVX.
*
CALL ZLACPY( ' ', 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 190 J = 1, N
DO 180 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
180 CONTINUE
190 CONTINUE
ELSE
INDX = 1
DO 210 J = 1, N
DO 200 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
200 CONTINUE
210 CONTINUE
END IF
*
NTEST = NTEST + 1
INDWRK = N*( N+1 ) / 2 + 1
CALL ZHPEVD( 'V', UPLO, N, WORK, D1, Z, LDU,
$ WORK( INDWRK ), LWEDC, RWORK, LRWEDC, IWORK,
$ LIWEDC, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHPEVD(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 270
END IF
END IF
*
* Do tests 13 and 14.
*
CALL ZHET21( 1, UPLO, N, 0, A, LDA, D1, D2, Z, LDU, V,
$ LDU, TAU, WORK, RWORK, RESULT( NTEST ) )
*
IF( IUPLO.EQ.1 ) THEN
INDX = 1
DO 230 J = 1, N
DO 220 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
220 CONTINUE
230 CONTINUE
ELSE
INDX = 1
DO 250 J = 1, N
DO 240 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
240 CONTINUE
250 CONTINUE
END IF
*
NTEST = NTEST + 2
INDWRK = N*( N+1 ) / 2 + 1
CALL ZHPEVD( 'N', UPLO, N, WORK, D3, Z, LDU,
$ WORK( INDWRK ), LWEDC, RWORK, LRWEDC, IWORK,
$ LIWEDC, IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHPEVD(N,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 270
END IF
END IF
*
* Do test 15.
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 260 J = 1, N
TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D3( J ) ) )
TEMP2 = MAX( TEMP2, ABS( D1( J )-D3( J ) ) )
260 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.
*
270 CONTINUE
IF( IUPLO.EQ.1 ) THEN
INDX = 1
DO 290 J = 1, N
DO 280 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
280 CONTINUE
290 CONTINUE
ELSE
INDX = 1
DO 310 J = 1, N
DO 300 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
300 CONTINUE
310 CONTINUE
END IF
*
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
*
CALL ZHPEVX( 'V', 'A', UPLO, N, WORK, VL, VU, IL, IU,
$ ABSTOL, M, WA1, Z, LDU, V, RWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHPEVX(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 370
END IF
END IF
*
* Do tests 16 and 17.
*
CALL ZHET21( 1, UPLO, N, 0, A, LDU, WA1, D2, Z, LDU, V,
$ LDU, TAU, WORK, RWORK, RESULT( NTEST ) )
*
NTEST = NTEST + 2
*
IF( IUPLO.EQ.1 ) THEN
INDX = 1
DO 330 J = 1, N
DO 320 I = 1, J
WORK( INDX ) = A( I, J )
INDX = INDX + 1
320 CONTINUE
330 CONTINUE
ELSE
INDX = 1
DO 350 J = 1, N
DO 340 I = J, N
WORK( INDX ) = A( I, J )
INDX = INDX + 1
340 CONTINUE
350 CONTINUE
END IF
*
CALL ZHPEVX( 'N', 'A', UPLO, N, WORK, VL, VU, IL, IU,
$ ABSTOL, M2, WA2, Z, LDU, V, RWORK, IWORK,
$ IWORK( 5*N+1 ), IINFO )
IF( IINFO.NE.0 ) THEN
WRITE( NOUNIT, FMT = 9999 )'ZHPEVX(N,A,' // UPLO //
$ ')', IINFO, N, JTYPE, IOLDSD
INFO = ABS( IINFO )
IF( IINFO.LT.0 ) THEN
RETURN
ELSE
RESULT( NTEST ) = ULPINV
GO TO 370
END IF
END IF
*
* Do test 18.
*
TEMP1 = ZERO
TEMP2 = ZERO
DO 360 J = 1, N
TEMP1 = MAX( TEMP1, ABS( WA1( J ) ), ABS( WA2( J ) ) )
TEMP2 = MAX( TEMP2, ABS( WA1( J )-WA2( J ) ) )
360 CONTINUE
RESULT( NTEST ) = TEMP2 / MAX( UNFL,
$ ULP*MAX( TEMP1, TEMP2 ) )
*
370 CONTINUE
NTEST = NTEST + 1
IF( IUPLO.EQ.1 ) THEN
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?