MatrixCovers

This package computes covers of matrices. Given a matrix A, a cover (more specifically, a hard cover) is a matrix C that can be defined as C = a * b', where a and b are non-negative vectors. C must satisfy

\[C_{ij} \;\geq\; |A_{ij}| \quad \text{for all } i, j.\]

For a symmetric matrix the cover is symmetric (b = a), so a single vector suffices: a[i] * a[j] >= abs(A[i, j]).

An optimal cover is one for which C is "as tight as possible" in bounding A (equivalently, no larger than strictly necessary), by criteria that will be described below.

This package also supports soft covers: these satisfy C[i, j] ⪆ abs(A[i, j]), meaning that C matches or exceeds A at most indexes but not necessarily all; this intuitive notion will be made concrete below.

Why covers?

Covers provide a natural scale-covariant "summary" of a matrix. If you rescale rows by a positive diagonal factor D_r and columns by D_c, the optimal cover transforms as a → D_r * a, b → D_c * b, exactly mirroring how the matrix entries change. Scalar summaries like norm(A) or maximum(abs, A) do not have this property and therefore implicitly encode an arbitrary choice of units.

Most users will employ matrices that store pure numbers, and this package works well with such matrices. But to emphasize the scale-covariance, we'll start with an example of a 3×3 matrix whose rows and columns correspond to physical variables with different units — position in meters, velocity in m/s, force in Newtons. Loading Unitful lets the matrix carry those units itself:

julia> using MatrixCovers, Unitful

julia> L, V, F = u"m", u"m/s", u"N";  # length, velocity, force

julia> A = [1e6/L^2   1e3/(L*V)  1.0/(L*F)
            1e3/(L*V) 1.0/V^2    1e-3/(V*F)
            1.0/(L*F) 1e-3/(V*F) 1e-6/F^2];

julia> a = symcover(A);

julia> round.(typeof.(a), a; digits=6)
3-element Vector{Quantity{Float64}}:
 1000.0 m^-1
    1.0 s m^-1
    0.001 N^-1

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.

Normalizing by the cover cancels the units along with the magnitudes, leaving a matrix that is all-ones, dimensionless, and scale-invariant:

julia> round.(A ./ (a .* a'); digits=6)
3×3 Matrix{Float64}:
 1.0  1.0  1.0
 1.0  1.0  1.0
 1.0  1.0  1.0

It is worth noting that this yields 1 only for entries where the cover bound is tight; had A been, say, diagonal, then A ./ (a .* a') would also be diagonal.

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.

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

\[\sum_{i,j} \phi\!\left(\frac{|A_{ij}|}{a_i\, b_j}\right),\]

which the solvers minimize. Two penalty families are provided:

AbsLog{p} — ϕ(r) = |log r|^p (and ϕ(0) = 0). Convex in log 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.
AbsLinear{p} — ϕ(r) = |1 - r|^p. Non-convex, but unlike AbsLog these are 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.

cover_objective evaluates either penalty for a given cover:

julia> using MatrixCovers

julia> A = [4.0 2.0; 2.0 16.0];

julia> a = symcover(A)
2-element Vector{Float64}:
 2.0
 4.0

julia> cover_objective(AbsLog{1}(), a, A)   # sum of log-excesses (L1)
2.772588722239781

julia> cover_objective(AbsLog{2}(), a, A)   # sum of squared log-excesses (L2)
3.843624111345611

Both objectives are zero if and only if every constraint is exactly tight.

You can override the default penalty by supplying it as an argument to the solvers.

Choosing a cover algorithm

Function	Symmetric	Constraint	Default (or alternative) objective	Requires
`symcover`	yes	hard (`r ≤ 1`)	heuristic	—
`cover`	no	hard (`r ≤ 1`)	heuristic	—
`symcover_min`	yes	hard (`r ≤ 1`)	`AbsLog{2}` (or `AbsLog{1}`, `AbsLinear`)	native for `AbsLog{2}`; else JuMP
`cover_min`	no	hard (`r ≤ 1`)	`AbsLog{2}` (or `AbsLog{1}`, `AbsLinear`)	native for `AbsLog{2}`; else JuMP
`soft_symcover`	yes	soft (penalized)	`AbsLinear{2}` (or `AbsLog`, `AbsLinear{1}`)	native for `AbsLog`; else —
`soft_cover`	no	soft (penalized)	`AbsLinear{2}` (or `AbsLog`, `AbsLinear{1}`)	native for `AbsLog`; else —
`soft_symcover_min`	yes	soft (penalized)	`AbsLog{2}`, `AbsLinear`	native for `AbsLog{2}`; else JuMP
`soft_cover_min`	no	soft (penalized)	`AbsLog{2}`, `AbsLinear`	native for `AbsLog{2}`; else JuMP

The two soft tiers are separated by a different axis than the two hard ones. For hard covers, symcover trades optimality for speed while guaranteeing feasibility, and symcover_min is optimal. Both soft tiers minimize the same unconstrained objective, and differ instead in what they promise about reaching its minimum:

soft_symcover and soft_cover are always native and best-effort. They require no extension for any penalty, and they own their multistart — but what they return is a coordinate-descent fixed point, which for the non-convex and nonsmooth penalties need not be a minimizer.
soft_symcover_min and soft_cover_min return a true minimizer of the basin they start in, and may require an extension. AbsLog{2} is native; the AbsLinear penalties need JuMP and Ipopt; AbsLog{1} is not implemented.

Both reduce to the same trade: cheap and always available, against best quality and possibly an extra dependency.

Under AbsLog{2} the objective is convex with a single minimizer, so the tiers coincide — soft_symcover is soft_symcover_min there, and likewise for the asymmetric pair. That is the degenerate case of the contract rather than an exception to it: with one minimizer there is nothing for a best-effort descent and a minimizer to disagree about. Under AbsLog{1} they part company most sharply — the soft AbsLog{1} covers are coordinate descents that reach a deterministic fixed point rather than a minimizer, and soft_symcover_min/soft_cover_min do not accept AbsLog{1} at all.

symcover, cover, and any native implementation can be recommended for production use, possibly with relaxed convergence bounds. The heuristic solvers are particularly fast: they run in $O(mn)$ time for an $m\times n$ matrix and often land within a few percent of the objective-minimal cover (see the quality tests involving test/testmatrices.jl). Native solvers (both hard and soft) are intermediate, still roughly $O(mn)$ but requiring many iterations (and for non-convex cases, multiple start points by default) for convergence; still, they are much faster than their JuMP-counterparts, which are provided mainly as a reference.

Objective-minimal covers

symcover_min and cover_min return a cover that minimizes the chosen penalty subject to the hard constraint. For the default AbsLog{2} penalty they are solved natively (no external solver) by penalty-continuation with a damped semismooth Newton iteration:

julia> using MatrixCovers

julia> A = [1 2 3; 6 5 4];

julia> a, b = cover(A);          # fast heuristic

julia> aq, bq = cover_min(AbsLog{2}(), A);   # AbsLog{2}-minimal, native

julia> a * b'
2×3 Matrix{Float64}:
 2.16541  2.03444  3.0
 6.0      5.63709  8.31251

julia> aq * bq'
2×3 Matrix{Float64}:
 2.21042  2.0      3.0
 6.0      5.42884  8.14325

julia> round(cover_objective(AbsLog{2}(), a, b, A); digits=6)
1.146646

julia> round(cover_objective(AbsLog{2}(), aq, bq, A); digits=6)
1.141281

The native solver is near-exact (relative objective excess typically a few $\times 10^{-7}$, growing slowly with problem size) and orders of magnitude faster than a general-purpose convex solver. The other penalties — AbsLog{1} (a linear program) and the non-convex AbsLinear variants — are solved through JuMP and are loaded on demand as a package extension:

using JuMP, HiGHS   # HiGHS for AbsLog penalties
using MatrixCovers

a    = symcover_min(AbsLog{1}(), A)   # L1-minimal symmetric hard cover
a, b = cover_min(AbsLog{1}(), A)      # L1-minimal general hard cover

The soft *_min solvers divide along the same line, but not at the same place: soft_symcover_min and soft_cover_min solve AbsLog{2} natively and reach for JuMP with Ipopt only for the AbsLinear penalties. They do not accept AbsLog{1}; the soft AbsLog{1} covers are available through soft_symcover and soft_cover, which are native.

Uniqueness

The AbsLog{2}() penalty generally has a unique minimum, with one exception: row/column scaling a → γ*a, b → b/γ does not affect C and thus invisible to the objective function. For non-symmetric (i.e., not symcover) problems, the scaling of each is pinned separately by the balance convention ∑ n_i log a[i] = ∑ m_j log b[j], where n_i, m_j are the the nonzero counts of row i and column j, respectively. This convention is not scale-invariant but has no impact on the cover itself.

Other penalties may be more degenerate. The AbsLog{1}() penalty is identical over a whole face of the feasible polytope, and its members are genuinely different covers — the products a[i]*b[j] differ — that merely happen to score the same objective. To make the result deterministic, we select the one that additionally minimizes the AbsLog{2} objective.

AbsLinear penalties typically have isolated minima, so are not as degenerate as AbsLog{1}(), but these minima occur in separate basins. There is no guarantee of global optimality.

Starting points: initialize and refine

For objectives with multiple minima, the solver starts from a specified point and descends. At a lower level, this package's interface is organized in three layers:

Initializers name the starting points. initialize_symcover and initialize_cover take a strategy — :geomean, :leaveout, :diagfeasible, or :hardcover — and return that point. Each is a property of A alone; no objective is involved, so an initializer takes no penalty. A second keyword, feasible, says how the point is brought up to covering A: :inflate (the default) scales it bodily by one common factor, :boost raises only the rows touching a violated entry, and :none leaves it as it is. The hard-cover solvers need a cover, so they take one of the first two; the soft covers want :none, since forcing the geometric mean to cover A would spoil the very property that makes it the soft AbsLog{2} optimum.
The two feasible routes land on the boundary at different points, hence in different basins — which is why the choice is a named part of the start rather than an internal detail. The heuristic cover is itself a composition of these: the geometric mean, boosted, then tightened.
Refiners improve a starting point in place, and are the !-suffixed forms of the solvers: symcover_min!, cover_min!, soft_symcover!, soft_cover!, soft_symcover_min!, and soft_cover_min! validate the start, then optimize from it. Which basin they reach is the caller's choice, by construction, and supplying the start is the caller's job. The hard refiners require a start that covers A; the soft ones do not, since their objective constrains nothing — build theirs with feasible=:none.
Solvers bundle the two. symcover_min, cover_min, soft_symcover, soft_cover, soft_symcover_min, and soft_cover_min refine every start on a menu (the strategies keyword, or the multistart's own list) and return the best cover by cover_objective, so their result depends on A and not on an initialization the caller never chose.

That is the rule for the whole grid: the plain form owns the menu, so its result is a property of A; the ! form refines the one start you give it, so its result is a property of A and that start. symcover! and cover! are the exception that proves it — they are initializers, not refiners, and construct their cover from scratch rather than reading the vector passed in.

For finer control, you can run these manually:

using JuMP, Ipopt   # Ipopt for the AbsLinear penalties
using MatrixCovers

a = symcover_min(AbsLinear{2}(), A)                      # multistart over the whole menu
a = symcover_min(AbsLinear{2}(), A; strategies=(:geomean,))   # or commit to one start

a0 = initialize_symcover(A; strategy=:geomean)           # or drive it yourself
symcover_min!(AbsLinear{2}(), a0, A)

The same menu supplies the starting points of the soft_symcover and soft_cover multistarts, adding (by default) a few randomized perturbations of a base point up to a user-controllable number of starts.

For the convex AbsLog penalties the start cannot change the result, and the refiners accept one only so that the two families share an interface.

Worked example: roundoff in `A \ b`

Because a cover names each variable's natural scale, it also says how to measure a solution in units that do not depend on how the problem was parameterized.

Solving x = A \ b is contravariant: rescaling A → D*A*D and b → D*b sends x → x ./ d, while the cover is covariant, a → d .* a. The products x .* a are therefore unchanged, and ∑ᵢ |xᵢ * aᵢ| is a measure of the solution's size that is the same in every frame.

That quantity can be estimated from the magnitudes of A and b alone, without forming x at all:

julia> using MatrixCovers, LinearAlgebra

julia> A = [1e6 1e3; 1e3 4.0];

julia> b = [1.5e3, 6.0];

julia> a = symcover(A);

julia> round.(a; digits=6)
2-element Vector{Float64}:
 1000.0
    2.0

julia> mag = sum(abs(bi / ai) for (bi, ai) in zip(b, a))
4.5

The cover reports natural scales of 1000 and 2, and mag estimates the size of the solution measured against them — here within a factor of 1.5 of the truth:

julia> x = A \ b;

julia> sum(abs.(x .* a))
3.0

Both numbers are scale-invariant, so the estimate is unchanged by any diagonal rescaling of the problem:

julia> d = [0.05, 3.0];

julia> Ad, bd = d .* A .* d', d .* b;

julia> ad = symcover(Ad);

julia> sum(abs(bi / ai) for (bi, ai) in zip(bd, ad))
4.5

This makes eps(mag) a scale-invariant estimate of the roundoff floor of the sum. For a well-conditioned A, the error of the Float64 solve meets that floor:

julia> xbig = big.(A) \ big.(b);

julia> abs(sum(abs.(x .* a)) - sum(abs.(Float64.(xbig) .* a))) <= 2 * eps(mag)
true

The estimate is built from magnitudes only, so it knows nothing about the conditioning of A or about cancellation during the solve. When A is ill-conditioned the true error sits far above the floor:

julia> Aill = [1.0 -0.9999; -0.9999 1.0];

julia> bill = [0.75, 7.0];

julia> aill = symcover(Aill);

julia> magill = sum(abs(bi / ai) for (bi, ai) in zip(bill, aill));

julia> xill = Aill \ bill;

julia> xbigill = big.(Aill) \ big.(bill);

julia> err = abs(sum(abs.(xill .* aill)) - sum(abs.(Float64.(xbigill) .* aill)));

julia> err > 1e6 * eps(magill)
true

Folding in the condition number of the normalized matrix A ./ (a .* a') — itself scale-invariant, since normalizing cancels the frame — restores a usable bound:

julia> κ = cond(Aill ./ (aill .* aill'));

julia> err <= 1e3 * eps(κ * magill)
true

Index of available tools

MatrixCovers.AbsLinear
MatrixCovers.AbsLog
MatrixCovers.AbstractCoverPenalty
MatrixCovers.cover
MatrixCovers.cover!
MatrixCovers.cover_min
MatrixCovers.cover_min!
MatrixCovers.cover_objective
MatrixCovers.foreach_support
MatrixCovers.foreach_support_sym
MatrixCovers.initialize_cover
MatrixCovers.initialize_cover!
MatrixCovers.initialize_symcover
MatrixCovers.initialize_symcover!
MatrixCovers.iscover
MatrixCovers.soft_cover
MatrixCovers.soft_cover!
MatrixCovers.soft_cover_min
MatrixCovers.soft_cover_min!
MatrixCovers.soft_symcover
MatrixCovers.soft_symcover!
MatrixCovers.soft_symcover_min
MatrixCovers.soft_symcover_min!
MatrixCovers.symcover
MatrixCovers.symcover!
MatrixCovers.symcover_min
MatrixCovers.symcover_min!

Reference documentation

MatrixCovers.AbsLinear — Type

AbsLinear{p}

Penalty type for φ(r) = |1 - r|^p. Unlike AbsLog, this penalty is continuous at r = 0 (φ(0) = 1), so zero entries in A naturally contribute a constant penalty.

The resulting optimization problems are non-convex and may have multiple local minima.

source

MatrixCovers.AbsLog — Type

AbsLog{p}

Penalty type for

φ(r) = |log(r)|^p  if r > 0
       0           if r = 0

The discontinuity at r=0 prevents zero entries in A from sending the objective value to infinity.

This leads to convex optimization problems in log space. AbsLog{1} typically has a flat minimum-basin in which members of an entire family of solutions are equally good. AbsLog{2}, except in degenerate cases like [0 1; 1 0], has a unique minimum.

See also: initialize_symcover.

source

MatrixCovers.initialize_symcover — Method

a = initialize_symcover(A; strategy=:hardcover, feasible=:inflate, kwargs...)

Build a starting point for the symmetric cover of A, as consumed by symcover_min and by the soft_symcover multistart.

No penalty is taken: every strategy below is a property of A alone, so the starting point does not depend on the objective it will be refined against.

strategy names the point:

:geomean — the geometric mean of the nonzero entries of each row, and not a cover. It minimizes the soft AbsLog{2} objective exactly when every entry of A is nonzero; on a sparse support it approximates that minimum, which soft_symcover_min(AbsLog{2}(), A) returns exactly.
:leaveout — the geometric mean recomputed with the most-underweighted support entry dropped, which lands in the basin that treats that entry as effectively zero. Raises an ArgumentError when no entry can be dropped (empty support, or dropping it would empty a row). Not a cover.
:diagfeasible — a cover grown from the diagonal by nearest-neighbor propagation. Feasible by construction.
:hardcover — the tightened hard cover of symcover, which is :geomean boosted to feasibility and then tightened. Forwards maxiter to the tightening pass. Feasible by construction, so feasible has no effect on it.

feasible names how the point is brought up to covering A — that is, to a[i]*a[j] >= abs(A[i,j]), up to the roundoff of the log-domain arithmetic:

:inflate (the default) multiplies every scale by the smallest common factor that achieves coverage, moving the point bodily and leaving its shape intact.
:boost raises only the rows that touch a violated entry, so it changes the shape of the point. This is the route symcover itself takes.
:none returns the strategy's own point, with no coverage guarantee. This is what the soft covers want: forcing the geometric mean to cover A would destroy the very property that makes it the soft AbsLog{2} optimum.

The two feasible routes land on the boundary at different points, and so in different basins of the non-convex AbsLinear objectives — which is exactly why a menu of starts is worth having, and why the choice is exposed rather than fixed.

Under every setting the result is strictly positive on every row that carries support and exactly zero on every row that carries none.

An unrecognized strategy or feasible raises an ArgumentError.

source

MatrixCovers.iscover — Method

iscover(a, b, A; rtol=0, atol=0)
iscover(a, A; rtol=0, atol=0)

Test whether a and b cover A, that is, whether a[i]*b[j] >= abs(A[i,j]) for every entry. The two-argument form tests the symmetric cover a*a', and requires A to be square.

This is the inequality the whole package is organized around. cover, symcover, and the *_min solvers all satisfy it by construction — up to the roundoff noted below — so the predicate earns its keep mainly on covers that carry no such guarantee: those from soft_cover and soft_symcover, which penalize under-coverage rather than forbidding it, and covers a caller has adjusted by hand.

rtol and atol supply the slack the producing algorithm warrants, testing a[i]*b[j] >= abs(A[i,j])*(1 - rtol) - atol. Neither defaults to any slack at all, so the bare call is the exact inequality. Use rtol for roundoff that scales with entry magnitude (the log-domain arithmetic behind the heuristics warrants a few multiples of eps), and atol for the convergence tolerance of an iterative solver. atol is subtracted only when it is nonzero: abs(A[i,j]) - atol is undefined when the two carry different units, so an rtol-only check never forms it. That makes rtol the only one of the two meaningful for a dimensional A, whose entries need not share units — no single scalar atol is commensurate with all of them.

a and b must be nonnegative; a negative scale raises an ArgumentError. Zero is allowed, and is what every solver here returns for a row or column with no support. The requirement is not cosmetic: only nonnegativity makes it sound to skip the zero entries of A, which is what lets this run in time proportional to the support rather than to length(A).

eachindex(a) must match axes(A, 1) and eachindex(b) must match axes(A, 2).

See also: cover_objective, cover, symcover.

Examples

julia> A = [1.0 2.0; 3.0 4.0];

julia> a, b = cover(A);

julia> iscover(a, b, A; rtol=8eps())
true

julia> iscover([1.0, 1.0], [1.0, 1.0], A)   # a*b' = ones, which does not reach A[2,2]
false

source

MatrixCovers.soft_cover! — Function

a, b = soft_cover!(ϕ, a, b, A; maxiter=...)
a, b = soft_cover!(a, b, A; maxiter=...)

Refine the starting point (a, b) into a soft cover of A, in place, and return it. This is the asymmetric counterpart of soft_symcover! and the refiner half of soft_cover, carrying the same contract on the start: finite and strictly positive on every supported row and column, inert (and zero on output) elsewhere, and under no obligation to cover A. Build one with initialize_cover and feasible=:none. The no-ϕ form defaults to AbsLinear{2}(), matching soft_cover.

The product a[i]*b[j] is unchanged by a -> c*a, b -> b/c, so the start is read only up to that gauge: (a, b) and (2a, b/2) give the same result. The result itself is pinned to the balance convention of cover_min, as every asymmetric cover in this package is.

source

MatrixCovers.soft_cover — Method

a, b = soft_cover(ϕ, A; maxiter=200, starts=4, σ=2.0, rng=MersenneTwister(0))
a, b = soft_cover(A; maxiter=200, starts=4, σ=2.0, rng=MersenneTwister(0))

Given a matrix A, return vectors a and b approximately minimizing the soft-cover objective ∑_{i,j} ϕ(|A[i,j]| / (a[i]*b[j])). This is the asymmetric analog of soft_symcover.

Unlike cover, there is no hard coverage constraint: a[i]*b[j] may be less than |A[i,j]|, with violations penalized by ϕ.

Supported penalty functions:

AbsLog{2}(): convex, and returns its exact unconstrained minimum from a single linear solve. Identical to soft_cover_min(AbsLog{2}(), A) — with one minimizer there is nothing for a heuristic and a minimizer to disagree about.
AbsLog{1}(): initializes from the AbsLog{2}() minimum, then refines by alternating weighted-median row and column updates, reaching a deterministic and scale-covariant fixed point. As in soft_symcover, that point is not in general a minimizer: each half-sweep minimizes exactly, but the objective's nonsmoothness couples a[i] with b[j], so the descent can settle where no such sweep improves and the objective still sits materially above its minimum. soft_cover_min does not yet offer an exact AbsLog{1} alternative.
AbsLinear{2}() (default): in the inverse-scale variables u = 1 ./ a, v = 1 ./ b, the objective ∑_{i,j∈S} (1 - |A[i,j]| u[i] v[j])² (sum over the nonzero support S) is biconvex, so alternating least squares with the closed-form half-sweeps
```
u[i] = ∑_j |A[i,j]| v[j] / ∑_j (|A[i,j]| v[j])²   (dually for v[j])
```
is monotone, stopping when the relative objective decrease drops below 1e-14 or after maxiter sweeps.
AbsLinear{1}(): initializes from the AbsLinear{2}() result, then refines by alternating weighted-median updates — each row/column block is minimized exactly, so the descent is monotone. Its flat basins are broken by a deterministic lower-median tie-break, giving a scale-covariant representative.

Rows or columns of A that are entirely zero receive scale 0. As with cover, only the products a[i] * b[j] are determined by the problem; the split is fixed by the balance convention ∑ nzaᵢ log a[i] = ∑ nzbⱼ log b[j].

The objective is non-convex, so starts starting points are tried and the lowest-objective result kept: the geometric mean boosted until it covers A, the tightened hard cover cover, and — for the remaining starts — multiplicative log-normal perturbations a .* exp.(σ .* ξ), with spread σ, of that boosted point, ξ drawn from rng. Every start co-varies with an independent row/column rescaling of A and the objective is scale-invariant, so the selection is scale-covariant. The default rng is a fresh MersenneTwister(0) per call, making repeated calls (and the two frames of a covariance check) agree; pass your own rng for reproducibility you control, since default RNG streams are not stable across Julia versions. sigma is accepted as an ASCII alias for σ.

Examples

julia> A = [1 2 3; 6 5 4];

julia> a, b = soft_cover(A);

julia> a * b'
2×3 Matrix{Float64}:
 1.93288  1.97239  2.50673
 4.97144  5.07307  6.44741

source

MatrixCovers.soft_cover_min — Function

a, b = soft_cover_min(ϕ, A)
a, b = soft_cover_min(A)

Return the ϕ-minimal asymmetric soft cover of A: minimizes ∑_{i,j} ϕ(|A[i,j]|/(a[i]*b[j])) with no coverage constraints. This is the asymmetric analog of soft_symcover_min. The no-ϕ form defaults to AbsLinear{2}(), matching soft_cover.

The row/column scales are pinned to the balance convention ∑ nzaᵢ log a[i] = ∑ nzbⱼ log b[j] (nzaᵢ, nzbⱼ = nonzero counts of row i, column j), as in cover_min: the objective depends on a and b only through the products a[i]*b[j], so without a convention the split between them would be arbitrary.

Supported ϕ values and required extensions:

AbsLog{2}(): solved natively (no external solver) — the same analytic geometric-mean minimum cover computes as its initial point. Convex, so the minimizer is unique.
AbsLinear{1}(), AbsLinear{2}(): requires JuMP and Ipopt. These objectives are non-convex, so the solver returns the minimum of the basin it starts in. Rather than commit to one start, these methods refine each of strategies — the initialize_cover menu, by default (:hardcover, :geomean), taken raw (feasible=:none, since the soft objective constrains nothing) — and return the best cover found, at a cost of one solve per start. The result is the best local minimum on that menu: the multistart is a hedge against a poor basin, not a certificate of global optimality.
AbsLog{1}(): not yet implemented. The objective is an LP in log space, but its optimum is a face, and the lexicographic AbsLog{2} selection that symcover_min uses to pin one member of the corresponding hard face does not carry over: the hard face is bounded by the coverage constraints, while this one is a level set of an unconstrained piecewise- linear objective, across which the quadratic pulls far enough to cost most of the exactly tight residuals that make AbsLog{1} worth choosing.

Every start on the menu co-varies with a rescaling of A and the objective is scale-invariant, so the selection — and hence the result — is scale-covariant.

source

MatrixCovers.soft_cover_min! — Function

a, b = soft_cover_min!(ϕ, a, b, A)
a, b = soft_cover_min!(a, b, A)

Refine the starting point (a, b) into the ϕ-minimal asymmetric soft cover of A, in place. This is the asymmetric counterpart of soft_symcover_min!. The no-ϕ form defaults to AbsLinear{2}(), matching soft_cover_min, whose supported ϕ values these methods share.

a and b must be strictly positive on every supported row and column; scales on unsupported rows and columns are inert, and are zero on output. As with soft_symcover_min! — and unlike cover_min! — the start need not cover A. Build one with feasible=:none.

The product a[i]*b[j] is unchanged by a -> c*a, b -> b/c, so the start is read only up to that gauge, and the result is pinned to the balance convention of soft_cover_min.

source

MatrixCovers.soft_symcover! — Function

a = soft_symcover!(ϕ, a, A; maxiter=...)
a = soft_symcover!(a, A; maxiter=...)

Refine the starting point a into a symmetric soft cover of A, in place, and return it. The no-ϕ form defaults to AbsLinear{2}(), matching soft_symcover, whose supported ϕ values these methods share.

This is the refiner half of soft_symcover: the non-mutating form owns a multistart menu, so its result is a property of A, while this one descends from the single start you hand it, so its result is a property of A and that start. Building the start is the caller's job — see initialize_symcover, and pass feasible=:none, since a soft cover is under no obligation to cover A.

a must be finite and strictly positive on every row of A that carries support; scales on rows carrying no support are inert, and are zero on output. Unlike symcover_min!, a need not cover A — the soft objective imposes no coverage constraint.

maxiter bounds the descent sweeps; its default matches the corresponding soft_symcover method. Under AbsLog{2} the objective is convex with a unique minimizer, so the start is honored but not visible in the result.

source

MatrixCovers.soft_symcover — Method

a = soft_symcover(ϕ, A; maxiter=32, starts=5, σ=2.0, rng=MersenneTwister(0))
a = soft_symcover(A; maxiter=32, starts=5, σ=2.0, rng=MersenneTwister(0))

Given a square matrix A assumed to be symmetric, return a vector a approximately minimizing the soft-cover objective ∑_{i,j} ϕ(|A[i,j]| / (a[i]*a[j])).

Unlike symcover, there is no hard coverage constraint: a[i]*a[j] may be less than |A[i,j]|, with violations penalized by ϕ.

Supported penalty functions:

AbsLog{2}(): convex, and returns its exact unconstrained minimum from a single linear solve. Identical to soft_symcover_min(AbsLog{2}(), A) — with one minimizer there is nothing for a heuristic and a minimizer to disagree about.
AbsLog{1}(): initializes from the AbsLog{2} minimum, then refines by coordinate descent with a log-space weighted-median step, reaching a deterministic and scale-covariant fixed point. That point is not in general a minimizer: each step minimizes exactly over one coordinate, but the objective's nonsmoothness couples a[i] with a[j], so the descent can settle where no single-coordinate move improves and the objective still sits materially above its minimum. soft_symcover_min does not yet offer an exact AbsLog{1} alternative.
AbsLinear{2}() (default): non-convex; refined by coordinate descent from starts scale-covariant starting points, keeping the lowest-objective result (see below).
AbsLinear{1}(): initializes from the AbsLinear{2}() result, coordinate descent uses a weighted-median step.

For the AbsLinear penalties the objective is non-convex, so starts starting points are tried and the best kept, taken in this order: the geometric-mean minimum, the tightened hard cover, the geometric-mean minimum inflated uniformly until it covers A, a leave-one-out geometric mean that drops the support entry with the most negative log-residual (this start keeps the result continuous as an entry |A[i,j]| approaches zero), and — only when A has a zero entry — a greedy feasible cover. Any slots left over are multiplicative log-normal perturbations a .* exp.(σ .* ξ) of the geometric-mean point with spread σ, ξ drawn from rng; at the default starts=5 there is at most one such perturbation. Every start co-varies with a diagonal rescaling of A and the objective is scale-invariant, so the selection is scale-covariant. The default rng is a fresh MersenneTwister(0) per call, making repeated calls (and the two frames of a covariance check) agree; pass your own rng for reproducibility you control, since default RNG streams are not stable across Julia versions. sigma is accepted as an ASCII alias for σ.

Examples

The multistart converges to the covariant minimizer to within its objective tolerance; round to compare against exact values.

julia> A = [4 -1; -1 0];

julia> round.(soft_symcover(A); digits=4)
2-element Vector{Float64}:
 2.0
 0.5

julia> round.(soft_symcover([0 1; 1 0]); digits=4)
2-element Vector{Float64}:
 1.0
 1.0

source

MatrixCovers.soft_symcover_min — Function

a = soft_symcover_min(ϕ, A)
a = soft_symcover_min(A)

Return the ϕ-minimal symmetric soft cover of A: minimizes ∑_{i,j} ϕ(|A[i,j]|/(a[i]*a[j])) with no coverage constraints. The no-ϕ form defaults to AbsLinear{2}(), matching soft_symcover.

Supported ϕ values and required extensions:

AbsLog{2}(): solved natively (no external solver). In log space the objective is a linear least-squares, so one solve settles it, and being convex it has a unique minimizer that no start can influence. linsolve selects the inner solve, exactly as in symcover_min.
AbsLinear{1}(), AbsLinear{2}(): requires JuMP and Ipopt. These objectives are non-convex, so the solver returns the minimum of the basin it starts in. Rather than commit to one start, these methods refine each of strategies — the initialize_symcover menu, by default (:hardcover, :geomean, :leaveout), without forcing feasibility — and return the best cover found, at a cost of one solve per start.
AbsLog{1}(): not yet implemented. The objective is an LP in log space, but its optimum is a face, and the lexicographic AbsLog{2} selection that symcover_min uses to pin one member of the corresponding hard face does not carry over: the hard face is bounded by the coverage constraints, while this one is a level set of an unconstrained piecewise- linear objective, across which the quadratic pulls far enough to cost most of the exactly tight residuals that make AbsLog{1} worth choosing.

source

MatrixCovers.soft_symcover_min! — Function

a = soft_symcover_min!(ϕ, a, A)
a = soft_symcover_min!(a, A)

Refine the starting point a into the ϕ-minimal symmetric soft cover of A, in place. This is the soft counterpart of symcover_min!, and the second half of the initialize/refine pair whose first half is initialize_symcover. The no-ϕ form defaults to AbsLinear{2}(), matching soft_symcover_min, whose supported ϕ values these methods share.

a must be strictly positive on every row of A that carries support; scales on rows carrying no support are inert, and are zero on output. Unlike symcover_min!, a need not cover A — the soft objective imposes no coverage constraint, and the natural starts do not satisfy one. Pass feasible=:none when building a start with initialize_symcover.

The AbsLinear penalties are non-convex, so the start selects the local minimum the solver descends into; that is why soft_symcover_min tries several rather than committing to one. Under AbsLog{2} the objective is convex with a unique minimizer, so the start is honored but not visible in the result.

source

MatrixCovers.symcover! — Method

a = symcover!(ϕ, a, A; maxiter=3)
a = symcover!(a, A; maxiter=3)

Mutating counterpart of symcover: writes the symmetric hard cover into a and returns it, rather than allocating a new vector. eachindex(a) must match axes(A, 1) (and A must be square). ϕ has the same meaning as in symcover, and is likewise ignored by the current heuristic covers.