Fix Sparspak symbolic reuse invalidated by per-solve dropzeros#1036
Merged
ChrisRackauckas merged 2 commits intoJun 13, 2026
Merged
Conversation
When the DefaultLinearSolver's nonstructural_zeros feature operates in per-solve dropzeros mode (cache_union=false), successive calls to reduce_operand! can return sparse matrices with different nnz counts. Sparspak's sparspaklu! reuses the symbolic factorization from the previous call without checking whether the sparsity structure matches. This produces an invalid factorization when nnz changes, causing wrong linear-solve results and stalled nonlinear (BVP) iterations. Fix: track `reduced_nnz` and `structure_changed` on SparseReduction. The Sparspak solve! reads `structure_changed` and skips symbolic reuse (falls back to a full sparspaklu) when the structure is unstable. Fixes regression introduced by the nonstructural_zeros feature (v3.85) that caused BigFloat BVP solves to stall (SciML/BoundaryValueDiffEq.jl#509). Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
PureUMFPACKFactorization (added in SciML#1004) is an AbstractSparseFactorization but was missing from resolve.jl's sparse-conversion lists, so the Re-solve testset exercised it with a dense matrix and errored in do_factorization. Add it alongside the other sparse LU factorizations. Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
5 tasks
Contributor
Author
|
Root-cause companion: #1037 makes PureKLU the default sparse LU for generic (non-BLAS) element types, so the default path for BigFloat no longer selects Sparspak at all. This PR remains valuable as defense-in-depth for users who explicitly request |
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
2 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.
Summary
The
nonstructural_zerosfeature introduced in v3.85 (per-solvedropzerosmode,cache_union=false) can return sparse matrices with differentnnzcounts across successive Newton iterations. Sparspak'ssparspaklu!reuses the symbolic factorization from the previous call without checking whether the sparsity structure has changed — unlike KLU/UMFPACK, which both guard against this. Whennnzchanges,sparspaklu!produces a corrupt factorization, leading to wrong linear-solve results and stalled nonlinear (e.g. BVP) iterations.ext/LinearSolveSparseArraysExt.jl: Addreduced_nnz::Intandstructure_changed::BooltoSparseReduction. Updatereduce_operand!to track whether the per-solve dropzeros changed the matrix structure.ext/LinearSolveSparspakExt.jl: Before callingsparspaklu!, checkstructure_changedvia theDefaultLinearSolverInitcache. If the structure is unstable, skip symbolic reuse and fall back to a fullsparspaklu.Fixes the BigFloat BVP stalling regression reported in SciML/BoundaryValueDiffEq.jl#509 (comment from @ChrisRackauckas identifying FIRK test failures).
Test plan
test/core/nonstructural_zeros.jl— 122/122 passBigFloat+SparspakFactorizationdirect solve — correct resultsBigFloat+DefaultLinearSolver(Sparspak path) — correct results