Propagate a NaN through xLAE2 and xLAEV2 instead of returning finite eigenvalues - #1405
Open
rmlarsen wants to merge 2 commits into
Open
Propagate a NaN through xLAE2 and xLAEV2 instead of returning finite eigenvalues#1405rmlarsen wants to merge 2 commits into
rmlarsen wants to merge 2 commits into
Conversation
…eigenvalues
xLAE2 and xLAEV2 compute the eigenvalues of the symmetric 2 by 2
matrix [a b; b c] from sm = a + c, adf = |a - c| and ab = |2b|:
IF( ADF.GT.AB ) THEN
RT = ADF*SQRT( ONE+( AB / ADF )**2 )
ELSE IF( ADF.LT.AB ) THEN
RT = AB*SQRT( ONE+( ADF / AB )**2 )
ELSE
RT = AB*SQRT( TWO )
END IF
The last branch is meant for adf = ab, but a NaN diagonal entry makes
adf a NaN, both comparisons false, and rt = |2b| sqrt(2), finite. The
sign test on sm that follows is false for a NaN as well, so the
routines returned rt1 = -rt2 = |b| sqrt(2) with no trace of the NaN.
Every 2 by 2 tridiagonal eigenproblem goes through this code: xSTEQR,
xSTERF, xSTEV, xSTEVD, xSTEDC and xSTEMR returned INFO = 0 and finite
eigenvalues for a 2 by 2 matrix with a NaN on the diagonal, while a
NaN off the diagonal was propagated.
Take the last branch only for adf = ab and let a NaN in either operand
propagate, rt = adf + ab, so that the eigenvalues, and in xLAEV2 the
eigenvector, come out as NaN. Finite arguments are not affected: the
branch they take is unchanged and so is the value it computes.
xCHKST gets the case as a regression test: it calls xSTEQR, xSTERF
and xSTEMR on the 2 by 2 matrix [1 1; 1 2] with a NaN in place of
either diagonal entry and reports a finite eigenvalue as a failure.
On the parent commit all six calls return INFO = 0 and finite
eigenvalues in every precision; here every eigenvalue is a NaN. Over
the 1 by 1 and 2 by 2 NaN cases of a sweep of every tridiagonal solver
the parent returned an all-finite result in 62 of 111 cases and this
branch in 16, all of them calls that return no eigenvalue at all
(RANGE = 'V' or 'I' with M = 0) or bisection paths, which do the same
for larger n; every finite case is bit-identical.
The full LAPACK test suite passes: 0 numerical errors, 0 other errors,
120 tests more than the parent from the new calls.
The regression test fails on the parent and passes with the fix, and the
reproducer prints the same before and after output, with gfortran 13
(x86-64 Release and Debug with -fcheck=all, and under QEMU on aarch64,
ppc64le, s390x and riscv64), flang-19 and Intel ifx 2025.3.
Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## master #1405 +/- ##
========================================
Coverage 69.36% 69.37%
========================================
Files 6122 6122
Lines 486337 486485 +148
Branches 23268 23268
========================================
+ Hits 337330 337476 +146
- Misses 148569 148571 +2
Partials 438 438
Continue to review full report in Codecov by Harness.
|
Add STEQR(I) and STEMR(V,A) cases for both diagonal NaN positions in all four precisions. Include the job option in diagnostics and check the real workspace capacity for the complex drivers. Validation: 8 sep/se2 driver runs and 8 AddressSanitizer cases passed. All four precisions detect reverting only SLAEV2/DLAEV2, and all four also retain the expected failures against the parent kernels.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
1 participant
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Disclaimer: The initial changes were prepared using Claude Code; review and subsequent test fixes used Codex.
Summary
xLAE2andxLAEV2compute the eigenvalues of the symmetric 2 by 2 matrix[a b; b c]through a three-way comparison of|a-c|with|2b|, whose last branch is meant for the equal case. For a NaN diagonal entry both comparisons are false, the last branch computesrt = |2b| sqrt(2), finite, and the sign test ona+cthat follows is false too, so the routines returnedrt1 = -rt2 = |b| sqrt(2)with no trace of the NaN. Every 2 by 2 tridiagonal eigenproblem goes through this code:xSTEQR,xSTERF,xSTEV,xSTEVD,xSTEDCandxSTEMRreturnedINFO = 0and finite eigenvalues for a 2 by 2 matrix with a NaN on its diagonal (a NaN off the diagonal was propagated). This PR takes the equal-case branch only when the operands are equal and propagates a NaN otherwise. The numerical changes are in{s,d}lae2.fand{s,d}laev2.f, with regression coverage in{s,d,c,z}chkst.f;cLAEV2andzLAEV2call the real routine.Description
With
aa NaN,sm = a + candadf = |a - c|are NaN,ab = |2b|is finite, both branches of the first test are skipped,rt = ab sqrt(2), and thesmtest lands in its last branch as well.Fix. The last branch of the first test is taken only for
ADF.EQ.AB; otherwise one of the two is a NaN andRT = ADF + ABpropagates it, after whichsmis a NaN too (a NaN diagonal entry makes bothsmandadfNaN), soRT1andRT2come out as NaN. Finite arguments take exactly the branch they took before and get the same value.Minimal reproducer
+/- sqrt(2) = +/- |b| sqrt(2): the diagonal, one entry of which is a NaN, played no part in the result.Regression test.
xCHKSTcallsxSTEQR('N'),xSTERF,xSTEMR('N','A'),xSTEQR('I'), andxSTEMR('V','A')on the 2 by 2 matrix[1 1; 1 2]with a NaN in either diagonal entry. The eigenvalue-only calls exercisexLAE2, and the eigenvector-producing calls exercisexLAEV2. A successful call must return NaN eigenvalues; diagnostics identify the routine and job option. Reverting onlySLAEV2andDLAEV2makes the new vector cases fail in every precision, while the complete fix passes.Validation
-O2 -fcheck=all; all 8 focusedsep.in/se2.indriver runs pass.SLAEV2andDLAEV2are reverted.Found while auditing the symmetric tridiagonal eigensolvers for the NaN and overflow handling of #1377-#1391.
Update: Codex review identified that the eigenvalue-only regression calls did not reach LAEV2. STEQR(I) and STEMR(V,A) now exercise that path, and reverting only LAEV2 makes the added cases fail.