Route the structured-sparse default LU to SupernodalLUFactorization (stacked on #1102)

[ChrisRackauckas-Claude]

⚠️ Draft — please ignore until reviewed by @ChrisRackauckas.

Stacked on #1102 — merge that first. Split out of #1102 so the additive algorithm and this behavior change can be judged separately. GitHub shows this against main, so the diff currently also contains #1102's commits (the vendored src/SupernodalLU solver); the change unique to this PR is the last commit only — src/default.jl, the sparse defaultalg in LinearSolveSparseArraysExt, and the docs wording (3 files, ~40 lines). Once #1102 merges, this diff reduces to exactly that.

What

The default polyalgorithm's :UMFPACKFactorization slot resolves to SupernodalLUFactorization() — the same pattern as the :KLUFactorization slot resolving to PureKLUFactorization(). Explicit UMFPACKFactorization() requests are unchanged.

Occupying the existing slot rather than adding a new one means the established sparse-LU → column-pivoted-QR singular fallback chain applies unchanged, and the DefaultAlgorithmChoice enum is untouched (so default_algs.jl passes as-is).

Side effect worth calling out: both sparse LU default slots now resolve to pure-Julia solvers, so the LU side of the sparse default no longer depends on Base.USE_GPL_LIBS. The non-GPL build previously sent structured systems to the scalar KLU path (131–195 s on poisson3d_48 vs 3.8 s here). The QR side keeps its SPQR gating.

Benchmark basis

SciMLBenchmarks SparsePDE fdmatrix sweep (1D/2D/3D to n≈164k) plus 18 SuiteSparse-collection matrices via MatrixDepot, accuracy-aware (residual ≤ 1e-6 required to count as a win — a fast wrong answer doesn't):

• UMFPACK won 0 of 45 PDE-sweep rows and 1 of 18 collection matrices. Its one win is the bayer01 class: fully unsymmetric pattern (structural symmetry ≈ 0) but numerically tame, where COLAMD-style unsymmetric ordering avoids the A+Aᵀ symmetrization penalty. SupernodalLU still won the refactorization on that same matrix.
• KLU keeps banded/1D at every size, circuits at every size tested (to n=171k), and small n (≲1k 3D-like, ≲5–10k 2D-like) — that is the use_klulike_sparse_structure side of the split and is unaffected by this PR.
• Robustness: on goodwin, SupernodalLU is the only accurate solver of the five tested (KLU res 6e2, MKL PARDISO 0.31, UMFPACK 6e-4, here 1.2e-9).

Known boundary case (not fixed here)

Large circuit matrices (scircuit, n=171k) fall on the structured side of use_klulike_sparse_structure, but KLU is ~1.9× faster there. Density alone can't separate circuits from 2D meshes — torso2 is sparser than scircuit yet SupernodalLU wins it 3×. A predicted-fill probe (KLU-like iff estimated nnz(L)/nnz(A) ≲ 4, available from the AMD analysis both libraries already run) would fix the routing; left as a follow-up so this PR stays a single decision.

Verification

GROUP=Core passes locally. Default routing verified end to end: structured sparse selects SupernodalLUFactorization, the cacheval is a SupernodalLUFactor, and the refactorization and singular-fallback paths were exercised (consistent singular system falls back to sparse QR and solves to residual 0.0).

🤖 Generated with Claude Code

https://claude.ai/code/session_017rr42T1vFRRT8D7DMAmJzf

[ChrisRackauckas-Claude]

Rebased onto merged main (#1102 is in), so this PR's diff is now exactly the defaults change: src/default.jl, the sparse defaultalg in LinearSolveSparseArraysExt, and the docs wording — 3 files, ~40 lines, one commit. Verified against the merged main: structured sparse resolves to SupernodalLUFactorization and solves at 3.5e-17 residual.