Keep the divide and conquer workspace of cLAED0 and zLAED0 out of the caller's Z - #1406
Open
rmlarsen wants to merge 2 commits into
Open
Keep the divide and conquer workspace of cLAED0 and zLAED0 out of the caller's Z#1406rmlarsen wants to merge 2 commits into
rmlarsen wants to merge 2 commits into
Conversation
… caller's Z cLAED0 and zLAED0 solve the leaf problems into QSTORE, the N by N workspace with leading dimension LDQS, and then merge them with xLAED7 while passing the caller's Q, which has leading dimension LDQ, as the QSIZ*N packed workspace of xLAED7: CALL ZLAED7( ..., D( SUBMAT ), QSTORE( 1, SUBMAT ), LDQS, ..., $ Q( 1, SUBMAT ), RWORK( IWREM ), ... ) xLAED7 and xLAED8 address that workspace with leading dimension QSIZ, so whenever LDQ > QSIZ the merge writes into rows QSIZ+1 to LDQ of Q, which belong to the caller. For cSTEDC and zSTEDC with COMPZ = 'V' and N > SMLSIZ that is every entry of Z(N+1:LDZ, :) for a user whose Z is a block of a larger array: with LDZ = N + 8 the rows N+1 to N+8 of the first columns come back overwritten for any matrix. The real routines keep their workspace in WORK and are not affected. Keep the eigenvectors in Q, with its own leading dimension, copying each leaf product back from QSTORE, and hand QSTORE, which is contiguous and holds LDQS*N >= QSIZ*MATSIZ entries, to xLAED7 as its workspace. The final permutation of the columns goes through QSTORE and back into Q. This costs one copy of QSIZ by MATSIZ per leaf and one of QSIZ by N at the end, O(N^2) against the O(N^3) merges, and changes no arithmetic: the eigenvalues and eigenvectors are bit-identical to the parent commit. The documentation of QSTORE and of its leading dimension, LDQS >= max(1,QSIZ), which the leaf products already required, is updated with it. cCHKST and zCHKST get the case as a regression test for the divide and conquer path, which the sizes in sep.in (N <= 20, below SMLSIZ = 25) never reach: after the size and type loops they call xSTEDC with COMPZ = 'V' on an N = LDU - 1 tridiagonal matrix with the row below N of Z set to a marker, and report an overwritten marker, a nonzero INFO, or a residual or orthogonality ratio above the threshold as a failure. On the parent commit the marker row is overwritten in every column in both precisions; here it is intact and the ratios pass. The full LAPACK test suite passes: 0 numerical errors, 0 other errors, 40 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>
❌ 4 Tests Failed:
View the top 3 failed test(s) by shortest run time
View the full list of 76 ❄️ flaky test(s)
To view more test analytics, go to the Test Analytics Dashboard |
Allocate local eigenvalue arrays and a matrix with one padding row, using a problem size above SMLSIZ. LDU is only a row stride and does not establish the caller's array capacities. Preserve padding, residual, and orthogonality checks. Validation: 4 sep/se2 driver runs and 16 AddressSanitizer capacity cases passed. Both precisions retain the expected padding 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
cLAED0andzLAED0solve the leaf problems of the divide and conquer recursion intoQSTORE, theNbyNworkspace with leading dimensionLDQS, and then merge them withxLAED7while passing the caller'sQ, which has leading dimensionLDQ, as theQSIZ*Npacked workspace ofxLAED7.xLAED7andxLAED8address that workspace with leading dimensionQSIZ, so wheneverLDQ > QSIZthe merge writes into rowsQSIZ+1toLDQofQ, which belong to the caller. ForcSTEDCandzSTEDCwithCOMPZ = 'V'andN > SMLSIZthat isZ(N+1:LDZ, :)wheneverZis stored in a larger array, for any matrix: the memory is inside the array the caller passed, so nothing crashes, but the rows belowNcome back overwritten. This PR keeps the eigenvectors inQwith its own leading dimension and givesQSTORE, which is contiguous, toxLAED7as its workspace. The numerical changes are in{c,z}laed0.f, with regression coverage in{c,z}chkst.f. The real routines keep their merge workspace inWORKand are not affected.Description
The
Q(1,SUBMAT)actual argument isxLAED7'sWORK, documented as "dimension (QSIZ*N)".xLAED7passes it toxLAED8asQ2withLDQ2 = QSIZand toxLACRMas aQSIZbyKmatrix with leading dimensionQSIZ, so the workspace is laid out as ifQhad leading dimensionQSIZ; withLDQ > QSIZcolumnjof the workspace starts at offset(j-1)*QSIZofQ, which is inside the padding of an earlier column. The comment is right thatQis free at that point, since the leaf products have moved its contents intoQSTORE; the leading dimension was not accounted for.Fix. After each leaf product
xLACRMintoQSTORE, the block is copied back intoQ( 1, SUBMAT )withLDQ, so the eigenvectors accumulate inQ;xLAED7getsQ( 1, SUBMAT )withLDQas its eigenvector matrix andQSTOREas its workspace, which holdsLDQS*N >= QSIZ*MATSIZcontiguous entries. The final permutation of the columns goes throughQSTOREand back intoQ. This costs one copy ofQSIZbyMATSIZper leaf and one ofQSIZbyNat the end,O(N^2)against theO(N^3)merges, and changes no arithmetic: the eigenvalues and eigenvectors are bit-identical to the parent commit. The documentation ofQSTOREand ofLDQS >= max(1,QSIZ), which the leaf products already required, is updated with it.Minimal reproducer
Regression test.
cCHKSTandzCHKSTcallxSTEDC('V')on a tridiagonal matrix of orderN = max(2, SMLSIZ+1)(26 by default). The case owns its diagonal and off-diagonal arrays, eigenvalues, and a(N+1)byNeigenvector matrix. The extra row contains a marker. This reaches divide and conquer independently of the input sizes and respects the caller's array capacities. An overwritten marker, nonzeroINFO, or excessive residual/orthogonality ratio fromxSTT21is a failure. The parent overwrites 25 markers in both precisions; the fix preserves all markers and passes both ratios.Validation
-O2 -fcheck=all; all 4 focusedsep.in/se2.indriver runs pass.SMLSIZ = 25and40, and verify that the dedicated regression leaves caller matrix and eigenvalue arrays untouched.Found while auditing the symmetric tridiagonal eigensolvers for the NaN and overflow handling of #1377-#1391: the sweep's canaries caught the write for a NaN input, and the reproducer shows it is independent of the input.
Update: Codex review identified the regression test's reliance on LDU as array capacity. Dedicated arrays now preserve the padding, residual, and orthogonality checks without accessing caller storage beyond its documented bounds.