Skip to content

Keep a NaN inside its block in xSTEDC and report it through INFO #1403

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
rmlarsen wants to merge 2 commits into
base: master
Choose a base branch
from
Open

Keep a NaN inside its block in xSTEDC and report it through INFO #1403

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
b9ed114
Keep a NaN inside its block in xSTEDC and report it through INFO
rmlarsen Sep 10, 2026
8cca9b5
TESTING: Bound the STEDC NaN regression storage
rmlarsen Sep 11, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
13 changes: 10 additions & 3 deletions SRC/cstedc.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -230,10 +230,10 @@ SUBROUTINE CSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
REAL EPS, ORGNRM, P, TINY
* ..
* .. External Functions ..
LOGICAL LSAME
LOGICAL SISNAN, LSAME
INTEGER ILAENV
REAL SLAMCH, SLANST, SROUNDUP_LWORK
EXTERNAL ILAENV, LSAME, SLAMCH, SLANST,
EXTERNAL ILAENV, LSAME, SISNAN, SLAMCH, SLANST,
$ SROUNDUP_LWORK
* ..
* .. External Subroutines ..
Expand Down Expand Up @@ -395,7 +395,10 @@ SUBROUTINE CSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
IF( FINISH.LT.N ) THEN
TINY = EPS*SQRT( ABS( D( FINISH ) ) )*
$ SQRT( ABS( D( FINISH+1 ) ) )
IF( ABS( E( FINISH ) ).GT.TINY ) THEN
* A NaN in D or E must not split the matrix: keep it in the
* block so that it is reported below instead of being
* dropped or isolated.
IF( .NOT.( ABS( E( FINISH ) ).LE.TINY ) ) THEN
FINISH = FINISH + 1
GO TO 40
END IF
Expand All @@ -409,6 +412,10 @@ SUBROUTINE CSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
* Scale.
*
ORGNRM = SLANST( 'M', M, D( START ), E( START ) )
IF( SISNAN( ORGNRM ) ) THEN
INFO = START*( N+1 ) + FINISH
GO TO 70
END IF
CALL SLASCL( 'G', 0, 0, ORGNRM, ONE, M, 1, D( START ),
$ M,
$ INFO )
Expand Down
13 changes: 10 additions & 3 deletions SRC/dstedc.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -205,10 +205,10 @@ SUBROUTINE DSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
DOUBLE PRECISION EPS, ORGNRM, P, TINY
* ..
* .. External Functions ..
LOGICAL LSAME
LOGICAL DISNAN, LSAME
INTEGER ILAENV
DOUBLE PRECISION DLAMCH, DLANST
EXTERNAL LSAME, ILAENV, DLAMCH, DLANST
EXTERNAL DISNAN, LSAME, ILAENV, DLAMCH, DLANST
* ..
* .. External Subroutines ..
EXTERNAL DGEMM, DLACPY, DLAED0, DLASCL, DLASET,
Expand Down Expand Up @@ -359,7 +359,10 @@ SUBROUTINE DSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
IF( FINISH.LT.N ) THEN
TINY = EPS*SQRT( ABS( D( FINISH ) ) )*
$ SQRT( ABS( D( FINISH+1 ) ) )
IF( ABS( E( FINISH ) ).GT.TINY ) THEN
* A NaN in D or E must not split the matrix: keep it in the
* block so that it is reported below instead of being
* dropped or isolated.
IF( .NOT.( ABS( E( FINISH ) ).LE.TINY ) ) THEN
FINISH = FINISH + 1
GO TO 20
END IF
Expand All @@ -377,6 +380,10 @@ SUBROUTINE DSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
* Scale.
*
ORGNRM = DLANST( 'M', M, D( START ), E( START ) )
IF( DISNAN( ORGNRM ) ) THEN
INFO = START*( N+1 ) + FINISH
GO TO 50
END IF
CALL DLASCL( 'G', 0, 0, ORGNRM, ONE, M, 1, D( START ),
$ M,
$ INFO )
Expand Down
13 changes: 10 additions & 3 deletions SRC/sstedc.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -205,10 +205,10 @@ SUBROUTINE SSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
REAL EPS, ORGNRM, P, TINY
* ..
* .. External Functions ..
LOGICAL LSAME
LOGICAL SISNAN, LSAME
INTEGER ILAENV
REAL SLAMCH, SLANST, SROUNDUP_LWORK
EXTERNAL ILAENV, LSAME, SLAMCH, SLANST,
EXTERNAL ILAENV, LSAME, SISNAN, SLAMCH, SLANST,
$ SROUNDUP_LWORK
* ..
* .. External Subroutines ..
Expand Down Expand Up @@ -360,7 +360,10 @@ SUBROUTINE SSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
IF( FINISH.LT.N ) THEN
TINY = EPS*SQRT( ABS( D( FINISH ) ) )*
$ SQRT( ABS( D( FINISH+1 ) ) )
IF( ABS( E( FINISH ) ).GT.TINY ) THEN
* A NaN in D or E must not split the matrix: keep it in the
* block so that it is reported below instead of being
* dropped or isolated.
IF( .NOT.( ABS( E( FINISH ) ).LE.TINY ) ) THEN
FINISH = FINISH + 1
GO TO 20
END IF
Expand All @@ -378,6 +381,10 @@ SUBROUTINE SSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
* Scale.
*
ORGNRM = SLANST( 'M', M, D( START ), E( START ) )
IF( SISNAN( ORGNRM ) ) THEN
INFO = START*( N+1 ) + FINISH
GO TO 50
END IF
CALL SLASCL( 'G', 0, 0, ORGNRM, ONE, M, 1, D( START ),
$ M,
$ INFO )
Expand Down
13 changes: 10 additions & 3 deletions SRC/zstedc.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -230,10 +230,10 @@ SUBROUTINE ZSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
DOUBLE PRECISION EPS, ORGNRM, P, TINY
* ..
* .. External Functions ..
LOGICAL LSAME
LOGICAL DISNAN, LSAME
INTEGER ILAENV
DOUBLE PRECISION DLAMCH, DLANST
EXTERNAL LSAME, ILAENV, DLAMCH, DLANST
EXTERNAL DISNAN, LSAME, ILAENV, DLAMCH, DLANST
* ..
* .. External Subroutines ..
EXTERNAL DLASCL, DLASET, DSTEDC, DSTEQR, DSTERF,
Expand Down Expand Up @@ -394,7 +394,10 @@ SUBROUTINE ZSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
IF( FINISH.LT.N ) THEN
TINY = EPS*SQRT( ABS( D( FINISH ) ) )*
$ SQRT( ABS( D( FINISH+1 ) ) )
IF( ABS( E( FINISH ) ).GT.TINY ) THEN
* A NaN in D or E must not split the matrix: keep it in the
* block so that it is reported below instead of being
* dropped or isolated.
IF( .NOT.( ABS( E( FINISH ) ).LE.TINY ) ) THEN
FINISH = FINISH + 1
GO TO 40
END IF
Expand All @@ -408,6 +411,10 @@ SUBROUTINE ZSTEDC( COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK,
* Scale.
*
ORGNRM = DLANST( 'M', M, D( START ), E( START ) )
IF( DISNAN( ORGNRM ) ) THEN
INFO = START*( N+1 ) + FINISH
GO TO 70
END IF
CALL DLASCL( 'G', 0, 0, ORGNRM, ONE, M, 1, D( START ),
$ M,
$ INFO )
Expand Down
59 changes: 59 additions & 0 deletions TESTING/EIG/cchkst.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -651,11 +651,16 @@ SUBROUTINE CCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
$ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP,
$ ULPINV, UNFL, VL, VU
* ..
CHARACTER COMPZ
INTEGER JCOMPZ, JNAN, NSMLSZ
REAL RNAN, RONE
* .. Local Arrays ..
INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
$ KMAGN( MAXTYP ), KMODE( MAXTYP ),
$ KTYPE( MAXTYP )
REAL DUMMA( 1 )
REAL, ALLOCATABLE :: DREG( : ), EREG( : )
COMPLEX, ALLOCATABLE :: ZREG( :, : )
* ..
* .. External Functions ..
INTEGER ILAENV
Expand Down Expand Up @@ -1953,11 +1958,65 @@ SUBROUTINE CCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
300 CONTINUE
310 CONTINUE
*
*
* CSTEDC must report a NaN in the matrix through INFO. A NaN
* off-diagonal entry used to split the matrix and was dropped, and
* a NaN diagonal entry was isolated as a 1 by 1 block, both with
* INFO = 0, whenever the block left over was large enough for the
* divide and conquer recursion, which the sizes in the input file
* do not reach. Use local arrays because LDU is a row stride,
* not the capacity of the caller's eigenvalue arrays or Z columns.
*
NSMLSZ = ILAENV( 9, 'CSTEDC', ' ', 0, 0, 0, 0 )
IF( ILAENV( 10, 'CSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 .AND.
$ ILAENV( 11, 'CSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 ) THEN
N = MAX( 2, NSMLSZ+1 )
ALLOCATE( DREG( N ), EREG( N ), ZREG( N, N ) )
RONE = ONE
RNAN = SQRT( -RONE )
CALL CSTEDC( 'V', N, DREG, EREG, ZREG, N, WORK, -1, RWORK, -1,
$ IWORK, -1, IINFO )
IF( IINFO.EQ.0 .AND. INT( REAL( WORK( 1 ) ) ).LE.LWORK .AND.
$ INT( RWORK( 1 ) ).LE.LRWORK .AND. IWORK( 1 ).LE.LIWORK )
$ THEN
DO 340 JNAN = 1, 2
DO 330 JCOMPZ = 1, 2
IF( JCOMPZ.EQ.1 ) THEN
COMPZ = 'I'
ELSE
COMPZ = 'V'
END IF
DO 320 J = 1, N
DREG( J ) = REAL( J )
EREG( J ) = ONE / REAL( J+1 )
320 CONTINUE
IF( JNAN.EQ.1 ) THEN
EREG( N / 2 ) = RNAN
ELSE
DREG( N ) = RNAN
END IF
CALL CLASET( 'Full', N, N, CZERO, CONE, ZREG, N )
CALL CSTEDC( COMPZ, N, DREG, EREG, ZREG, N,
$ WORK, LWORK, RWORK, LRWORK, IWORK,
$ LIWORK, IINFO )
IF( IINFO.EQ.0 ) THEN
WRITE( NOUNIT, FMT = 9985 )COMPZ, JNAN
NERRS = NERRS + 1
END IF
NTESTT = NTESTT + 1
330 CONTINUE
340 CONTINUE
END IF
DEALLOCATE( DREG, EREG, ZREG )
END IF
*
* Summary
*
CALL SLASUM( 'CST', NOUNIT, NERRS, NTESTT )
RETURN
*
9985 FORMAT( ' CCHKST: CSTEDC( ', A1, ' ) returned INFO=0 for a',
$ ' matrix with a NaN, case ', I1 )
9999 FORMAT( ' CCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
$ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
*
Expand Down
57 changes: 57 additions & 0 deletions TESTING/EIG/dchkst.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -634,11 +634,16 @@ SUBROUTINE DCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
$ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP,
$ ULPINV, UNFL, VL, VU
* ..
CHARACTER COMPZ
INTEGER JCOMPZ, JNAN, NSMLSZ
DOUBLE PRECISION RNAN, RONE
* .. Local Arrays ..
INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
$ KMAGN( MAXTYP ), KMODE( MAXTYP ),
$ KTYPE( MAXTYP )
DOUBLE PRECISION DUMMA( 1 )
DOUBLE PRECISION, ALLOCATABLE :: DREG( : ), EREG( : )
DOUBLE PRECISION, ALLOCATABLE :: ZREG( :, : )
* ..
* .. External Functions ..
INTEGER ILAENV
Expand Down Expand Up @@ -1931,11 +1936,63 @@ SUBROUTINE DCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
300 CONTINUE
310 CONTINUE
*
*
* DSTEDC must report a NaN in the matrix through INFO. A NaN
* off-diagonal entry used to split the matrix and was dropped, and
* a NaN diagonal entry was isolated as a 1 by 1 block, both with
* INFO = 0, whenever the block left over was large enough for the
* divide and conquer recursion, which the sizes in the input file
* do not reach. Use local arrays because LDU is a row stride,
* not the capacity of the caller's eigenvalue arrays or Z columns.
*
NSMLSZ = ILAENV( 9, 'DSTEDC', ' ', 0, 0, 0, 0 )
IF( ILAENV( 10, 'DSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 .AND.
$ ILAENV( 11, 'DSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 ) THEN
N = MAX( 2, NSMLSZ+1 )
ALLOCATE( DREG( N ), EREG( N ), ZREG( N, N ) )
RONE = ONE
RNAN = SQRT( -RONE )
CALL DSTEDC( 'V', N, DREG, EREG, ZREG, N, WORK, -1, IWORK, -1,
$ IINFO )
IF( IINFO.EQ.0 .AND. INT( WORK( 1 ) ).LE.LWORK .AND.
$ IWORK( 1 ).LE.LIWORK ) THEN
DO 340 JNAN = 1, 2
DO 330 JCOMPZ = 1, 2
IF( JCOMPZ.EQ.1 ) THEN
COMPZ = 'I'
ELSE
COMPZ = 'V'
END IF
DO 320 J = 1, N
DREG( J ) = DBLE( J )
EREG( J ) = ONE / DBLE( J+1 )
320 CONTINUE
IF( JNAN.EQ.1 ) THEN
EREG( N / 2 ) = RNAN
ELSE
DREG( N ) = RNAN
END IF
CALL DLASET( 'Full', N, N, ZERO, ONE, ZREG, N )
CALL DSTEDC( COMPZ, N, DREG, EREG, ZREG, N,
$ WORK, LWORK, IWORK, LIWORK, IINFO )
IF( IINFO.EQ.0 ) THEN
WRITE( NOUNIT, FMT = 9985 )COMPZ, JNAN
NERRS = NERRS + 1
END IF
NTESTT = NTESTT + 1
330 CONTINUE
340 CONTINUE
END IF
DEALLOCATE( DREG, EREG, ZREG )
END IF
*
* Summary
*
CALL DLASUM( 'DST', NOUNIT, NERRS, NTESTT )
RETURN
*
9985 FORMAT( ' DCHKST: DSTEDC( ', A1, ' ) returned INFO=0 for a',
$ ' matrix with a NaN, case ', I1 )
9999 FORMAT( ' DCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
$ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
*
Expand Down
57 changes: 57 additions & 0 deletions TESTING/EIG/schkst.f
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -634,11 +634,16 @@ SUBROUTINE SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
$ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP,
$ ULPINV, UNFL, VL, VU
* ..
CHARACTER COMPZ
INTEGER JCOMPZ, JNAN, NSMLSZ
REAL RNAN, RONE
* .. Local Arrays ..
INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
$ KMAGN( MAXTYP ), KMODE( MAXTYP ),
$ KTYPE( MAXTYP )
REAL DUMMA( 1 )
REAL, ALLOCATABLE :: DREG( : ), EREG( : )
REAL, ALLOCATABLE :: ZREG( :, : )
* ..
* .. External Functions ..
INTEGER ILAENV
Expand Down Expand Up @@ -1931,11 +1936,63 @@ SUBROUTINE SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
300 CONTINUE
310 CONTINUE
*
*
* SSTEDC must report a NaN in the matrix through INFO. A NaN
* off-diagonal entry used to split the matrix and was dropped, and
* a NaN diagonal entry was isolated as a 1 by 1 block, both with
* INFO = 0, whenever the block left over was large enough for the
* divide and conquer recursion, which the sizes in the input file
* do not reach. Use local arrays because LDU is a row stride,
* not the capacity of the caller's eigenvalue arrays or Z columns.
*
NSMLSZ = ILAENV( 9, 'SSTEDC', ' ', 0, 0, 0, 0 )
IF( ILAENV( 10, 'SSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 .AND.
$ ILAENV( 11, 'SSTEDC', 'V', 1, 0, 0, 0 ).EQ.1 ) THEN
N = MAX( 2, NSMLSZ+1 )
ALLOCATE( DREG( N ), EREG( N ), ZREG( N, N ) )
RONE = ONE
RNAN = SQRT( -RONE )
CALL SSTEDC( 'V', N, DREG, EREG, ZREG, N, WORK, -1, IWORK, -1,
$ IINFO )
IF( IINFO.EQ.0 .AND. INT( WORK( 1 ) ).LE.LWORK .AND.
$ IWORK( 1 ).LE.LIWORK ) THEN
DO 340 JNAN = 1, 2
DO 330 JCOMPZ = 1, 2
IF( JCOMPZ.EQ.1 ) THEN
COMPZ = 'I'
ELSE
COMPZ = 'V'
END IF
DO 320 J = 1, N
DREG( J ) = REAL( J )
EREG( J ) = ONE / REAL( J+1 )
320 CONTINUE
IF( JNAN.EQ.1 ) THEN
EREG( N / 2 ) = RNAN
ELSE
DREG( N ) = RNAN
END IF
CALL SLASET( 'Full', N, N, ZERO, ONE, ZREG, N )
CALL SSTEDC( COMPZ, N, DREG, EREG, ZREG, N,
$ WORK, LWORK, IWORK, LIWORK, IINFO )
IF( IINFO.EQ.0 ) THEN
WRITE( NOUNIT, FMT = 9985 )COMPZ, JNAN
NERRS = NERRS + 1
END IF
NTESTT = NTESTT + 1
330 CONTINUE
340 CONTINUE
END IF
DEALLOCATE( DREG, EREG, ZREG )
END IF
*
* Summary
*
CALL SLASUM( 'SST', NOUNIT, NERRS, NTESTT )
RETURN
*
9985 FORMAT( ' SCHKST: SSTEDC( ', A1, ' ) returned INFO=0 for a',
$ ' matrix with a NaN, case ', I1 )
9999 FORMAT( ' SCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
$ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
*
Expand Down
Loading
Loading