From 544c67666347e3d2dc8f8f20a5208b6ac7cce6ab Mon Sep 17 00:00:00 2001 From: Tim Holy Date: Sat, 1 Aug 2026 13:55:45 -0500 Subject: [PATCH 1/2] Pin what gramcover returns when there is one component A single support component leaves the per-component construction nothing to accumulate separately, so it returns exactly the global pair bound norm(a)*b. Only the block-diagonal case was covered, where the two differ by a wide margin and an inequality suffices; with one component the bound is attained, and an inequality alone would also pass for a construction that gave up slack it did not have to. The new testset states the agreement two-sided. It also records why the bare `all(s .<= norm(a) .* b)` is not the assertion to write: each sqrt is inflated by 1 + (n+3)*eps so the cover holds despite naive summation, while the right-hand side carries no such margin, so the inequality is false by a few ulps wherever it is tight. The tolerance is that inflation plus the reference's own roundoff, and coverage is asserted separately, since coverage is what the margin exists to deliver. Assisted-by: Claude Opus 5 --- test/gram_covers.jl | 24 ++++++++++++++++++++++++ test/runtests.jl | 1 + 2 files changed, 25 insertions(+) diff --git a/test/gram_covers.jl b/test/gram_covers.jl index dd38a20..864fe01 100644 --- a/test/gram_covers.jl +++ b/test/gram_covers.jl @@ -17,6 +17,30 @@ end end + @testset "one component: the global bound is attained" begin + # With a single support component there is nothing to accumulate separately, so the + # per-component construction returns exactly the global pair bound `norm(a)*b`. The + # agreement is two-sided: an inequality alone would also pass for a construction that + # gave up slack it did not have to. + # + # It is not bitwise, and the deviation is one-sided by design. Each `sqrt` is inflated by + # `1 + (n+3)*eps` so the cover holds despite naive summation; `norm` accumulates its own + # `n` terms with no such margin. So the tolerance is `(2n+3)*eps` — the declared inflation + # plus the reference's own roundoff — and `all(s .<= norm(a) .* b)` is *false* by a few + # ulps. That inequality is exposition, not a test. + for (seed, m, k) in ((3, 5, 4), (7, 4, 4), (11, 12, 3)) + rng = StableRNG(seed) + J = randn(rng, m, k) + a, b = cover(J) + s = gramcover(a, b, J) + @test ncomponents(support_components(J)) == 1 + @test s ≈ norm(a) .* b rtol = (2 * length(a) + 3) * eps(Float64) + @test !all(s .<= norm(a) .* b) + # Coverage, which the margin exists to deliver, is unconditional. + @test all(s * s' .>= abs.(J' * J)) + end + end + @testset "block-diagonal: per-component structure" begin rng = StableRNG(1) B = randn(rng, 4, 3) diff --git a/test/runtests.jl b/test/runtests.jl index e2bbdf0..00412e0 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,5 +1,6 @@ using MatrixCovers using MatrixCovers: foreach_support, foreach_support_sym, unconstrained_min!, tighten_cover! +using MatrixCovers: ncomponents, support_components using JuMP, HiGHS, Ipopt # triggers MatrixCoversJuMPExt and MatrixCoversIpoptExt extensions using SparseArrays # triggers MatrixCoversSparseArraysExt extension using Unitful # triggers MatrixCoversUnitfulExt extension From 5f835ae5b950393021dae0d1864d9bccd9480127 Mon Sep 17 00:00:00 2001 From: Tim Holy Date: Thu, 6 Aug 2026 03:04:27 -0500 Subject: [PATCH 2/2] Wording tweaks to documentation Adds: - hints on potential applications to the README - polish on documentation phrasing --- README.md | 9 +++++++-- docs/src/index.md | 25 +++++++++++++------------ 2 files changed, 20 insertions(+), 14 deletions(-) diff --git a/README.md b/README.md index 0ab7b0c..554c233 100644 --- a/README.md +++ b/README.md @@ -6,11 +6,16 @@ [![Coverage](https://codecov.io/gh/HolyLab/MatrixCovers.jl/graph/badge.svg?token=trG4HXo9N4)](https://codecov.io/gh/HolyLab/MatrixCovers.jl) [![Aqua QA](https://juliatesting.github.io/Aqua.jl/dev/assets/badge.svg)](https://github.com/JuliaTesting/Aqua.jl) -This package computes **covers** of matrices: non-negative vectors `a` (and `b`) +This package computes **covers** of matrices: non-negative vectors `a` and `b` such that `a[i] * b[j] >= abs(A[i, j])` for all `i`, `j`. Covers are the natural scale-covariant representation of a matrix — under row/column diagonal scaling they transform exactly as the matrix entries do — making them a useful -building block for scale-invariant numerical analysis. +building block for scale-invariant numerical analysis. In particular, +``\hat A = A ./ (a * b^T)`` is scale-invariant, and because ``|\hat A[i, j]| <= 1`` +for all `i` and `j`, this simple construct finds applications that range from +[statistical normalization](https://en.wikipedia.org/wiki/Normalization_(statistics)) +of data to the design of well-behaved numerical algorithms (thanks, e.g., +to [bounds on ``\hat A``'s eigenvalues](https://en.wikipedia.org/wiki/Gershgorin_circle_theorem)). Fast O(mn) heuristics (`symcover`, `cover`) are provided for everyday use, along with *soft* covers (`soft_symcover`, `soft_cover`) that penalize under-coverage diff --git a/docs/src/index.md b/docs/src/index.md index 65440c0..0cb24cd 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -60,9 +60,10 @@ julia> round.(typeof.(a), a; digits=6) `A[i,j]` has units `1/(u[i]*u[j])` (modeling a Hessian matrix for functions of parameter vectors with units `u[i]`), so `a[i]` comes back with gradient-like -units of `1/u[i]`: the cover names each variable's natural scale outright, here -1 mm, 1 m/s, and 1 kN. Had we expressed the original matrix in those units -directly, we would have gotten the equivalent cover stated in those units. +units of `1/u[i]`: the cover provides each variable's natural scale inferred +from `A`, here 1 mm, 1 m/s, and 1 kN. Had we expressed the original matrix in +those units directly, we would have gotten the equivalent cover stated in those +units. Normalizing by the cover cancels the units along with the magnitudes, leaving a matrix that is all-ones, dimensionless, and scale-invariant: @@ -80,18 +81,18 @@ It is worth noting that this yields 1 only for entries where the cover bound is A cover exists only when the units of `A` factor as `unit(A[i,j]) == unit(a[i])*unit(b[j])`, and a matrix that fails this is rejected with a -`DimensionMismatch` naming the entries that conflict. The requirement is not one -this package adds: without it the terms in a row of `A*x` carry different units, -so `A*x` is undefined for every `x`. If a matrix can be used in matrix-vector -multiplication, it has a cover. +`DimensionMismatch`. This requirement is not exhorbitant: without it, the terms +in a row of `A*x` carry different units, so `A*x` is undefined for every `x`. +If a matrix can be used in matrix-vector multiplication, it has a cover. ## Penalty functions A cover is valid as long as every constraint is satisfied, but tighter covers better capture the scaling of `A`. Cover quality is measured through the ratios -`r_{ij} = |A[i,j]| / (a[i] * b[j])`: a hard cover has every `r ≤ 1`, and `r = 1` -means the constraint is exactly tight. A **penalty function** `ϕ` combines those -ratios into a scalar objective +`r[i, j] = |A[i,j]| / (a[i] * b[j])`: a hard cover has every `0 ≤ r[i, j] ≤ 1`, +and `r[i, j] == 1` means the constraint is exactly tight. + +A **penalty function** `ϕ` combines those ratios into a scalar objective ```math \sum_{i,j} \phi\!\left(\frac{|A_{ij}|}{a_i\, b_j}\right), @@ -103,9 +104,9 @@ which the solvers minimize. Two penalty families are provided: space, which makes them a favorable (and therefore default) penalty for *hard* covers, where `r ≤ 1` and `|log r|` is the log-excess of a constraint. `AbsLog{1}` sums the log-excesses (L1), `AbsLog{2}` sums their squares (L2). - Their principal disadvantage is the discontinuity at `r = 0`. + Their principal disadvantage is the divergence and discontinuity at `r = 0`. - [`AbsLinear`](@ref)`{p}` — `ϕ(r) = |1 - r|^p`. Non-convex, but unlike - `AbsLog` these are continuous at `r = 0` (`ϕ(0) = 1`), so zero entries of `A` + `AbsLog` these are finite and continuous at `r = 0` (`ϕ(0) = 1`), so zero entries of `A` contribute a bounded penalty. This is the penalty used by default for the *soft* covers, where `r > 1` (an uncovered entry) is allowed but penalized.