TESTING: Improve eigenvalue consistency tests (chkhs) - #1410
Merged
langou merged 2 commits intoSep 12, 2026
Merged
Conversation
Test 8 is described as | W(Z computed) - W(Z not computed) | over |W| ulp, and in the real drivers that is what it is: WR1 comes from xHSEQR( 'S', 'V' ) and WR2 from xHSEQR( 'S', 'N' ), so only COMPZ differs. The complex drivers have no W2. They compute the ( 'S', 'N' ) eigenvalues into W1, overwrite them with ( 'S', 'V' ), and compare against W3, which comes from ( 'E', 'N' ) -- so JOB differs too, and the check is a stronger one than it claims to be. That matters because JOB is what makes the two paths diverge. xLAQR5 sets JTOP to 1 when the Schur form is wanted and to KTOP otherwise, so once a deflation moves KTOP the block update covers the same rows of the active block but partitions them differently between xGEMM calls. Where the compiler's xGEMM is sensitive to a row's index within the call, the last bits of H differ, and on the graded matrix types an O(ulp) difference in H is an O(1) relative difference in W. LLVM Flang 21 and later on aarch64 fail the complex test for that reason, reaching 1/ulp. Give the complex drivers a W2 of their own so test 8 compares the two JOB = 'S' calls, as the description says and as the real drivers already do. COMPZ decides only whether Z is accumulated and never feeds back into the Hessenberg data, so the two lists agree exactly: measured over 200 seeds and every matrix type, 126000 ratios per precision, all zero. The JOB = 'E' path is left without a cross-check by this commit; the next one adds it back with a normalization that suits it. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
With test 8 restored to the COMPZ comparison it claims to be, nothing compares the eigenvalues xHSEQR returns for JOB = 'E' against the ones it returns for JOB = 'S'. The real drivers never did; the complex ones did, and it was the only cross-check on the cheap eigenvalues-only path. Add it back as test 11, with a normalization that suits it. The two lists cannot be required to agree to within ulp times the largest eigenvalue. Matrix types 9 to 18 are graded and non-normal, with eigenvalue condition numbers reaching 1e15, so an O(ulp) difference in H is an O(1) relative difference in W: measured on type 10 at n = 16, a relative O(ulp) perturbation of H moves the spectrum by 1e12 to 2e14 in those units, 200 trials out of 200 over the threshold of 20. Weight each difference by s(j), the reciprocal condition number of that eigenvalue, which the eigenvectors computed for tests 9 and 10 give for free. A backward error eps in H moves the eigenvalue by eps / s(j), so |dW(j)| s(j) is the backward-error-scale quantity to test. Eigenvalues with s(j) < sqrt(ulp) are left out: those are the only ones the two paths have been seen to return in a different order, and the comparison is index-wise, so including them would mean differencing eigenvalues that do not correspond. Over 200 seeds and four precisions every reordered eigenvalue had s(j) <= 1.7e-12 in double and <= 4.3e-5 in single, against cutoffs of 1.5e-8 and 4.9e-4. The check needs both eigenvector matrices, and the xTREVC call a few lines further on reuses EVECTR for its own output, so it goes directly after test 10 and the tests that followed are renumbered 12 to 17, comments included. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Codecov Report✅ All modified and coverable lines are covered by tests. Additional details and impacted files@@ Coverage Diff @@
## master #1410 +/- ##
=======================================
Coverage 69.36% 69.36%
=======================================
Files 6122 6122
Lines 486337 486413 +76
Branches 23268 23268
=======================================
+ Hits 337330 337406 +76
Misses 148569 148569
Partials 438 438
Continue to review full report in Codecov by Harness.
|
langou
approved these changes
Sep 12, 2026
Contributor
|
Hi @ACSimon33, Thanks for cleaning this up. This is a great fix.
Very good idea. Julien. |
langou
approved these changes
Sep 12, 2026
Contributor
|
The CI is failing on I don't know how to fix that, I'm afraid. |
Collaborator
Author
Hi, this will be fixed once #1371 is merged (soon). You can ignore that for now. |
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
3 participants
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.
The failing test
ATfL 22.1.0, which is what the
ubuntu-24.04-armjob uses, on the extended_64API shows the following error:Summary
xCHKHStest 8 is documented as| W(Z computed) - W(Z not computed) | / ( |W| ulp ), and in the real drivers that is exactly what it is. The complex drivers have noW2array and compared againstW3, which comes fromJOB = 'E', so they were varyingJOBas well asCOMPZ-- a stronger check than the description, and the one that fails on every LLVM Flang from 21 onwards. This PR restores test 8 to theCOMPZcomparison in all four drivers, then adds theJOBcomparison back as a new test 11 with a normalization that suits it.Commit 1: compare the test 8 eigenvalues at equal JOB
schkhs,dchkhsWR1from('S','V')vsWR2from('S','N')cchkhs,zchkhsW1from('S','V')vsW3from('E','N')W1vs a newW2from('S','N')COMPZdecides only whetherZis accumulated and never feeds back into the Hessenberg data, so the two lists agree exactly. Measured over 200 seeds and every matrix type, 126,000 ratios per precision: all zero, andmax = 0.00with no scaling of any kind.Commit 2: check JOB = 'E' against JOB = 'S' as test 11
With test 8 restored, nothing compared the eigenvalues-only path against the full Schur path. The real drivers never did and the complex ones no longer do, so the new test adds that check to all four -- the first time the real drivers have had it.
It cannot be normalized by
ulptimes the largest eigenvalue. On ZHS type 10 at n=16 a relative O(ulp) perturbation of the Hessenberg matrix moves the spectrum by 1e12 to 2e14 in those units, with 200 of 200 trials over the threshold of 20, and the eigenvalue that moves most is not a denormal straggler: it sits at 6 to 9% of the spectral radius with its own relative error 4.6% median, up to 57%.So each difference is weighted by
s(j) = |y(j)**H x(j)| / ( ||y(j)|| ||x(j)|| ), the reciprocal condition number of eigenvaluej, built from the eigenvectors tests 9 and 10 already compute. A backward errorepsinHmoves that eigenvalue byeps / s(j), so|dW(j)| s(j)is the backward-error-scale quantity to test.s(j) = 1for a normal matrix, so the diagonal and Jordan types get the unscaled comparison.Eigenvalues with
s(j) < sqrt(ulp)are skipped. Those are the only ones the two paths have been seen to return in a different order, and the comparison is index-wise, so including them would mean comparing eigenvalues that do not correspond.Results
200 seeds x 21 types x 5 parameter sets x 7 orders = 126,000 ratios per precision,
THRESHlowered to 0 so every value is recorded, flang 21.1.8 with default contraction:(nonzero / max / at-or-over-20)
Notes
sqrt(ulp)cutoff rests on measurement and on a mechanism, not on a bound. There is no theorem of the form "s(j) >= cimplies the two paths agree on the ordering": the ordering follows from deflation decisions taken across the whole sweep, and a reordering among admitted eigenvalues would report a spurious difference. Every reordered eigenvalue found in 200 seeds and four precisions lay below the cutoff: by a factor of 8,900 in double but only 11 in single precision.