Add PureUMFPACKFactorization extension (PureUMFPACK.jl)#1004
Merged
ChrisRackauckas merged 3 commits intoJun 11, 2026
Conversation
Member
|
This should also get an opt-out #1023 |
Member
|
No opt out since it uses blas |
Adds an opt-in `PureUMFPACKFactorization <: AbstractSparseFactorization` backed by the pure-Julia UMFPACK port PureUMFPACK.jl, mirroring the existing SuiteSparse `UMFPACKFactorization` and `PureKLUFactorization`. PureUMFPACK is a weak dependency loaded via the new `LinearSolvePureUMFPACKExt` extension; the SuiteSparse `UMFPACKFactorization` and the default polyalgorithm are unchanged. PureUMFPACK has no in-place numeric-refactorization (`lu!`-style) API, so each fresh solve refactorizes via `splu`; `reuse_symbolic`/`check_pattern` are accepted for API parity but cannot skip the symbolic step. Singular factors (zero on U's diagonal under `check=false`) map to `ReturnCode.Infeasible`. PureUMFPACK.jl is not yet registered in General, so the test adds it by URL. Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
ChrisRackauckas-Claude
force-pushed
the
add-pureumfpack-ext
branch
from
c904868 to
9a9db91
Compare
ChrisRackauckas
approved these changes
Jun 11, 2026
ChrisRackauckas-Claude
pushed a commit
to ChrisRackauckas-Claude/LinearSolve.jl
that referenced
this pull request
Jun 13, 2026
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>
ChrisRackauckas
added a commit
that referenced
this pull request
Jun 13, 2026
* Test PureUMFPACKFactorization with a sparse matrix in resolve.jl PureUMFPACKFactorization (added in #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> * Default to PureKLU for generic-eltype sparse LU instead of Sparspak For non-BLAS sparse element types (e.g. BigFloat), the default sparse LU was SparspakFactorization, which requires `using Sparspak` and whose symbolic-factorization reuse produces invalid factorizations when the sparsity pattern changes across solves (the per-solve dropzeros of the nonstructural_zeros path). That manifested as stalled BVP Newton iterations for BigFloat problems (SciML/BoundaryValueDiffEq.jl#509). PureKLU is a pure-Julia, hard dependency that factors any `Number` element type and correctly re-analyzes on pattern changes, so make it the default sparse LU for generic eltypes. Benchmarks (BigFloat, n=10..200) show PureKLU is competitive-to-faster than Sparspak (e.g. n=200: 0.219s vs 0.307s single solve) with comparable accuracy, and it needs no extra `using`. - Generic-eltype `init_cacheval` for PureKLU and SparseColumnPivotedQR so the default polyalgorithm's KLU and SPQR slots are concretely typed for non-BLAS eltypes (the latter is required for the singular-LU QR fallback to `setfield!` a real factorization). - `defaultalg(::AbstractSparseMatrixCSC)` for square systems now selects the KLU (PureKLU) slot. Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com> --------- Co-authored-by: Chris Rackauckas (Claude) <accounts@chrisrackauckas.com>
ChrisRackauckas
added a commit
that referenced
this pull request
Jun 13, 2026
* Fix Sparspak symbolic reuse with per-solve dropzeros 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> * Test PureUMFPACKFactorization with a sparse matrix in resolve.jl PureUMFPACKFactorization (added in #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> --------- Co-authored-by: Chris Rackauckas (Claude) <accounts@chrisrackauckas.com>
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
Adds an opt-in
PureUMFPACKFactorization <: AbstractSparseFactorizationthat wraps the pure-Julia UMFPACK port PureUMFPACK.jl, mirroring the existing SuiteSparseUMFPACKFactorizationand the pure-JuliaPureKLUFactorization.PureUMFPACKis a weak dependency; the algorithm is provided by a newLinearSolvePureUMFPACKExtextension and is only available afterusing PureUMFPACK.UMFPACKFactorizationand the default polyalgorithm are unchanged — this is purely additive/opt-in.Details
src/factorization.jl:PureUMFPACKFactorization(; reuse_symbolic = true, check_pattern = true)struct + docstring + genericinit_cacheval(...) = nothingstub, placed next toKLUFactorization/PureKLUFactorization.src/LinearSolve.jl: exportPureUMFPACKFactorization.Project.toml:PureUMFPACKadded to[weakdeps],[compat],[extras], thetesttarget, and the new extensionLinearSolvePureUMFPACKExt = ["PureUMFPACK", "SparseArrays"].ext/LinearSolvePureUMFPACKExt.jl:init_cacheval(returns an empty typedPureLUfor sparse inputs so the cache field type can hold the real factor) andSciMLBase.solve!modeled on theUMFPACKFactorization/PureKLUFactorizationmethods.Caveat on caching semantics
PureUMFPACK exposes
splu(full factorization) andsolve/\; it has no in-place numeric-refactorization (lu!-style) API and no symbolic-reuse API. So each fresh solve (cache.isfresh) refactorizes viasplu.reuse_symbolic/check_patternare accepted for API parity withUMFPACKFactorizationbut cannot skip the symbolic step here; this is documented in the docstring and code comments. Cached factorizations are still reused across solves that only change the right-hand side. Singular factors (a zero onU's diagonal undercheck = false) are mapped toReturnCode.Infeasible.Registration
PureUMFPACK.jl is not yet registered in the General registry, so
test/pureumfpack.jladds it viaPkg.add(url = "https://github.com/SciML/PureUMFPACK.jl.git")before running.Tests
test/pureumfpack.jl(wired into theCore/Allgroup ofruntests.jl) verifies the extension loads (Base.get_extension(...) !== nothing), a square sparse solve, caching + value-only refactor,reuse_symbolic = false, and singular ->Infeasible. Run locally on Julia 1.10, all passing.Please ignore until reviewed by @ChrisRackauckas.
🤖 Generated with Claude Code