eliminate over-approximates the projection by up to ~2.5e-4 on some instances

[cl445]

Summary

eliminate can return a projected polytope whose bounding box is noticeably larger than the true projection, by up to ~2.5e-4 on some instances. Since a projection preserves the extent in the kept dimensions, bounding_box(eliminate(p, ids)) should equal bounding_box(p) in those dimensions. On degenerate random polytopes it does not, and the gap is far above the solver tolerance.

Reproduction

using PolyDAQP, Random
Random.seed!(1305)
A, b = randn(6,40), rand(40)
p  = Polyhedron(A, b)
ub, lb = bounding_box(p)          # exact bounding box of p (all 6 dims)
pn = eliminate(p, [4,5,6])        # project onto dims 1,2,3
ubn, lbn = bounding_box(pn)
ubn[3] - ub[3]                    # 0.00025336583177260064, should be ~0

The eliminated polytope is a strict over-approximation (ubn > ub, never smaller): a facet that should bound the projection is missing.

Root cause

The over-approximation comes from the combinatorial pruning in eliminate, in two parts.

1. The sum(Hcand) > level+1 Chernikov filter drops real facets. That cardinality bound (a combination's history has at most level+1 original constraints) only holds under non-degeneracy. Random polytopes are degenerate, so a genuine facet of the projection can have a larger history and is discarded. On seed 1305 (error in the dim-3 upper bound):

variant	error	#constraints
current (cardinality filter)	2.5e-4	116
subset-adjacency filter	5.9e-8	129
no combinatorial filter	6.7e-16	4658

Disabling the filter entirely makes it exact, which confirms the filter is what drops the facet, but that blows the intermediate constraint count up about 40x, so it is not a usable fix.

A cheap adjacency filter is not sufficient either. Replacing the cardinality test with the standard subset-adjacency test (skip a combination only if some other constraint's history is a subset of the combined one) is never worse than the current code over 2000 random seeds and fixes the single worst case, but it still leaves 8/2000 failures and a 5.5e-5 worst case, because the subset test also has false positives on degenerate instances.

2. A second, smaller source is minrep. Even with no combinatorial filter, a residual around 1e-7 remains, from the redundancy check running at DAQP's default primal_tol = 1e-6: a constraint protruding by about 1e-7 is classified as redundant and dropped.

Bottom line

The accuracy eliminate can currently deliver on degenerate inputs is around 1e-4 with the fast filter. The elimination testset in test/runtests.jl asserts this to δ = 1e-8, so it fails whenever CI lands on one of these instances (8 of 2000 random seeds here). The follow-up comment identifies the fix: replacing the cardinality filter with history-dedup removes the macroscopic errors (machine epsilon on most seeds). A residual around 1e-7 remains from LP-based redundancy removal on the degenerate intermediate polytopes, so the test's δ should also be relaxed to about 1e-6.

Relation to darnstrom/daqp#160

This is what turns the PolyDAQP.jl check red on darnstrom/daqp#160, but the two are unrelated: that PR is a division-to-reciprocal rewrite in the DAQP LDP transform, and on arm64 it produces bit-identical bounding boxes to master across all 2000 seeds here. The over-approximation is present on unmodified master.

[cl445]

The correct filter is history-dedup

I kept digging and found the actual fix. The sum(Hcand) > level+1 cardinality test is the wrong pruning rule. The right one is to keep one constraint per unique history bitset and never drop by cardinality. A facet of the projection corresponds to a unique set of original constraints (its history), so two combinations with the same history are the same facet, and different histories are different or redundant facets that a cardinality bound must not drop.

Replacing the cardinality filter with history-dedup makes the projection exact on 4 of the 5 failing seeds I checked, down to machine epsilon:

seed	cardinality (current)	subset-adjacency	history-dedup
269	1.6e-7	1.6e-7	3.6e-15
1380	4.5e-7	4.5e-7	4.4e-15
1849	5.5e-5	5.5e-5	3.8e-15
1305	2.5e-4	5.9e-8	1.2e-14
1114	3.4e-7	3.4e-7	3.4e-7

So the macroscopic over-approximation is entirely the cardinality filter dropping real facets on degenerate inputs, and history-dedup removes it. (The subset-adjacency test I mentioned earlier only fixes the single worst case, confirming it is not the right criterion.)

The cost is constraint count, and it is tamable

History-dedup keeps every distinct-history constraint, and the current minrep leaves many geometrically redundant ones in. Its is_redundant[AS] .= 0 step (marking every active-set constraint non-redundant after an optimal solve) is unsound for completeness on the near-degenerate intermediate polytopes, so the peak count inflates to tens of thousands.

The count comes back down with a complete redundancy pass between levels (for each constraint, maximize its normal over the others, redundant iff the max is within tolerance of its bound). The design space I measured:

filter + minrep	max error	final count	peak
cardinality (current)	2.5e-4	~120	~250
history-dedup + current minrep	3.4e-7 (4/5 exact)	up to 805	up to 38000
history-dedup + complete LP minrep	~4.5e-7	~100	~5000

There is a tension between the two minrep variants: the current one is exact but leaves redundant constraints (bloated), while a complete LP pass is compact but occasionally drops a thin facet, reintroducing ~1e-7 errors. Solving the LP redundancy checks at a tight primal_tol (1e-12) recovers exactness on some seeds while staying compact, but not uniformly.

Residual floor

One seed (1114) stays at 3.4e-7 across every configuration, including LP solves at primal_tol = 1e-14. It comes from redundancy decisions on the near-degenerate intermediate polytopes, where a facet thinner than the LP resolution cannot be distinguished from redundant. Getting below this floor needs exact or rational arithmetic, or a different projection method such as double description, not tolerance tuning.

Suggested path

Replace the cardinality filter with history-dedup. That removes the macroscopic errors, which are the real bug. Then relax the elimination test's δ to about 1e-6 to match the residual LP floor. A projection exact to 1e-8 on degenerate inputs is a larger, separate effort.