Speed up [sym]cover_min - #45
Merged
Merged
Conversation
The sparse support traversal and the sparse `linsolve` defaults live in `src/sparse_support.jl`, and the Unitful extension carries the sparse methods that resolve its ambiguity with them. Assisted-by: Claude Fable 5 <noreply@anthropic.com> Assisted-by: Claude Opus 5 <noreply@anthropic.com>
Each AbsLog{2} Newton step solves reweighted normal equations whose
matrix, on a near-complete support, is the complete-support matrix
corrected by the zero set and by the currently violated entries. That
split is `C + U*Uᵀ` with `C` sparse and positive definite and `U` of rank
one (symmetric) or two (asymmetric), so a sparse Cholesky of `C` plus a
Woodbury update replaces the dense factorization. `linsolve=:woodbury`
selects it, `:auto` takes it wherever its requirements hold, and
`stats.linsolve` reports which path ran.
Assisted-by: Claude Fable 5 <noreply@anthropic.com>
Assisted-by: Claude Opus 5 <noreply@anthropic.com>
The dense and Woodbury Newton steps minimize the quadratic model of the current stage exactly. When such a step is taken whole and the violated set at the new iterate is the one the step froze its weights on, the gradient of the convex stage objective there is the model's, which vanishes: the iterate is the stage minimizer and the stage ends without a confirmation solve. The inexact `:lsqr` steps keep the decrease test as their sole criterion. Assisted-by: Claude Fable 5 <noreply@anthropic.com> Assisted-by: Claude Opus 5 <noreply@anthropic.com>
`C + U*Uᵀ` applies at O(n + |Z| + |V|) without being formed, and Gershgorin bounds its condition number by `1 + (κ−1)*2*maxdeg(V)/n`. Where that estimate stays under 1000, Jacobi-preconditioned conjugate gradients reach the same answer to rounding in a few hundred such applications; the sparse Cholesky serves the ill-conditioned stages and any run that exhausts the iteration cap. `stats.cgiters` counts the iterations. The objective sweep now reports whether the violated set still matches the one the step froze its weights on, so the line search and the stage's stopping test share one pass over the support. Assisted-by: Claude Fable 5 <noreply@anthropic.com> Assisted-by: Claude Opus 5 <noreply@anthropic.com>
The penalty strength enters the least-squares matrix only through the
rows it reweights, so `M = diag(RᵀR) + (κ−1)*Σ_{e∈V} rₑ*rₑᵀ` carries all
of it. Every generalized eigenvalue of `(RᵀWR, M)` is a mediant of
eigenvalues of `(RᵀR, diag(RᵀR))` and so lies in their range, which
holds the iteration count fixed as the continuation raises κ. LSQR runs
on `√W*R*K⁻ᵀ`, `K` the permuted sparse Cholesky factor of `M`, applied
through CHOLMOD's factor components; the violated rows join `M` only
once diagonal scaling alone would leave the system ill conditioned.
CHOLMOD is Float64-only, so other working types keep the plain
matrix-free iteration.
Assisted-by: Claude Fable 5 <noreply@anthropic.com>
Assisted-by: Claude Opus 5 <noreply@anthropic.com>
The conjugate-gradient residual recurrence drifts from `f − B x`, so convergence is confirmed against a freshly computed residual and a disagreement restarts the iteration; only a confirmed result is reported exact, which is what the sign-stability stopping test needs. A support thin enough per row can still carry a quadratic number of zeros, which would make the sparse correction dense work, so the Woodbury path also requires the total zero count to be O(n). The LSQR preconditioner is a plain scaling while its condition-number estimate stays low, forming and factorizing nothing. `stats.cholsolves` counts the Woodbury solves that reached the factorization, the two condition-number gates are named constants, and the Woodbury workspace is allocated only on that path. Assisted-by: Claude Fable 5 <noreply@anthropic.com> Assisted-by: Claude Opus 5 <noreply@anthropic.com>
`C` is assembled sparsely on every Woodbury solve, in full rather than in one triangle, so the conjugate-gradient sub-path applies it with `mul!` and adds the low-rank term directly instead of walking the edge lists. The symmetric Sherman-Morrison and the asymmetric 2x2 Woodbury update become one `_woodbury_solve!` over a low-rank block of k columns, taking `[f U]` through the factor in a single multi-right-hand-side solve. Assisted-by: Claude Fable 5 <noreply@anthropic.com> Assisted-by: Claude Opus 5 <noreply@anthropic.com>
Codecov Report❌ Patch coverage is
Additional details and impacted files@@ Coverage Diff @@
## main #45 +/- ##
==========================================
- Coverage 98.84% 98.41% -0.44%
==========================================
Files 14 13 -1
Lines 2417 2958 +541
==========================================
+ Hits 2389 2911 +522
- Misses 28 47 +19 ☔ View full report in Codecov by Harness. 🚀 New features to boost your workflow:
|
Merged
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.
The analysis of the unconstrained problem established a fast (
O(n^3) -> O(n)) (though sometimes approximate) approach based on the Woodbury matrix identity and factorizations of the support. This brings the same architecture to the[sym]cover_minsolvers. It's more complex here because it acts over violations, which change from iteration to iteration.Nevertheless the speedups are quite dramatic. Timings (σ=1 lognormal): sym n=1000 5.2 s → 1.35 s; asym 600×500 5.6 s → 0.67 s; sparse n=10⁵ LSQR 747 s → 29 s