Background: deciding support
From the exact cone test to linear programs: how the hold gate decides whether contacts can hold the box, proves its answer either way, and measures by how much.
What this explains. Before every edit, the hold gate asks whether a set of contacts can hold the box against gravity. This page derives the answer from first principles in two stages. First an exact test that works in the plane: it turns the question into membership in a cone in three dimensions and decides it with a witness, a certificate and a margin. Then the linear programs (LPs) that answer the same question in any dimension, with limits on the manipulator’s force and margins against a set of disturbances. The LPs are the gate today; the exact test is its cross-check. You need Coulomb friction and a little linear algebra. The code is in cmg/src/cmg/physics/planar_hold.py (cone_test) and hold_lp.py on the branch phase1/planar.
1. The question
A rigid box rests in the x–z plane on some contacts. Each contact pushes on the box with a force that Coulomb friction confines to a cone. The box is supported if some choice of admissible contact forces cancels gravity exactly: the net force and the net moment. That is statics only: no dynamics and no time.
2. From forces to one cone
Friction cones have two edges. At a contact with inward normal \(n\), tangent \(t\) and friction coefficient \(\mu\), the admissible forces are those with \(|f\cdot t| \le \mu\, f\cdot n\). In the plane that set is spanned by its two edge rays: \[ f = \lambda_+\,(n+\mu t) + \lambda_-\,(n-\mu t), \qquad \lambda_\pm \ge 0 . \]
A segment contact acts through its endpoints. A contact set can be a segment \([p_0,p_1]\) (an edge lying on the floor) with the same \(n\) and \(\mu\) along it. A force \(f\) applied at \(p = s\,p_0 + (1-s)\,p_1\) has the same force and moment as \(s f\) at \(p_0\) plus \((1-s) f\) at \(p_1\), because the moment \((p-c)\times f\) is linear in \(p\). So two points per contact suffice. A point contact is the case \(p_0 = p_1\).
Wrenches. Write each ray \(r\) at each endpoint \(p\) as a wrench: force and moment about the center of mass \(c\), made dimensionless, \[ w = \Big(r_x,\; r_z,\; \tfrac{(p-c)\times r}{\ell}\Big) \in \R^3, \qquad (p-c)\times r = (p-c)_z\, r_x - (p-c)_x\, r_z , \] with forces in units of the weight \(mg\) and moments in units of \(mg\,\ell\), where \(\ell = a/2\) is the box’s half-size. Gravity applies the wrench \((0,-1,0)\) about the center of mass, so the contacts must supply \[ \text{need} = (0,\, 1,\, 0). \]
The question becomes geometry. Let \(w_1,\dots,w_k\) be all the ray wrenches (at most 4 contacts × 2 endpoints × 2 rays = 16). Then \[ \boxed{\ \text{supported} \iff \text{need} \in C = \Big\{\textstyle\sum_i \lambda_i w_i : \lambda_i \ge 0\Big\}\ } \] where \(C\) is the convex cone the generators span.
Running example. A square box flat on the floor, \(\mu = 0.2\), touching at its two bottom corners A (left) and B (right). In units of \(\ell\), the corners sit at \((\mp1,-1)\) from the center of mass. The four generators are:
| Ray | Force \((f_x, f_z)\) | Moment / ℓ | Wrench \(w\) |
|---|---|---|---|
| A₊ | (+0.2, 1) | +0.8 | (0.2, 1, 0.8) |
| A₋ | (−0.2, 1) | +1.2 | (−0.2, 1, 1.2) |
| B₊ | (+0.2, 1) | −1.2 | (0.2, 1, −1.2) |
| B₋ | (−0.2, 1) | −0.8 | (−0.2, 1, −0.8) |
3. Seeing the cone
A cone in \(\R^3\) is hard to draw, but when every generator pushes upward (\(f_z > 0\)) it is fixed by where its rays cross the plane \(f_z = 1\). Each generator becomes a point \((f_x/f_z,\ \tau/f_z)\), the cone becomes the convex polygon of those points, and need becomes the origin. Supported then means the origin lies in the polygon.
scripts/bg_cone_test.py, which runs cone_test on both cases.The slice is only a picture. A disk pushing sideways can have \(f_z \le 0\), and every test below works in \(\R^3\).
Part I: the exact test. Sections 4–8 decide support by enumeration, in the plane.
4. “Yes”, with a witness: Carathéodory
Carathéodory’s theorem (conic form). If \(x \in \operatorname{cone}\{w_1,\dots,w_k\} \subset \R^d\), then \(x\) is a nonnegative combination of at most \(d\) linearly independent generators.
Why. Suppose \(x = \sum_i \lambda_i w_i\) with \(\lambda > 0\) on generators that are linearly dependent, so \(\sum_i \alpha_i w_i = 0\) for some \(\alpha \ne 0\), which we may take to have a positive entry. Replace \(\lambda\) by \(\lambda - s\alpha\) and increase \(s\) until the first coefficient hits zero. The sum still equals \(x\), the coefficients stay nonnegative, and one generator has dropped out. Repeat until the generators are independent; there are at most \(d\) of them.
Algorithm. Here \(d = 3\), so the test solves every subset of 1, 2 and 3 generators: \(16 + 120 + 560 = 696\) tiny systems (a projection, a 2×2 system, Cramer’s rule). A subset whose coefficients are \(\ge 0\) and that reproduces need is a witness: actual contact forces that hold the box.
Example. \(\tfrac12 w_{A_+} + \tfrac12 w_{B_-} = (0, 1, 0)\): half the weight through each corner, both forces leaning slightly inward (panels a, c).
5. “No”, with a certificate: Farkas
Farkas’s lemma (cone form). \(x \notin \operatorname{cone}\{w_i\}\) if and only if some \(y\) has \(y\cdot w_i \ge 0\) for every \(i\) and \(y \cdot x < 0\).
Geometrically, \(y\) is the normal of a plane through the origin with every generator on one side and need strictly on the other. Such a \(y\) is a certificate: anyone can check the \(k+1\) inequalities.
Which \(y\) to try. A full-dimensional polyhedral cone is the intersection of the half-spaces of its facets, and each facet contains at least two independent generators. So the normals \(\pm\, w_i \times w_j\) of the planes through two generators contain a certificate. Flat cones and single rays need a few more candidates (\(\pm w_i\), the negated part of need perpendicular to \(w_i\), normals of planes containing one generator, and \(-\text{need}\)): 385 candidates for 16 generators in all.
Example. On corner A alone the generators \(w_{A_\pm}\) span a flat cone. The plane through both has normal \(y = (1,-1,1)/\sqrt3\), with \(y\cdot w_{A_\pm} = 0\) and \(y\cdot \text{need} = -0.577 < 0\): certified unsupported (panel d). Physically, every force at A has a moment about the center of mass unless it points straight at it, which is outside A’s friction cone.
6. How much: the margin is a signed distance
The test also returns \[ m \;=\; \min_{y \text{ supporting},\ \|y\| = 1}\; y\cdot \text{need}, \] the minimum over every candidate \(y\) with all generators on its side.
Claim. \(m\) is the signed Euclidean distance from need to the boundary of \(C\): inside \(C\) the distance to the boundary, outside minus the distance to \(C\), and \(0\) when need lies on a flat cone.
Inside. \(C\) is the intersection of its facet half-spaces \(\{x : y_f\cdot x \ge 0\}\). A ball of radius \(\rho\) around need fits in \(C\) exactly when \(\rho \le y_f\cdot\text{need}\) for every facet, so the distance to the boundary is \(\min_f y_f\cdot\text{need}\); any other supporting plane bounds a larger half-space and gives a larger value.
Outside. The nearest point of \(C\) lies inside a two-dimensional face, on a generator’s ray, or at the apex. The unit vector from need toward it is then that face’s normal, the negated perpendicular part of need, or \(-\text{need}\); all are candidates and support \(C\). By duality, \(\operatorname{dist}(x, C) = \max\{-y\cdot x : y \text{ supporting},\ \|y\|\le1\}\).
Checked. On random cones the margin agrees with minus the NNLS distance outside (to 3·10⁻⁶) and with the shortest exit along hundreds of directions inside (to 4·10⁻⁷); this is a test in the suite.
The flat box has four facets:
| Facet | Through | Normal y | y · need | Meaning |
|---|---|---|---|---|
| Sliding (friction toward +x) | A₊, B₊ | (−0.981, 0.196, 0) | 0.196 | friction gives out |
| Sliding (friction toward −x) | A₋, B₋ | (0.981, 0.196, 0) | 0.196 | friction gives out |
| Tipping about A | A₊, A₋ | (−1, 1, −1)/√3 | 0.577 | B lifts off |
| Tipping about B | B₊, B₋ | (1, 1, 1)/√3 | 0.577 | A lifts off |
So \(m = \mu/\sqrt{1+\mu^2} = 0.196\), and sliding binds first. A purely sideways push of up to \(\mu\,mg = 0.2\,mg\) is absorbed; the worst direction also lifts the box a little, which unloads friction, hence 0.196. The tipping value does not depend on \(\mu\), so with \(\mu > 1/\sqrt2\) tipping binds: \(m = \min(\mu/\sqrt{1+\mu^2},\ 1/\sqrt3)\). For corner A alone, \(m = -0.636\): need lies 0.636 from the cone. The certificate of §5 gave −0.577: any certificate proves “no”, but only the nearest one measures the gap. The distance mixes forces and moments, so it depends on \(\ell\); the yes/no answer does not.
7. Robustness comes cheaply
- Friction. A smaller \(\mu\) narrows every cone, so for an interval of \(\mu\) the worst case is its lower end.
- Mass and center of mass enter only through need, which is linear in \((m,\ m c)\). Over a box of masses and centers the needs form a polytope, membership in a convex cone holds over a polytope iff it holds at its corners, and the margin’s worst case is also at a corner.
8. Where the exact test stops
The enumeration has two limits.
- Force limits. The disk’s push is bounded by its motors. With bounded forces the supportable wrenches are no longer a cone, and the candidate planes of §5 no longer describe them. Without limits, the disk’s unbounded push often makes \(C\) all of \(\R^3\), and the margin is infinite.
- Three dimensions. A spatial wrench has 6 components, so Carathéodory needs subsets of up to 6 generators. With 4 contacts × 4 points × 8 rays of a pyramidal friction cone (128 generators), that is \(\binom{128}{6} \approx 5.4\times10^{9}\) subsets.
Both limits disappear when the same question is posed as linear programs.
Part II: the linear programs. Sections 9–14 answer the same question in any dimension, with force limits and margins.
9. One backbone, any dimension
The LPs take four pieces of data:
| Data | Meaning |
|---|---|
| \(W \in \R^{d\times k}\) | generators: one column per friction-cone edge at each support point, as a wrench |
| \(\text{need} \in \R^d\) | the wrench the contacts must supply |
| \(B\,\lambda \le f_{\max}\) (optional) | force limits |
| \(u_1, \dots, u_J\) (optional) | disturbance directions for the margin |
The unknowns are the edge weights \(\lambda \ge 0\). Only the builder of \(W\) and need knows the dimension: planar wrenches have \(d = 3\) and exactly two edges per friction cone; spatial ones have \(d = 6\) and an inner pyramid of \(m\) edges \(n + \mu(\cos\theta_j t_1 + \sin\theta_j t_2)\) at the vertices of each contact polygon. The pyramid lies inside the true cone, so it can only reject wrongly, never accept, and 8 edges keep at least \(\cos(\pi/8) = 92\%\) of the friction in every direction.
10. LP-1: feasibility, with a witness and a certificate
“Find \(\lambda \ge 0\) with \(W\lambda = \text{need}\)” has no solution exactly when the answer is no, and some solvers then report nothing more. Adding slack that pays for any mismatch makes it always solvable: \[ \text{LP-1:}\qquad \delta \;=\; \min_{\lambda,\, s^+,\, s^- \,\ge\, 0}\; \mathbf 1^\top s^+ + \mathbf 1^\top s^- \quad\text{s.t.}\quad W\lambda + s^+ - s^- = \text{need}. \]
- \(\lambda = 0\) with \(s^\pm\) the positive and negative parts of need is feasible, and \(\delta \ge 0\).
- \(\delta = \min_{\lambda \ge 0}\|W\lambda - \text{need}\|_1\), the ℓ₁ distance from need to \(C\); \(\delta = 0\) means supported.
- Yes comes with a witness. The simplex method returns a basic optimal solution with at most \(d\) nonzero variables: Carathéodory’s theorem again, now produced by the solver.
Dual of LP-1. \(\displaystyle \delta = \max_{y}\ \text{need}^\top y \quad\text{s.t.}\quad W^\top y \le 0,\ \ -\mathbf 1 \le y \le \mathbf 1.\)
No comes with a certificate. If \(\delta > 0\), then \(z = -y\) has \(z\cdot w_i \ge 0\) for every generator and \(z\cdot\text{need} < 0\): the Farkas certificate of §5, at no extra cost.
Example. Both corners: \(\delta = 0\) with witness \(\lambda_{A_+} = \lambda_{B_-} = \tfrac12\), as in §4. Corner A alone: \(\delta = 1.0\) with dual \(y = (-1, 1, -1)\), so \(z = (1, -1, 1)\): the certificate of §5, unscaled.
11. LP-2: the margin, as exit distances
A Euclidean ball is not a polytope, so the margin of §6 is not an LP. A polytope of disturbances is. For a direction \(u\), the exit distance is how far need can move along \(u\) before leaving the supportable set: \[ \text{LP-2:}\qquad t(u) \;=\; \max_{\lambda \ge 0,\ 0 \le t \le T}\; t \quad\text{s.t.}\quad W\lambda = \text{need} + t\,u , \] capped at \(T = 10\) where the set never ends in that direction.
Claim (set inclusion). For a convex polytope \(D\) of disturbances, the largest \(r\) with \(\text{need} + rD\) supported is \(\min_j t(d_j)\) over the corners \(d_j\) of \(D\).
Why. \(\text{need} + rD\) is the convex hull of its corners, so a convex set contains it exactly when it contains the corners; along each corner’s direction the supportable points form an interval \([0, t(d_j)]\).
The gate uses \(D\) = the ℓ₁ unit ball, whose corners are the pure disturbances \(\pm e_i\) (a push along x, a lift, a twist). With LP-1 this gives a signed ℓ₁ margin: \(\min_i t(\pm e_i)\) inside and \(-\delta\) outside. It is always within a factor \(\sqrt3\) of the Euclidean margin of §6 (\(\sqrt6\) in 3D).
scripts/bg_hold_lp.py (HiGHS, float64).The flat box’s ℓ₁ margin is 0.2, exactly μ: the smallest pure push that breaks support. Its Euclidean margin was 0.196, because the worst Euclidean direction mixes a push with a slight lift.
12. LP-3: force limits
Every edge has unit normal component, so a contact’s normal force is the sum of its edge weights, and a limit is one more row. The gate uses the simulator’s motors instead: held still, the disk’s two slide motors (±1 N each, read from the MuJoCo model) must supply its push on the box plus its own weight, \(u = A(F - f_M)\), which gives two rows per motor. The feasible wrenches are then the environment’s cone plus a polytope from the disk, and directions the disk helps have finite exits: in panel (b), \(t(+f_x) = 0.72 = 0.5\) from the disk plus \(0.2 \times 1.1\) of floor friction, since the disk’s friction also presses the box onto the floor. The dual gains one multiplier per row; the certificate stays a few inequalities anyone can check.
13. Solving it inside JAX
The LPs are tiny and dense (\(d = 3\) or \(6\) rows, tens to a hundred-odd columns), many are solved at once, and the interesting ones sit on the support boundary, where they are degenerate.
| Method | Returns | Fit |
|---|---|---|
| Simplex | a vertex (witness with ≤ d generators) and the dual (certificate) | best: a handful of pivots, exact up to rounding |
| Interior point | an interior point; duals | workable; dense witness, delicate certificate near the boundary |
| First-order (PDHG, e.g. MPAX) | an approximate point | poor: slow near the boundary, and an early stop cannot certify “no”. MPAX falsely rejected about 193 of 416 supported near-boundary cases |
linrax (a two-phase tableau simplex in JAX, pinned to 0.3.0) is the solver; HiGHS (via SciPy, float64) checks every answer in the tests. hold_lp replaces linrax’s pivot loop at import with three changes: a cap of 300 pivots per phase (under vmap one LP that never finishes hangs the batch), Bland’s leaving rule (linrax’s rule let a 3D LP cycle), and a corrected test for artificial rows (which had returned a slack of −1 as “optimal” on the flat box). linrax returns no dual, so the dual is recovered from the final basis by solving \(B^\top y = c_B\).
Float32. All computations stay in float32. Four measures, each fixing a failure the tests found:
- the equality rows are rotated by a Householder reflection that makes every right-hand side positive (the planar need \((0,1,0)\) has two zero rows);
- the primal is recomputed from the final basis rather than read off the tableau;
- the pivot tolerance is \(10^{-4}\) (at linrax’s \(10^{-6}\) it pivoted on noise: 19 wrong answers on 600 random LPs, none at \(10^{-4}\));
- exits are capped at \(T = 10\), not \(10^3\) (the large cap’s coefficient fell to the pivot threshold and produced negative weights).
Every answer is then checked independently, with elementwise products (no reduced-precision TF32 on GPUs): a witness by its residual, a certificate by its signs, an exit by its weights. A failed exit check sets that exit to 0, so a solver error can shrink the margin but not enlarge it.
Benchmark, before force limits (record at commit e29d96f; 1,000 random cases and 416 within 2 mm of disk height of the boundary; HiGHS reference; float32):
| Device | Exact test, near boundary | LP, near boundary | ms per 256: exact / LP-1 / LP with margins |
|---|---|---|---|
| CPU | 390 of 416 | 386 of 416 | 4.1 / 2.7 / 6.3 |
| GPU | 389 of 416 | 383 of 416 | 1.5 / 2.0 / 2.5 |
Both get all 1,000 random cases right; every other disagreement is a pose reported uncertain and rejected (the safe direction), where a witness needs weights of 60–190 times the box’s weight, beyond what float32 can check. The gate’s own benchmark, with motor limits, is on the hold gate.
In 3D. On a cube touched by up to four contacts (8-edge pyramids, \(d = 6\), up to 128 generators, 1,024 cases), nothing hangs, nothing is certified wrongly, 763 of 768 random and 255 of 256 near-boundary cases are right, and the rest are “uncertain” (safe) except one false acceptance 0.1% outside the boundary. The uncertain cases come from float32 rounding over the longer pivot sequences; the cure, rebuilding the tableau from the final basis (reinversion), is planned for the 3D build.
14. How the gate uses them
- The gate is LP-1 and LP-2 with the motor rows, on the edit’s less-constrained end, at χ₀ and at the transition sweep’s poses (the hold gate). It accepts at margin ≥ 0.
- The margin is set inclusion with the ℓ₁ ball; a scene can pass its own polytope \(D\).
- The exact test (§4–6) stays as a cone-only cross-check.
| Exact test (plane only) | LPs (any dimension) | |
|---|---|---|
| Yes | witness, ≤ 3 generators | witness, ≤ d generators (basic solution) |
| No | Farkas plane from a candidate list | Farkas plane from the dual |
| Margin | signed Euclidean distance | signed ℓ₁ distance, or set inclusion of any polytope |
| Force limits | no | yes (rows) |
| 3D | no (about 5×10⁹ subsets) | yes (6 rows) |
Check yourself
\(0.5/\sqrt{1.25} = 0.447 < 1/\sqrt3 = 0.577\): sliding still binds. cone_test returns 0.4472.
No. Every wrench a force at A can produce lies in one plane, \(\tau/\ell = f_z - f_x\) (arm \((-1,-1)\)), whatever \(\mu\) is. Need \(= (0,1,0)\) gives \(1 \ne 0\), so it is off that plane. A single point contact holds a body only if the vertical through the center of mass passes through that point and the vertical lies in its friction cone.
A plane through the origin in \(\R^3\) is fixed by two independent directions. Every facet of a 3D cone is such a plane and contains at least two generators; three on one plane add nothing new.
\(\lambda = 0\) with \(s^\pm\) equal to the positive and negative parts of need satisfies every constraint, and the objective is at least 0. So the LP is feasible and bounded, and the solver always ends with an optimal basis and its dual.
No. The pyramid lies inside the true friction cone, so every force it allows is a real force. The error is one-sided: it may reject a pose that needs more than 92% of the available friction in an unlucky direction.
The margin is the smallest exit. The disk can only push in +x, so it cannot help resist a disturbance that needs −x force from the contacts; floor friction alone sets \(t(-f_x) = 0.2\). A second disk on the other face would.
Reproduce
cd ~/repos/research/c-gnn && git switch phase1/planar && git pull && cd cmg
JAX_PLATFORMS=cpu .venv314/bin/python <website>/research/nonprehensile/planar-build/scripts/bg_cone_test.py <out dir> # §2–6, cone figure
JAX_PLATFORMS=cpu .venv314/bin/python <website>/research/nonprehensile/planar-build/scripts/bg_hold_lp.py <out dir> # §10–12, exit figure
# the benchmark of §13 is a record: bg_hold_lp_bench.py runs at commit e29d96f