diff --git a/SRC/cbdsqr.f b/SRC/cbdsqr.f index fd6d4731b..df0ee5076 100644 --- a/SRC/cbdsqr.f +++ b/SRC/cbdsqr.f @@ -455,7 +455,11 @@ SUBROUTINE CBDSQR( UPLO, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, ABSE = ABS( E( LL ) ) IF( TOL.LT.ZERO .AND. ABSS.LE.THRESH ) $ D( LL ) = ZERO - IF( ABSE.LE.THRESH ) +* An exact zero always marks a split, also when THRESH is NaN +* because the input contains an Inf; otherwise the convergence +* tests below zero E( M-1 ) and return here without the split +* ever being found. + IF( ABSE.LE.THRESH .OR. ABSE.EQ.ZERO ) $ GO TO 80 SMAX = MAX( SMAX, ABSS, ABSE ) 70 CONTINUE diff --git a/SRC/dbdsqr.f b/SRC/dbdsqr.f index 99a12aa79..e9dd2ae2b 100644 --- a/SRC/dbdsqr.f +++ b/SRC/dbdsqr.f @@ -462,7 +462,11 @@ SUBROUTINE DBDSQR( UPLO, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, ABSE = ABS( E( LL ) ) IF( TOL.LT.ZERO .AND. ABSS.LE.THRESH ) $ D( LL ) = ZERO - IF( ABSE.LE.THRESH ) +* An exact zero always marks a split, also when THRESH is NaN +* because the input contains an Inf; otherwise the convergence +* tests below zero E( M-1 ) and return here without the split +* ever being found. + IF( ABSE.LE.THRESH .OR. ABSE.EQ.ZERO ) $ GO TO 80 SMAX = MAX( SMAX, ABSS, ABSE ) 70 CONTINUE diff --git a/SRC/sbdsqr.f b/SRC/sbdsqr.f index 866b13295..9c06037a1 100644 --- a/SRC/sbdsqr.f +++ b/SRC/sbdsqr.f @@ -462,7 +462,11 @@ SUBROUTINE SBDSQR( UPLO, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, ABSE = ABS( E( LL ) ) IF( TOL.LT.ZERO .AND. ABSS.LE.THRESH ) $ D( LL ) = ZERO - IF( ABSE.LE.THRESH ) +* An exact zero always marks a split, also when THRESH is NaN +* because the input contains an Inf; otherwise the convergence +* tests below zero E( M-1 ) and return here without the split +* ever being found. + IF( ABSE.LE.THRESH .OR. ABSE.EQ.ZERO ) $ GO TO 80 SMAX = MAX( SMAX, ABSS, ABSE ) 70 CONTINUE diff --git a/SRC/zbdsqr.f b/SRC/zbdsqr.f index ec67efeb0..74eb8776a 100644 --- a/SRC/zbdsqr.f +++ b/SRC/zbdsqr.f @@ -453,7 +453,11 @@ SUBROUTINE ZBDSQR( UPLO, N, NCVT, NRU, NCC, D, E, VT, LDVT, U, ABSE = ABS( E( LL ) ) IF( TOL.LT.ZERO .AND. ABSS.LE.THRESH ) $ D( LL ) = ZERO - IF( ABSE.LE.THRESH ) +* An exact zero always marks a split, also when THRESH is NaN +* because the input contains an Inf; otherwise the convergence +* tests below zero E( M-1 ) and return here without the split +* ever being found. + IF( ABSE.LE.THRESH .OR. ABSE.EQ.ZERO ) $ GO TO 80 SMAX = MAX( SMAX, ABSS, ABSE ) 70 CONTINUE diff --git a/TESTING/EIG/cerrbd.f b/TESTING/EIG/cerrbd.f index 71f49e80a..ed4fa15db 100644 --- a/TESTING/EIG/cerrbd.f +++ b/TESTING/EIG/cerrbd.f @@ -68,12 +68,13 @@ SUBROUTINE CERRBD( PATH, NUNIT ) * .. Parameters .. INTEGER NMAX, LW PARAMETER ( NMAX = 4, LW = NMAX ) - REAL ONE - PARAMETER ( ONE = 1.0E+0 ) + REAL ZERO, ONE + PARAMETER ( ZERO = 0.0E+0, ONE = 1.0E+0 ) * .. * .. Local Scalars .. CHARACTER*2 C2 INTEGER I, INFO, J, NT + REAL RZERO * .. * .. Local Arrays .. REAL D( NMAX ), E( NMAX ), RW( 4*NMAX ) @@ -279,6 +280,34 @@ SUBROUTINE CERRBD( PATH, NUNIT ) $ INFO ) CALL CHKXER( 'CBDSQR', INFOT, NOUT, LERR, OK ) NT = NT + 8 +* +* CBDSQR with singular vectors must return when D contains an +* infinity, which makes its convergence threshold a NaN, +* instead of iterating forever. +* + RZERO = ZERO + DO 40 J = 1, NMAX + D( J ) = REAL( J ) + E( J ) = ONE / REAL( J+1 ) + DO 30 I = 1, NMAX + U( I, J ) = ZERO + V( I, J ) = ZERO + 30 CONTINUE + U( J, J ) = ONE + V( J, J ) = ONE + 40 CONTINUE + D( 1 ) = ONE / RZERO + E( NMAX ) = ZERO + SRNAMT = 'CBDSQR' + INFOT = 0 + LERR = .FALSE. + CALL CBDSQR( 'U', NMAX, NMAX, NMAX, 0, D, E, V, NMAX, U, + $ NMAX, A, NMAX, RW, INFO ) + IF( LERR ) THEN + WRITE( NOUT, FMT = 9997 )'CBDSQR' + OK = .FALSE. + END IF + NT = NT + 1 END IF * * Print a summary line. @@ -293,6 +322,8 @@ SUBROUTINE CERRBD( PATH, NUNIT ) $ I3, ' tests done)' ) 9998 FORMAT( ' *** ', A3, ' routines failed the tests of the error ', $ 'exits ***' ) + 9997 FORMAT( ' *** ', A6, ' called XERBLA for a matrix with an ', + $ 'infinity ***' ) * RETURN * diff --git a/TESTING/EIG/derrbd.f b/TESTING/EIG/derrbd.f index 8f2b6603b..af8a21b23 100644 --- a/TESTING/EIG/derrbd.f +++ b/TESTING/EIG/derrbd.f @@ -67,13 +67,14 @@ SUBROUTINE DERRBD( PATH, NUNIT ) * * .. Parameters .. INTEGER NMAX, LW - PARAMETER ( NMAX = 4, LW = NMAX ) + PARAMETER ( NMAX = 4, LW = 4*NMAX ) DOUBLE PRECISION ZERO, ONE PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) * .. * .. Local Scalars .. CHARACTER*2 C2 INTEGER I, INFO, J, NS, NT + DOUBLE PRECISION RZERO * .. * .. Local Arrays .. INTEGER IQ( NMAX, NMAX ), IW( NMAX ) @@ -278,6 +279,34 @@ SUBROUTINE DERRBD( PATH, NUNIT ) CALL CHKXER( 'DBDSQR', INFOT, NOUT, LERR, OK ) NT = NT + 8 * +* DBDSQR with singular vectors must return when D contains an +* infinity, which makes its convergence threshold a NaN, +* instead of iterating forever. +* + RZERO = ZERO + DO 40 J = 1, NMAX + D( J ) = DBLE( J ) + E( J ) = ONE / DBLE( J+1 ) + DO 30 I = 1, NMAX + U( I, J ) = ZERO + V( I, J ) = ZERO + 30 CONTINUE + U( J, J ) = ONE + V( J, J ) = ONE + 40 CONTINUE + D( 1 ) = ONE / RZERO + E( NMAX ) = ZERO + SRNAMT = 'DBDSQR' + INFOT = 0 + LERR = .FALSE. + CALL DBDSQR( 'U', NMAX, NMAX, NMAX, 0, D, E, V, NMAX, U, + $ NMAX, A, NMAX, W, INFO ) + IF( LERR ) THEN + WRITE( NOUT, FMT = 9997 )'DBDSQR' + OK = .FALSE. + END IF + NT = NT + 1 +* * DBDSDC * SRNAMT = 'DBDSDC' @@ -369,6 +398,8 @@ SUBROUTINE DERRBD( PATH, NUNIT ) $ ' (', I3, ' tests done)' ) 9998 FORMAT( ' *** ', A3, ' routines failed the tests of the error ', $ 'exits ***' ) + 9997 FORMAT( ' *** ', A6, ' called XERBLA for a matrix with an ', + $ 'infinity ***' ) * RETURN * diff --git a/TESTING/EIG/serrbd.f b/TESTING/EIG/serrbd.f index 683b8952e..0d0202e83 100644 --- a/TESTING/EIG/serrbd.f +++ b/TESTING/EIG/serrbd.f @@ -67,13 +67,14 @@ SUBROUTINE SERRBD( PATH, NUNIT ) * * .. Parameters .. INTEGER NMAX, LW - PARAMETER ( NMAX = 4, LW = NMAX ) + PARAMETER ( NMAX = 4, LW = 4*NMAX ) REAL ZERO, ONE PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 ) * .. * .. Local Scalars .. CHARACTER*2 C2 INTEGER I, INFO, J, NS, NT + REAL RZERO * .. * .. Local Arrays .. INTEGER IQ( NMAX, NMAX ), IW( NMAX ) @@ -278,6 +279,34 @@ SUBROUTINE SERRBD( PATH, NUNIT ) CALL CHKXER( 'SBDSQR', INFOT, NOUT, LERR, OK ) NT = NT + 8 * +* SBDSQR with singular vectors must return when D contains an +* infinity, which makes its convergence threshold a NaN, +* instead of iterating forever. +* + RZERO = ZERO + DO 40 J = 1, NMAX + D( J ) = REAL( J ) + E( J ) = ONE / REAL( J+1 ) + DO 30 I = 1, NMAX + U( I, J ) = ZERO + V( I, J ) = ZERO + 30 CONTINUE + U( J, J ) = ONE + V( J, J ) = ONE + 40 CONTINUE + D( 1 ) = ONE / RZERO + E( NMAX ) = ZERO + SRNAMT = 'SBDSQR' + INFOT = 0 + LERR = .FALSE. + CALL SBDSQR( 'U', NMAX, NMAX, NMAX, 0, D, E, V, NMAX, U, + $ NMAX, A, NMAX, W, INFO ) + IF( LERR ) THEN + WRITE( NOUT, FMT = 9997 )'SBDSQR' + OK = .FALSE. + END IF + NT = NT + 1 +* * SBDSDC * SRNAMT = 'SBDSDC' @@ -369,6 +398,8 @@ SUBROUTINE SERRBD( PATH, NUNIT ) $ ' (', I3, ' tests done)' ) 9998 FORMAT( ' *** ', A3, ' routines failed the tests of the error ', $ 'exits ***' ) + 9997 FORMAT( ' *** ', A6, ' called XERBLA for a matrix with an ', + $ 'infinity ***' ) * RETURN * diff --git a/TESTING/EIG/zerrbd.f b/TESTING/EIG/zerrbd.f index 08e7821f5..d18654d12 100644 --- a/TESTING/EIG/zerrbd.f +++ b/TESTING/EIG/zerrbd.f @@ -68,12 +68,13 @@ SUBROUTINE ZERRBD( PATH, NUNIT ) * .. Parameters .. INTEGER NMAX, LW PARAMETER ( NMAX = 4, LW = NMAX ) - DOUBLE PRECISION ONE - PARAMETER ( ONE = 1.0D+0 ) + DOUBLE PRECISION ZERO, ONE + PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) * .. * .. Local Scalars .. CHARACTER*2 C2 INTEGER I, INFO, J, NT + DOUBLE PRECISION RZERO * .. * .. Local Arrays .. DOUBLE PRECISION D( NMAX ), E( NMAX ), RW( 4*NMAX ) @@ -279,6 +280,34 @@ SUBROUTINE ZERRBD( PATH, NUNIT ) $ INFO ) CALL CHKXER( 'ZBDSQR', INFOT, NOUT, LERR, OK ) NT = NT + 8 +* +* ZBDSQR with singular vectors must return when D contains an +* infinity, which makes its convergence threshold a NaN, +* instead of iterating forever. +* + RZERO = ZERO + DO 40 J = 1, NMAX + D( J ) = DBLE( J ) + E( J ) = ONE / DBLE( J+1 ) + DO 30 I = 1, NMAX + U( I, J ) = ZERO + V( I, J ) = ZERO + 30 CONTINUE + U( J, J ) = ONE + V( J, J ) = ONE + 40 CONTINUE + D( 1 ) = ONE / RZERO + E( NMAX ) = ZERO + SRNAMT = 'ZBDSQR' + INFOT = 0 + LERR = .FALSE. + CALL ZBDSQR( 'U', NMAX, NMAX, NMAX, 0, D, E, V, NMAX, U, + $ NMAX, A, NMAX, RW, INFO ) + IF( LERR ) THEN + WRITE( NOUT, FMT = 9997 )'ZBDSQR' + OK = .FALSE. + END IF + NT = NT + 1 END IF * * Print a summary line. @@ -293,6 +322,8 @@ SUBROUTINE ZERRBD( PATH, NUNIT ) $ I3, ' tests done)' ) 9998 FORMAT( ' *** ', A3, ' routines failed the tests of the error ', $ 'exits ***' ) + 9997 FORMAT( ' *** ', A6, ' called XERBLA for a matrix with an ', + $ 'infinity ***' ) * RETURN *