Introduction
Linear programming is a foundational tool of operations research, and for decades the simplex and interior-point methods have dominated the space of LP solvers [10 ] . Both rely on matrix factorizations whose memory traffic and fill-in scale poorly on modern massively parallel hardware. First-order methods, which trade factorizations for matrix-vector products, have recently emerged as a practical alternative for large-scale LP [34 , 41 ] . PDLP, a development of PDHG [12 ] , demonstrated that a carefully engineered first-order method is competitive on broad LP benchmarks [2 ] . Its computational core is matrix-vector multiplication, so PDHG suits GPUs. cuPDLP.jl [37 ] removed the remaining sequential components of PDLP and showed that GPU implementations can achieve near-universal performance on par with mature CPU solvers, while cuPDLP-C [39 ] provided a mature C implementation of this approach. D-PDLP [32 ] scaled PDLP to multi-GPU systems, while HPR-LP [13 ] pursued an alternative first-order approach based on the Halpern Peaceman–Rachford method. PDHG-type LP solvers now ship in production codes, including Gurobi’s PDHG [21 ] , NVIDIA’s cuOpt [43 ] , COPT’s PDLP [18 ] , HiGHS’s bundled CuPDLP-C [39 , 26 ] , FICO Xpress’s PDHG [16 ] , and Knitro’s PDLP [5 ] .
This paper presents a lever that the LP community has largely written off: reformulation . Folklore among practitioners holds that reformulating an LP offers limited benefit to interior-point methods. Their iteration counts are remarkably stable, nearly independent of instance size or data. Formulation affects performance mainly through the sparsity of the factorized systems. Customized linear solve and automatic presolve already capture most of that benefit [20 , 49 , 10 ] . Under this folklore there is little reason for a modeler to hand-craft an equivalent formulation.
PDHG-type methods are different in two ways. First, their per-iteration cost is the cost of multiplying by the constraint matrix, so a reformulation that shrinks or restructures that matrix changes the price of every iteration directly. Second, their convergence speed is governed by the geometry of the formulation. Restarted PDHG converges linearly at a rate set by sharpness constants of the primal–dual formulation [3 , 38 ] . Those constants are closely tied to Hoffman constants and level-set geometry [51 , 50 ] , and equivalent reformulations change them. Reformulation therefore has two levers on a PDHG-type solver where it has almost none on an interior-point method.
We evaluate reformulation by absorbing constraints into projected sets. PDLP projects onto boxes; PDHG needs only closed convex feasible sets with computable projections. So any subset of the constraints whose induced feasible set has an efficient Euclidean projection can be absorbed , moved out of the constraint matrix and into the projection step, instead of dualized. The solver then runs on a different saddle-point problem for the same LP, with a smaller coupling matrix, fewer dual variables, and different geometry.
Our contributions are summarized as follows.
A systematic framework for absorption. We formalize the absorption procedure and identify classes of polyhedra that admit efficient projections (§ 2 ). We also develop a construction for problems involving auxiliary variables that select maxima over groups (§ 2.4 ), which serves as the basis for our applications in § 5 .
Algorithmic improvements and convergence guarantees. We analyze how absorption affects four key quantities: the coupling norm, the termination certificate, the distance scale, and the sharpness modulus, as well as its potential to reduce per-iteration cost (§ 4 ). We establish that sharpness holds automatically for polyhedral absorbed formulations and show how a lattice structure in the vertex values of both value polytopes yields an explicit sharpness bound (Corollary 1), enabling quantitative convergence guarantees.
Two case studies. We apply the framework to the QUBO relaxation (§ 5.1 ) and the joint replenishment LP relaxation (§ 5.3 ), and implement each reformulation as a GPU solver: PRAQ for QUBO and PRIA for joint replenishment. For QUBO, half-integrality on both sides yields explicit constants and a linear convergence rate (Theorem 1). For joint replenishment we establish a sublinear rate whose only data dependence is the instance size (Theorem 2). We further evaluate both reformulations against cuPDLPx and Gurobi on the original formulations.
Code. PRAQ, PRIA, the baselines they are measured against and the benchmark harness are at https://github.com/MIT-Lu-Lab/reformulation.
Notation. For a positive integer n n n [ n ] : = { 1 , … , n } [n] := \{1,\ldots,n\} [ n ] := { 1 , … , n } S S S proj S \operatorname{proj}_{S} proj S S S S M M M ∥ v ∥ M : = v T M v \lVert v\rVert_{M} := \sqrt{v^{\mathsf{T}}Mv} ∥ v ∥ M := v T M v diam M ( S ) \operatorname{diam}_{M}(S) diam M ( S ) dist M ( ⋅ , S ) \operatorname{dist}_{M}(\cdot,S) dist M ( ⋅ , S ) S S S ∥ ⋅ ∥ M \lVert\cdot\rVert_{M} ∥ ⋅ ∥ M vec ( X ) \operatorname{vec}(X) vec ( X ) X X X ⟨ ⋅ , ⋅ ⟩ \langle\cdot,\cdot\rangle ⟨ ⋅ , ⋅ ⟩ 1 \mathbf{1} 1 Δ r : = { v ≥ 0 : 1 T v = r } \Delta^{r} := \{v \ge0 : \mathbf{1}^{\mathsf{T}}v = r\} Δ r := { v ≥ 0 : 1 T v = r }
Reformulation for PDHG-type methods
PDHG for linear programming
PDHG solves convex–concave saddle-point problems
min x ∈ X max y ∈ Y L ( x , y ) : = ⟨ c , x ⟩ + ⟨ y , q − A x ⟩ (1) \min_{x\in\mathcal{X}}\max_{y\in\mathcal{Y}} L(x,y) := \langle c,x\rangle+ \langle y,q-Ax\rangle \tag*{(1)} x ∈ X min y ∈ Y max L ( x , y ) := ⟨ c , x ⟩ + ⟨ y , q − A x ⟩ (1) with X ⊆ R n \mathcal{X}\subseteq\mathbb{R}^{n} X ⊆ R n Y ⊆ R m \mathcal{Y}\subseteq\mathbb{R}^{m} Y ⊆ R m [12 ] under positive diagonal scalings Λ X , Λ Y \Lambda_{\mathcal{X}},\Lambda_{\mathcal{Y}} Λ X , Λ Y
x k + 1 : = proj X ( x k − τ Λ X − 1 ( c − A T y k ) ) , y k + 1 : = proj Y ( y k + σ Λ Y − 1 ( q − A ( 2 x k + 1 − x k ) ) ) (2) x^{k+1}:=\operatorname{proj}_{\mathcal{X}}\left(x^{k}-\tau\Lambda_{\mathcal{X}}^{-1}\left(c-A^{\mathsf{T}}y^{k}\right)\right),\qquad y^{k+1}:=\operatorname{proj}_{\mathcal{Y}}\left(y^{k}+\sigma\Lambda_{\mathcal{Y}}^{-1}\left(q-A(2x^{k+1}-x^{k})\right)\right) \tag*{(2)} x k + 1 := proj X ( x k − τ Λ X − 1 ( c − A T y k ) ) , y k + 1 := proj Y ( y k + σ Λ Y − 1 ( q − A ( 2 x k + 1 − x k ) ) ) (2) which converges whenever τ σ ∥ A ∥ Λ 2 < 1 \tau\sigma\lVert A\rVert_{\Lambda}^{2}<1 τ σ ∥ A ∥ Λ 2 < 1 effective coupling norm
∥ A ∥ Λ : = ∥ Λ Y − 1 / 2 A Λ X − 1 / 2 ∥ 2 . (3) \lVert A\rVert_{\Lambda}:=\left\lVert\Lambda_{\mathcal{Y}}^{-1/2}A\Lambda_{\mathcal{X}}^{-1/2}\right\rVert_{2}. \tag*{(3)} ∥ A ∥ Λ := Λ Y − 1/2 A Λ X − 1/2 2 . (3) Both projections are taken in ∥ ⋅ ∥ Λ X \lVert\cdot\rVert_{\Lambda_{\mathcal{X}}} ∥ ⋅ ∥ Λ X ∥ ⋅ ∥ Λ Y \lVert\cdot\rVert_{\Lambda_{\mathcal{Y}}} ∥ ⋅ ∥ Λ Y (2) a proximal-point step [24 ] . Λ X \Lambda_{\mathcal{X}} Λ X ∥ A ∥ Λ \lVert A\rVert_{\Lambda} ∥ A ∥ Λ Λ X = Λ Y = I \Lambda_{\mathcal{X}}=\Lambda_{\mathcal{Y}}=I Λ X = Λ Y = I ∥ A ∥ Λ = ∥ A ∥ 2 \lVert A\rVert_{\Lambda}=\lVert A\rVert_{2} ∥ A ∥ Λ = ∥ A ∥ 2 τ = η / ω \tau=\eta/\omega τ = η / ω σ = η ω \sigma=\eta\omega σ = η ω step size η \eta η primal weight ω \omega ω (2) by T ω T_{\omega} T ω
A general LP fits (1) by taking X \mathcal{X} X Y \mathcal{Y} Y A A A standard form , the formulation solved by PDLP and its successors [2 , 35 ] , and each of its iterations costs two sparse matrix-vector products plus two trivial projections.
Absorbing constraints into the projection
PDLP projects onto boxes, but nothing in (2) requires that. Beyond the closed convexity of X \mathcal{X} X Y \mathcal{Y} Y proj X \operatorname{proj}_{\mathcal{X}} proj X proj Y \operatorname{proj}_{\mathcal{Y}} proj Y
Consider an LP with the constraint system
min c T x s.t. A x ≥ q , x ∈ box (4) \min c^{\mathsf{T}}x\quad\text{s.t. }Ax\geq q,\ x\in\operatorname{box} \tag*{(4)} min c T x s.t. A x ≥ q , x ∈ box (4) A row of (4) can be dualized , in which case PDHG carries a multiplier for it, or absorbed , in which case proj X \operatorname{proj}_{\mathcal{X}} proj X
Definition 1 (Absorption split). An absorption split of (4) is a partition of its rows into an absorbed set E \mathcal{E} E and a dualized set R \mathcal{R} R such that
X : = { x ∈ box : A E x ≥ q E } \mathcal{X}:=\{x\in\operatorname{box}:A_{\mathcal{E}}x\geq q_{\mathcal{E}}\} X := { x ∈ box : A E x ≥ q E } has a computable projection.
Under the split, (4) reads
min c T x s.t. A E x ≥ q E , A R x ≥ q R , x ∈ box , (5) \min c^{\mathsf{T}}x\quad\text{s.t. }A_{\mathcal{E}}x\geq q_{\mathcal{E}},\ A_{\mathcal{R}}x\geq q_{\mathcal{R}},\ x\in\operatorname{box}, \tag*{(5)} min c T x s.t. A E x ≥ q E , A R x ≥ q R , x ∈ box , (5) where proj X \operatorname{proj}_{\mathcal{X}} proj X (1) takes A : = A R A:=A_{\mathcal{R}} A := A R q : = q R q:=q_{\mathcal{R}} q := q R Y \mathcal{Y} Y R \mathcal{R} R
min x ∈ X max y ≥ 0 ⟨ c , x ⟩ + ⟨ y , q R − A R x ⟩ . (6) \min_{x\in\mathcal{X}}\max_{y\geq0}\ \langle c,x\rangle+\langle y,q_{\mathcal{R}}-A_{\mathcal{R}}x\rangle. \tag*{(6)} x ∈ X min y ≥ 0 max ⟨ c , x ⟩ + ⟨ y , q R − A R x ⟩ . (6) The standard form is the split E = ∅ \mathcal{E}=\varnothing E = ∅ § 2.4 . We propose to reformulate structured LPs using this absorption split.
A useful subset of polyhedra with fast projections
Table 1 presents some examples of practical polyhedral structures with efficient Euclidean projections.
Table 1. Candidate structures for an absorption split, with the constraint family each one absorbs. The projection column is the best worst-case cost over n n n
Each row of Table 1 is a family of constraints together with the set it becomes. Absorption into X \mathcal{X} X (6) .
Assignment and partition rows , ∑ i ∈ G g x i = r g \sum_{i\in G_g} x_i = r_g ∑ i ∈ G g x i = r g x ≥ 0 x \ge0 x ≥ 0 g g g Δ r g \Delta^{r_g} Δ r g [11 , 44 ] ; sort-and-threshold [25 ] costs a logarithmic factor more, and expected linear-time pivoting [15 ] and the variants of [14 ] are the fast choices in practice. These rows back assignment, transportation, covering-by-choice and lot-sizing models.
A weighted partition , a ⊤ x = r a^\top x = r a ⊤ x = r a > 0 a > 0 a > 0 x ≥ 0 x \ge0 x ≥ 0 X = { x ≥ 0 : a ⊤ x = r } \mathcal{X} = \{x \ge0 : a^\top x = r\} X = { x ≥ 0 : a ⊤ x = r } [45 ] .
A sparsity budget , ∥ x ∥ 1 ≤ r \lVert x\rVert_1 \le r ∥ x ∥ 1 ≤ r X = { x : ∥ x ∥ 1 ≤ r } \mathcal{X} = \{x : \lVert x\rVert_1 \le r\} X = { x : ∥ x ∥ 1 ≤ r }
Monotonicity , x 1 ≤ ⋯ ≤ x n x_1 \le\cdots\le x_n x 1 ≤ ⋯ ≤ x n [8 ] .
Max-selection rows , z g ≥ x j z_g \ge x_j z g ≥ x j G g G_g G g z g z_g z g Y \mathcal{Y} Y λ g ≥ 0 \lambda_g \ge0 λ g ≥ 0 z g z_g z g ∑ j ∈ G g λ g , j = c g \sum_{j\in G_g} \lambda_{g,j} = c_g ∑ j ∈ G g λ g , j = c g z z z § 2.4 carries this out.
Max-selection
Max-selection is common in practical LPs: auxiliary variables that select maxima over groups of the x x x
min x ∈ X 0 , z ⟨ h , x ⟩ + ∑ g ∈ G c g z g s.t. A R x ≥ q R , z g ≥ x j , ∀ j ∈ G g , ∀ g ∈ G , (7) \min_{x\in\mathcal{X}_0,\ z}\ \langle h,x\rangle+ \sum_{g\in G} c_g z_g
\quad\text{s.t.}\quad
A_Rx \ge q_R,\qquad
z_g \ge x_j,\ \forall j\in G_g,\ \forall g\in G, \tag*{(7)} x ∈ X 0 , z min ⟨ h , x ⟩ + g ∈ G ∑ c g z g s.t. A R x ≥ q R , z g ≥ x j , ∀ j ∈ G g , ∀ g ∈ G , (7) where each group G g ⊆ [ n ] G_g \subseteq[n] G g ⊆ [ n ] c g ≥ 0 c_g \ge0 c g ≥ 0 A R x ≥ q R A_Rx \ge q_R A R x ≥ q R
Reformulate (7) as a saddle problem.
min x ∈ X 0 max λ g ∈ Δ c g , g ∈ G y R ≥ 0 ⟨ h , x ⟩ + ∑ g ∈ G λ g ⊤ x G g + ⟨ y R , q R − A R x ⟩ (8) \min_{x\in\mathcal{X}_0}\ \max_{\substack{\lambda_g\in\Delta^{c_g},\ g\in G\\ y_R\ge0}}
\left\langle h,x\right\rangle+ \sum_{g\in G}\lambda_g^\top x_{G_g} + \left\langle y_R,\ q_R-A_Rx\right\rangle \tag*{(8)} x ∈ X 0 min λ g ∈ Δ c g , g ∈ G y R ≥ 0 max ⟨ h , x ⟩ + g ∈ G ∑ λ g ⊤ x G g + ⟨ y R , q R − A R x ⟩ (8) which is (1) with X = X 0 \mathcal{X}=\mathcal{X}_0 X = X 0 z z z
The dual of a max-selection family is an assignment family. Each primal variable contributes one row to the dual. Price z g ≥ x j z_g \ge x_j z g ≥ x j λ g , j ≥ 0 \lambda_{g,j} \ge0 λ g , j ≥ 0 A R x ≥ q R A_{\mathcal{R}}x \ge q_{\mathcal{R}} A R x ≥ q R y R ≥ 0 y_{\mathcal{R}} \ge0 y R ≥ 0
L ( x , z ; λ , y R ) = ⟨ h , x ⟩ + ∑ g ∈ G z g ( c g − ∑ j ∈ G g λ g , j ) + ∑ g ∈ G λ g T x G g + ⟨ y R , q R − A R x ⟩ . L(x,z;\lambda,y_{\mathcal{R}})=\langle h,x\rangle+\sum_{g\in G}z_g\left(c_g-\sum_{j\in G_g}\lambda_{g,j}\right)+\sum_{g\in G}\lambda_g^{\mathsf T}x_{G_g}+\langle y_{\mathcal{R}},q_{\mathcal{R}}-A_{\mathcal{R}}x\rangle. L ( x , z ; λ , y R ) = ⟨ h , x ⟩ + g ∈ G ∑ z g c g − j ∈ G g ∑ λ g , j + g ∈ G ∑ λ g T x G g + ⟨ y R , q R − A R x ⟩ . Each z g z_g z g
∑ j ∈ G g λ g , j = c g , λ g ≥ 0. \sum_{j\in G_g}\lambda_{g,j}=c_g,\qquad\lambda_g\ge0. j ∈ G g ∑ λ g , j = c g , λ g ≥ 0. one per group, over disjoint groups. That is the assignment family of Table 1 , with simplices of total c g c_g c g
So absorb it. Definition 1 applied to those rows puts each λ g \lambda_g λ g Δ c g \Delta^{c_g} Δ c g z g z_g z g (7) , and what is left is (8) .
When R = ∅ \mathcal{R}=\varnothing R = ∅ X 0 \mathcal{X}_0 X 0 (7) and of its dual has been absorbed, and we call the split complete . Both case studies are complete splits.
Restarted Halpern PDHG and theoretical guarantees
This section presents the restarted reflected Halpern algorithm the implementations run [38 ] . Throughout, (1) is a feasible and bounded LP presented under an absorption split, with X ⊆ R n \mathcal{X}\subseteq\mathbb{R}^n X ⊆ R n Y ⊆ R m \mathcal{Y}\subseteq\mathbb{R}^m Y ⊆ R m Λ X , Λ Y \Lambda_{\mathcal{X}},\Lambda_{\mathcal{Y}} Λ X , Λ Y (2) come from the split. Write Z : = X × Y \mathcal{Z}:=\mathcal{X}\times\mathcal{Y} Z := X × Y v = ( x , y ) v=(x,y) v = ( x , y ) Λ : = diag ( Λ X , Λ Y ) \Lambda:=\operatorname{diag}(\Lambda_{\mathcal{X}},\Lambda_{\mathcal{Y}}) Λ := diag ( Λ X , Λ Y ) Z ⋆ \mathcal{Z}^{\star} Z ⋆ L L L Z \mathcal{Z} Z Z ⋆ = X ⋆ × Y ⋆ \mathcal{Z}^{\star}=\mathcal{X}^{\star}\times\mathcal{Y}^{\star} Z ⋆ = X ⋆ × Y ⋆ X ⋆ , Y ⋆ \mathcal{X}^{\star},\mathcal{Y}^{\star} X ⋆ , Y ⋆ L P ( x ) : = max y ∈ Y L ( x , y ) L_{\mathrm{P}}(x):=\max_{y\in\mathcal{Y}}L(x,y) L P ( x ) := max y ∈ Y L ( x , y ) L D ( y ) : = min x ∈ X L ( x , y ) L_{\mathrm{D}}(y):=\min_{x\in\mathcal{X}}L(x,y) L D ( y ) := min x ∈ X L ( x , y ) gap ( v ) : = L P ( x ) − L D ( y ) \operatorname{gap}(v):=L_{\mathrm{P}}(x)-L_{\mathrm{D}}(y) gap ( v ) := L P ( x ) − L D ( y )
The algorithm
One step of PDHG (2) is a proximal-point step in the canonical norm
P : = [ ω η Λ X A T A 1 η ω Λ Y ] , ∥ v ∥ P : = v T P v . (9) P:=\begin{bmatrix}\frac{\omega}{\eta}\Lambda_{\mathcal{X}} & A^{\mathsf T}\\ A & \frac{1}{\eta\omega}\Lambda_{\mathcal{Y}}\end{bmatrix},\qquad\lVert v\rVert_P:=\sqrt{v^{\mathsf T}Pv}. \tag*{(9)} P := [ η ω Λ X A A T η ω 1 Λ Y ] , ∥ v ∥ P := v T P v . (9) taken on the saddle operator F ( v ) : = ( c − A T y + N X ( x ) , A x − q + N Y ( y ) ) F(v):=(c-A^{\mathsf T}y+N_{\mathcal{X}}(x),Ax-q+N_{\mathcal{Y}}(y)) F ( v ) := ( c − A T y + N X ( x ) , A x − q + N Y ( y )) N S N_{\mathcal{S}} N S F F F F − 1 ( 0 ) = Z ⋆ F^{-1}(0)=\mathcal{Z}^{\star} F − 1 ( 0 ) = Z ⋆ T ω T_\omega T ω F F F [24 ] . It is firmly nonexpansive in ∥ ⋅ ∥ P \lVert\cdot\rVert_P ∥ ⋅ ∥ P 2 T ω − I 2T_\omega-I 2 T ω − I [7 Proposition 4.4] . Conjugating by Λ − 1 / 2 \Lambda^{-1/2} Λ − 1/2
Λ − 1 / 2 P Λ − 1 / 2 = [ ω η I Λ X − 1 / 2 A T Λ Y − 1 / 2 Λ Y − 1 / 2 A Λ X − 1 / 2 1 η ω I ] . \Lambda^{-1/2}P\Lambda^{-1/2}=
\begin{bmatrix}
\frac{\omega}{\eta}I & \Lambda_{\mathcal{X}}^{-1/2}A^{\mathsf T}\Lambda_{\mathcal{Y}}^{-1/2}\\
\Lambda_{\mathcal{Y}}^{-1/2}A\Lambda_{\mathcal{X}}^{-1/2} & \frac{1}{\eta\omega}I
\end{bmatrix}. Λ − 1/2 P Λ − 1/2 = [ η ω I Λ Y − 1/2 A Λ X − 1/2 Λ X − 1/2 A T Λ Y − 1/2 η ω 1 I ] . so P ≻ 0 P\succ0 P ≻ 0 η < 1 / ∥ A ∥ Λ \eta<1/\lVert A\rVert_{\Lambda} η < 1/ ∥ A ∥ Λ ω > 0 \omega>0 ω > 0
∥ v ∥ P ≤ ω + ω − 1 η ∥ v ∥ Λ , and at ω = 1 , ∥ v ∥ P < 2 η ∥ v ∥ Λ . (10) \lVert v\rVert_P\le\sqrt{\frac{\omega+\omega^{-1}}{\eta}}\lVert v\rVert_{\Lambda},\qquad\text{and at }\omega=1,\quad\lVert v\rVert_P<\sqrt{\frac{2}{\eta}}\lVert v\rVert_{\Lambda}. \tag*{(10)} ∥ v ∥ P ≤ η ω + ω − 1 ∥ v ∥ Λ , and at ω = 1 , ∥ v ∥ P < η 2 ∥ v ∥ Λ . (10) Halpern iteration. Given an anchor v a v^{a} v a [22 ] applied to T ω T_{\omega} T ω
v k + 1 = k + 1 k + 2 T ω ( v k ) + 1 k + 2 v a . (11) v^{k+1} = \frac{k+1}{k+2}T_{\omega}(v^{k}) + \frac{1}{k+2}v^{a}. \tag*{(11)} v k + 1 = k + 2 k + 1 T ω ( v k ) + k + 2 1 v a . (11) Since T ω T_{\omega} T ω ∥ ⋅ ∥ P \lVert\cdot\rVert_{P} ∥ ⋅ ∥ P [33 ]
∥ v k − T ω ( v k ) ∥ P ≤ 2 dist P ( v a , Z ⋆ ) k + 1 . \lVert v^{k}-T_{\omega}(v^{k})\rVert_{P} \le\frac{2\operatorname{dist}_{P}(v^{a},Z^{\star})}{k+1}. ∥ v k − T ω ( v k ) ∥ P ≤ k + 1 2 dist P ( v a , Z ⋆ ) . Reflection. The reflection 2 T ω − I 2T_{\omega}-I 2 T ω − I T ω T_{\omega} T ω [7 ]
v k + 1 = k + 1 k + 2 ( 2 T ω ( v k ) − v k ) + 1 k + 2 v a , (12) v^{k+1} = \frac{k+1}{k+2}\left(2T_{\omega}(v^{k})-v^{k}\right) + \frac{1}{k+2}v^{a}, \tag*{(12)} v k + 1 = k + 2 k + 1 ( 2 T ω ( v k ) − v k ) + k + 2 1 v a , (12) which halves the constant of (11) [38 Corollary 1] .
Restarting. The restarted reflected Halpern iteration [38 , 35 ] resets the anchor when a restart condition is met. Let j j j k k k j j j
v j , k + 1 = k + 1 k + 2 ( 2 T ω ( v j , k ) − v j , k ) + 1 k + 2 v j , 0 , (13) v^{j,k+1} = \frac{k+1}{k+2}\left(2T_{\omega}(v^{j,k})-v^{j,k}\right) + \frac{1}{k+2}v^{j,0}, \tag*{(13)} v j , k + 1 = k + 2 k + 1 ( 2 T ω ( v j , k ) − v j , k ) + k + 2 1 v j , 0 , (13) and when a restart is triggered after inner iteration k j k_j k j
v j + 1 , 0 : = T ω ( v j , k j ) . (14) v^{j+1,0} := T_{\omega}(v^{j,k_j}). \tag*{(14)} v j + 1 , 0 := T ω ( v j , k j ) . (14) Write v a : = v j , 0 v^{a}:=v^{j,0} v a := v j , 0 w : = T ω ( v ) w:=T_{\omega}(v) w := T ω ( v ) Z Z Z [35 ] . Sufficient decay fires at ∥ v − T ω ( v ) ∥ P ≤ β s u f f ∥ v a − T ω ( v a ) ∥ P \lVert v-T_{\omega}(v)\rVert_{P}\le\beta_{\mathrm{suff}}\lVert v^{a}-T_{\omega}(v^{a})\rVert_{P} ∥ v − T ω ( v ) ∥ P ≤ β suff ∥ v a − T ω ( v a ) ∥ P Necessary decay with stall fires at ∥ v − T ω ( v ) ∥ P ≤ β n e c ∥ v a − T ω ( v a ) ∥ P \lVert v-T_{\omega}(v)\rVert_{P}\le\beta_{\mathrm{nec}}\lVert v^{a}-T_{\omega}(v^{a})\rVert_{P} ∥ v − T ω ( v ) ∥ P ≤ β nec ∥ v a − T ω ( v a ) ∥ P ∥ v − T ω ( v ) ∥ P \lVert v-T_{\omega}(v)\rVert_{P} ∥ v − T ω ( v ) ∥ P Artificial fires once the epoch has run for a β a r t \beta_{\mathrm{art}} β art (14) is taken at the updated primal weight, where nonexpansiveness holds its residual at or below that of v j , k j v^{j,k_j} v j , k j
The primal weight is updated at restarts only, toward the ratio minimizing ∥ ⋅ ∥ P \lVert\cdot\rVert_{P} ∥ ⋅ ∥ P
ln ω ← ( 1 − θ ) ln ω + θ ln ∥ y − y a ∥ Λ Y ∥ x − x a ∥ Λ X , θ ← θ t g t + θ d e c ( θ − θ t g t ) , (15) \ln\omega\leftarrow(1-\theta)\ln\omega+\theta\ln\frac{\lVert y-y^{a}\rVert_{\Lambda_Y}}{\lVert x-x^{a}\rVert_{\Lambda_X}},\qquad\theta\leftarrow\theta_{\mathrm{tgt}}+\theta_{\mathrm{dec}}(\theta-\theta_{\mathrm{tgt}}), \tag*{(15)} ln ω ← ( 1 − θ ) ln ω + θ ln ∥ x − x a ∥ Λ X ∥ y − y a ∥ Λ Y , θ ← θ tgt + θ dec ( θ − θ tgt ) , (15) skipping the update if either norm vanishes. Early restarts therefore adapt ω \omega ω
Computational implementation. The method is (13) applied to the absorbed saddle problem (6) , run at the parameters of Table 14 in Appendix § 10 and started at
Table 14. Default parameters of the method of § 3.1 for the two instantiations, tuned by randomized search with Optuna [1 ] on instances disjoint from the benchmark sets. The step sizes are the two coupling norms of § 5.1 and § 5.3 . The sufficient restart parameters are the analytic optima for the ∥ ⋅ ∥ P \|\cdot\|_P ∥ ⋅ ∥ P
v 0 , 0 = ( proj X ( 0 ) , proj Y ( 0 ) ) . v^{0,0} = \left(\operatorname{proj}_{X}(0),\operatorname{proj}_{Y}(0)\right). v 0 , 0 = ( proj X ( 0 ) , proj Y ( 0 ) ) . The instance data is scaled to unit norm, the scale at which Table 14 is tuned, and unscaled at termination. Every ℓ \ell ℓ κ \kappa κ (4) at T ω ( v ) T_{\omega}(v) T ω ( v ) κ ≤ ε \kappa\le\varepsilon κ ≤ ε § 5 bounds the fixed-point residual and the distance to Z ⋆ Z^{\star} Z ⋆
Convergence guarantees of r2HPDHG
The reformulations below are measured against PDHG run on the standard form. The rates to improve on are the ones that method already has. Both are from [38 ] , which states them for the standard form. Appendix § 9 records the proof of Proposition 1(b) for T ω T_\omega T ω
Definition 2 (Sharpness modulus). α \alpha α α dist P ( v , Z ⋆ ) ≤ ∥ v − T ω ( v ) ∥ P \alpha\,\operatorname{dist}_{P}(v,\mathcal{Z}^{\star}) \le\lVert v-T_\omega(v)\rVert_{P} α dist P ( v , Z ⋆ ) ≤ ∥ v − T ω ( v ) ∥ P v v v
Proposition 1 ([38 ] , Corollary 1, and Theorem 2 with the reflection of Section 6).
(a) Run (13) from an anchor v a v^{\mathrm{a}} v a k ≥ 0 k \ge0 k ≥ 0
∥ v k − T ω ( v k ) ∥ P ≤ dist P ( v a , Z ⋆ ) k + 1 . \lVert v^{k}-T_\omega(v^{k})\rVert_{P} \le\frac{\operatorname{dist}_{P}(v^{\mathrm{a}},\mathcal{Z}^{\star})}{k+1}. ∥ v k − T ω ( v k ) ∥ P ≤ k + 1 dist P ( v a , Z ⋆ ) . (b) Let α > 0 \alpha> 0 α > 0 ω \omega ω j j j
∥ v j , k − T ω ( v j , k ) ∥ P ≤ e − 1 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P . (16) \lVert v^{j,k}-T_\omega(v^{j,k})\rVert_{P} \le e^{-1}\lVert v^{j,0}-T_\omega(v^{j,0})\rVert_{P}. \tag*{(16)} ∥ v j , k − T ω ( v j , k ) ∥ P ≤ e − 1 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P . (16) Writing τ j \tau^{j} τ j j j j τ j ≤ ⌈ e / α ⌉ \tau^{j} \le\lceil e/\alpha\rceil τ j ≤ ⌈ e / α ⌉ j ≥ 0 j \ge0 j ≥ 0
dist P ( v j + 1 , 0 , Z ⋆ ) ≤ e − j α ( τ 0 + 1 ) dist P ( v 0 , 0 , Z ⋆ ) . \operatorname{dist}_{P}(v^{j+1,0},\mathcal{Z}^{\star}) \le\frac{e^{-j}}{\alpha(\tau^{0}+1)}\operatorname{dist}_{P}(v^{0,0},\mathcal{Z}^{\star}). dist P ( v j + 1 , 0 , Z ⋆ ) ≤ α ( τ 0 + 1 ) e − j dist P ( v 0 , 0 , Z ⋆ ) . The rule (16) runs without α \alpha α β s u f f = e − 1 \beta_{\mathrm{suff}}=e^{-1} β suff = e − 1 ω \omega ω
Only α \alpha α
A positive α \alpha α X \mathcal{X} X Y \mathcal{Y} Y [3 , 36 , 17 ] is then finite on Z \mathcal{Z} Z v ∈ R n + m v \in\mathbb{R}^{n+m} v ∈ R n + m Z \mathcal{Z} Z α \alpha α
For the standard form the available constant is a Hoffman constant of the KKT system, H K K T H_{\mathrm{KKT}} H KKT [27 ] . Take the baseline at Λ = I \Lambda=I Λ = I η ≤ 1 / ( 2 ∥ A ∥ Λ ) \eta\le1/(2\lVert A\rVert_{\Lambda}) η ≤ 1/ ( 2 ∥ A ∥ Λ ) [38 ] analyzes. Write R : = 2 dist Λ ( v 0 , 0 , Z ⋆ ) R:=2\operatorname{dist}_{\Lambda}(v^{0,0},\mathcal{Z}^{\star}) R := 2 dist Λ ( v 0 , 0 , Z ⋆ ) z ⋆ z^{\star} z ⋆ v 0 , 0 v^{0,0} v 0 , 0 ρ : = R + 2 ∥ z ⋆ ∥ 2 \rho:=R+2\lVert z^{\star}\rVert_{2} ρ := R + 2 ∥ z ⋆ ∥ 2 [38 Lemmas 1 and 2(i), with the reflection of Section 6] , and on that ball the error bound of [3 Lemma 5] holds with constant H K K T 1 + 4 ρ 2 H_{\mathrm{KKT}}\sqrt{1+4\rho^{2}} H KKT 1 + 4 ρ 2 [38 Theorem 2, with the reflection of Section 6] with that constant put the standard form at
O ( ∥ A ∥ Λ R ε ) (17) O\left(\frac{\lVert A\rVert_{\Lambda}R}{\varepsilon}\right) \tag*{(17)} O ( ε ∥ A ∥ Λ R ) (17) applications of T ω T_\omega T ω ε \varepsilon ε
O ( ∥ A ∥ Λ H K K T 1 + 4 ρ 2 log 1 ε ) (18) O\left(\lVert A\rVert_{\Lambda}H_{\mathrm{KKT}}\sqrt{1+4\rho^{2}}\log\frac{1}{\varepsilon}\right) \tag*{(18)} O ( ∥ A ∥ Λ H KKT 1 + 4 ρ 2 log ε 1 ) (18) applications to contract dist Λ ( v , Z ⋆ ) \operatorname{dist}_{\Lambda}(v,\mathcal{Z}^{\star}) dist Λ ( v , Z ⋆ ) ε \varepsilon ε ∥ A ∥ Λ \lVert A\rVert_{\Lambda} ∥ A ∥ Λ X \mathcal{X} X Y \mathcal{Y} Y
Why reformulate
Absorption changes a solver’s performance through two independent channels: the cost of an iteration, and the number of iterations required. In this section we explain how reformulation affects these channels.
What changes numerically
Cost per iteration. An iteration of (2) in the standard form requires two sparse matrix–vector multiplications, one with A A A A T A^{\mathsf{T}} A T A T y A^{\mathsf{T}}y A T y A x Ax A x
The difference is important on GPUs. An SpMV kernel is inherently irregular, whereas an elementwise kernel is uniform. SpMV gathers entries of x x x A x Ax A x A T y A^{\mathsf{T}}y A T y
This speedup applies specifically to the interaction term. The tradeoff is that the additional projection cost introduced by the split can become dominant, as it does for the JRP example in Table 11 .
Table 11. Median microseconds per iteration. PRIA’s iterations are slower because of the simplex projections, and in particular the Δ 2 K 0 \Delta_{2}^{K_{0}} Δ 2 K 0 T T T N T NT N T
Iteration count. Reducing the number of dual variables lowers the cost of every vector operation, but the more consequential benefit is a reduction in the number of iterations. The constants of § 4.2 bound the saving; the measured counts realize it. Against cuPDLPx on the original formulation, PRAQ converges in 1.7–4.3 times fewer iterations and PRIA in 13.6–52.3 times fewer, over the sizes where cuPDLPx converges at all.
What changes theoretically
A formulation enters the convergence theory of § 3.2 through a short list of constants. Absorption acts on each separately, so its effect can be attributed term by term.
A data-free step size. The admissible step size satisfies η < 1 / ∥ A ∥ Λ \eta< 1 / \lVert A \rVert_{\Lambda} η < 1/ ∥ A ∥ Λ ∥ A ∥ Λ \lVert A \rVert_{\Lambda} ∥ A ∥ Λ § 3.1 ). In the standard formulation, the step size is therefore determined by the norm of the numerical constraint matrix. Diagonal rescaling can reduce this coupling norm to O ( 1 ) O(1) O ( 1 ) [47 ] , but does so by absorbing the data dependence into the scaling. Under an absorption split, by contrast, the coupling norm can be bounded purely from the structure of the split, yielding a data-free step size that is fixed for the entire problem family. In particular, A R A_{\mathcal{R}} A R ∥ A R ∥ 2 ≤ ∥ A ∥ 2 \lVert A_{\mathcal{R}} \rVert_{2} \le\lVert A \rVert_{2} ∥ A R ∥ 2 ≤ ∥ A ∥ 2 ∥ A ∥ Λ ≤ 1 \lVert A \rVert_{\Lambda} \le1 ∥ A ∥ Λ ≤ 1 ∥ A ∥ Λ = 2 \lVert A \rVert_{\Lambda} = \sqrt{2} ∥ A ∥ Λ = 2
A tighter termination certificate. The absorbed constraints and box constraints are satisfied exactly by construction, so the KKT residual needs to measure only the constraints in R \mathcal{R} R
An explicit distance scale. The sublinear rate of § 3.2 is measured against the distance scale R R R V \mathcal{V} V Z \mathcal{Z} Z R R R Z \mathcal{Z} Z R ≤ 2 diam Λ ( Z ) R \le2\operatorname{diam}_{\Lambda}(\mathcal{Z}) R ≤ 2 diam Λ ( Z ) diam Λ ( Z ) 2 ≤ 2 m + 2 \operatorname{diam}_{\Lambda}(\mathcal{Z})^2 \le2m+2 diam Λ ( Z ) 2 ≤ 2 m + 2 [33 ] and diam Λ ( Z ) 2 ≤ 2 N T + 2 T \operatorname{diam}_{\Lambda}(\mathcal{Z})^2 \le2NT+2T diam Λ ( Z ) 2 ≤ 2 N T + 2 T [40 ] .
Sharpness constant. Restarted PDHG converges linearly at a rate set by a sharpness constant, and every linear program supplies one in every formulation (§ 3.2 ). The standard-form analysis needs the Hoffman constant of the whole KKT system. Absorption removes the multipliers of E \mathcal{E} E explicit is a separate question, and for QUBO we bound it in Proposition 4.
Table 2 summarizes the four constants.
Table 2. The constants the convergence theory reads off a formulation, before and after a split.
Applications
In this section we apply the reformulation to two case studies. Both splits are complete, so X \mathcal{X} X Y \mathcal{Y} Y § 2.4 and instantiated against § 3 , and they differ in whether the sharpness modulus of Definition 2 can be made explicit. Tables 3 and 8 record the theoretical convergence rate improvements against (17) and (18) . Tables 4 and 10 give computational results.
Quadratic unconstrained binary optimization (QUBO)
Given a symmetric W ∈ R n × n W \in\mathbb{R}^{n \times n} W ∈ R n × n
min x ∈ { 0 , 1 } n x T W x . (19) \min_{x \in\{0,1\}^{n}} x^{\mathsf{T}} W x. \tag*{(19)} x ∈ { 0 , 1 } n min x T W x . (19) It is NP-hard and is the standard target format for Ising machines and quantum annealers. Its LP relaxation, the roof dual [23 ] , is the workhorse bound inside branch-and-bound and the source of variable fixings. This section reformulates that relaxation and instantiates Section 3 on it, with explicit constants.
The standard linearization
Using x i 2 = x i x_i^2=x_i x i 2 = x i
x T W x = ∑ i d i x i + ∑ i < j 2 W i j x i x j , d i : = W i i , (20) x^{\mathsf{T}}Wx=\sum_{i}d_ix_i+\sum_{i<j}2W_{ij}x_ix_j,\qquad d_i:=W_{ii}, \tag*{(20)} x T W x = i ∑ d i x i + i < j ∑ 2 W ij x i x j , d i := W ii , (20) and set E + : = { i < j : W i j > 0 } E_+:=\{i<j:W_{ij}>0\} E + := { i < j : W ij > 0 } E − : = { i < j : W i j < 0 } E_-:=\{i<j:W_{ij}<0\} E − := { i < j : W ij < 0 } E : = E + ∪ E − E:=E_+\cup E_- E := E + ∪ E − m : = ∣ E ∣ m:=|E| m := ∣ E ∣ w e : = 2 ∣ W i j ∣ w_e:=2|W_{ij}| w e := 2∣ W ij ∣ E E E [ n ] [n] [ n ] W W W (20) . Write deg ( i ) \operatorname{deg}(i) deg ( i ) i i i z e z_e z e [19 ] , whose relaxation we write (P):
m i n i m i z e x , z d T x + ∑ e ∈ E + w e z e − ∑ e ∈ E − w e z e s u b j e c t t o z e ≥ x i + x j − 1 , z e ≥ 0 , e = { i , j } ∈ E + , z e ≤ x i , z e ≤ x j , e = { i , j } ∈ E − , x ∈ [ 0 , 1 ] n , z ∈ [ 0 , 1 ] E . (21) \begin{aligned}
\mathop{\mathrm{minimize}}_{x,z}\quad&d^{\mathsf{T}}x+\sum_{e\in E_+}w_ez_e-\sum_{e\in E_-}w_ez_e\\
\mathrm{subject\ to}\quad&z_e\ge x_i+x_j-1,\quad z_e\ge0,\quad e=\{i,j\}\in E_+,\\
&z_e\le x_i,\quad z_e\le x_j,\quad e=\{i,j\}\in E_-,\\
&x\in[0,1]^n,\quad z\in[0,1]^E. \tag*{(21)}
\end{aligned} minimize x , z subject to d T x + e ∈ E + ∑ w e z e − e ∈ E − ∑ w e z e z e ≥ x i + x j − 1 , z e ≥ 0 , e = { i , j } ∈ E + , z e ≤ x i , z e ≤ x j , e = { i , j } ∈ E − , x ∈ [ 0 , 1 ] n , z ∈ [ 0 , 1 ] E . (21) Observations on the linearization. The following observations clarify the role of the variable bounds and diagonal terms in (21) .
Lower bounds on z e z_e z e For e ∈ E + e\in E_+ e ∈ E + w e w_e w e z e = max { 0 , x i + x j − 1 } z_e=\max\{0,x_i+x_j-1\} z e = max { 0 , x i + x j − 1 } z e ≥ 0 z_e\ge0 z e ≥ 0 x i = x j = 0 x_i=x_j=0 x i = x j = 0 z e = − 1 z_e=-1 z e = − 1 − w e -w_e − w e e ∈ E − e\in E_- e ∈ E − z e = min { x i , x j } ∈ [ 0 , 1 ] z_e=\min\{x_i,x_j\}\in[0,1] z e = min { x i , x j } ∈ [ 0 , 1 ]
Upper bounds on x x x The bounds x ≤ 1 x\le1 x ≤ 1 x ∈ { 0 , 1 } n x\in\{0,1\}^n x ∈ { 0 , 1 } n e = { i , j } ∈ E + e=\{i,j\}\in E_+ e = { i , j } ∈ E + z e ≥ x i + x j − 1 z_e\ge x_i+x_j-1 z e ≥ x i + x j − 1 z e ≤ 1 z_e\le1 z e ≤ 1 x i + x j ≤ 2 x_i+x_j\le2 x i + x j ≤ 2
Treatment of diagonal terms. The binary identity x i 2 = x i x_i^2=x_i x i 2 = x i W i i x i W_{ii}x_i W ii x i (20) . Introducing a separate product variable z i i z_{ii} z ii W i i > 0 W_{ii}>0 W ii > 0 W i i x i W_{ii}x_i W ii x i W i i max { 0 , 2 x i − 1 } W_{ii}\max\{0,2x_i-1\} W ii max { 0 , 2 x i − 1 } x i = 1 2 x_i=\frac{1}{2} x i = 2 1 W i i / 2 W_{ii}/2 W ii /2
The reformulation
The reformulation. QUBO has the max-selection structure of Section 2.4 on both edge families, with groups of size two. On E + E_+ E + + w e ≥ 0 +w_e\ge0 + w e ≥ 0 z e = max { 0 , x i + x j − 1 } z_e=\max\{0,x_i+x_j-1\} z e = max { 0 , x i + x j − 1 } ( 0 , x i + x j − 1 ) (0,x_i+x_j-1) ( 0 , x i + x j − 1 ) E − E_- E − − w e ≤ 0 -w_e\le0 − w e ≤ 0 z e = min { x i , x j } z_e=\min\{x_i,x_j\} z e = min { x i , x j } − w e min { x i , x j } = w e max { − x i , − x j } -w_e\min\{x_i,x_j\}=w_e\max\{-x_i,-x_j\} − w e min { x i , x j } = w e max { − x i , − x j } E − E_- E − ( − x i , − x j ) (-x_i,-x_j) ( − x i , − x j ) w e w_e w e z z z ( P ) = min x ∈ [ 0 , 1 ] n F ( x ) (P)=\min_{x\in[0,1]^n}F(x) ( P ) = min x ∈ [ 0 , 1 ] n F ( x )
F ( x ) : = d T x + ∑ e = { i , j } ∈ E + w e max { 0 , x i + x j − 1 } + ∑ e = { i , j } ∈ E − w e max { − x i , − x j } . (22) F(x):=d^{\mathsf{T}}x+\sum_{e=\{i,j\}\in E_+}w_e\max\{0,x_i+x_j-1\}+\sum_{e=\{i,j\}\in E_-}w_e\max\{-x_i,-x_j\}. \tag*{(22)} F ( x ) := d T x + e = { i , j } ∈ E + ∑ w e max { 0 , x i + x j − 1 } + e = { i , j } ∈ E − ∑ w e max { − x i , − x j } . (22) Write c k : = d k − ∑ e ∈ E − : k = max e w e c_k:=d_k-\sum_{e\in E_-:k=\max e}w_e c k := d k − ∑ e ∈ E − : k = m a x e w e q e : = − 1 q_e:=-1 q e := − 1 E + E_+ E + E − E_- E −
Y : = ∏ e ∈ E [ 0 , w e ] , (23) \mathcal{Y}:=\prod_{e\in E}[0,w_e], \tag*{(23)} Y := e ∈ E ∏ [ 0 , w e ] , (23) and let A ~ : R n → R E \widetilde{A}:\mathbb{R}^n\to\mathbb{R}^E A : R n → R E
( A ~ x ) e : = { x i − x j , e = { i , j } ∈ E − , − ( x i + x j ) , e = { i , j } ∈ E + , i < j . (24) (\widetilde{A}x)_e:=
\begin{cases}
x_i-x_j,&e=\{i,j\}\in E_-,\\
-(x_i+x_j),&e=\{i,j\}\in E_+,
\end{cases}
\qquad i<j. \tag*{(24)} ( A x ) e := { x i − x j , − ( x i + x j ) , e = { i , j } ∈ E − , e = { i , j } ∈ E + , i < j . (24) Proposition 2 (The QUBO reformulation). (P) is the saddle problem
min x ∈ [ 0 , 1 ] n max λ ∈ Y L ( x , λ ) : = ⟨ c , x ⟩ + ⟨ λ , q − A ~ x ⟩ (25) \min_{x\in[0,1]^n}\ \max_{\lambda\in\mathcal{Y}} L(x,\lambda) := \langle c,x\rangle+\langle\lambda,q-\widetilde{A}x\rangle \tag*{(25)} x ∈ [ 0 , 1 ] n min λ ∈ Y max L ( x , λ ) := ⟨ c , x ⟩ + ⟨ λ , q − A x ⟩ (25) which is (1) with X = [ 0 , 1 ] n \mathcal{X}=[0,1]^n X = [ 0 , 1 ] n F ( x ) F(x) F ( x ) x ∈ [ 0 , 1 ] n x\in[0,1]^n x ∈ [ 0 , 1 ] n x ∗ x^* x ∗
Proof. We prove Proposition 2 in three steps: price the two rows of each edge, absorb the simplex that pricing produces, and evaluate the inner maximum.
Step 1: the dual of each edge is a simplex. Write v e v_e v e e e e z e z_e z e E + E_+ E + − z e -z_e − z e E − E_- E − w e v e w_ev_e w e v e M e , 1 , M e , 2 ≥ 0 M_{e,1},M_{e,2}\ge0 M e , 1 , M e , 2 ≥ 0
e ∈ E − : v e ( w e − M e , 1 − M e , 2 ) − M e , 1 x i − M e , 2 x j , e ∈ E + : v e ( w e − M e , 1 − M e , 2 ) + M e , 2 ( x i + x j − 1 ) , e\in E_-:\quad v_e(w_e-M_{e,1}-M_{e,2})-M_{e,1}x_i-M_{e,2}x_j,\qquad e\in E_+:\quad v_e(w_e-M_{e,1}-M_{e,2})+M_{e,2}(x_i+x_j-1), e ∈ E − : v e ( w e − M e , 1 − M e , 2 ) − M e , 1 x i − M e , 2 x j , e ∈ E + : v e ( w e − M e , 1 − M e , 2 ) + M e , 2 ( x i + x j − 1 ) , and each v e v_e v e
M e , 1 + M e , 2 = w e , M e ≥ 0 (26) M_{e,1}+M_{e,2}=w_e,\qquad M_e\ge0 \tag*{(26)} M e , 1 + M e , 2 = w e , M e ≥ 0 (26) one per edge, over disjoint groups. This is the assignment family of Table 1 , with a simplex of total w e w_e w e
Step 2: absorb the simplices. Definition 1 applied to those rows puts each M e M_e M e Δ w e \Delta^{w_e} Δ w e E + E_+ E + x x x M e , 1 M_{e,1} M e , 1 λ e : = M e , 2 \lambda_e:=M_{e,2} λ e := M e , 2
min x ∈ [ 0 , 1 ] n max M e ∈ Δ w e ⟨ d , x ⟩ − ∑ e = { i , j } ∈ E − ( M e , 1 x i + M e , 2 x j ) + ∑ e = { i , j } ∈ E + λ e ( x i + x j − 1 ) . (27) \min_{x\in[0,1]^n}\ \max_{M_e\in\Delta^{w_e}}\ \langle d,x\rangle-\sum_{e=\{i,j\}\in E_-}(M_{e,1}x_i+M_{e,2}x_j)+\sum_{e=\{i,j\}\in E_+}\lambda_e(x_i+x_j-1). \tag*{(27)} x ∈ [ 0 , 1 ] n min M e ∈ Δ w e max ⟨ d , x ⟩ − e = { i , j } ∈ E − ∑ ( M e , 1 x i + M e , 2 x j ) + e = { i , j } ∈ E + ∑ λ e ( x i + x j − 1 ) . (27) Each Δ w e \Delta^{w_e} Δ w e
e ∈ E − : M e = ( λ e , w e − λ e ) , e ∈ E + : M e = ( λ e , λ e ) , λ e ∈ [ 0 , w e ] , (28) e\in E_-:\quad M_e=(\lambda_e,w_e-\lambda_e),\qquad e\in E_+:\quad M_e=(\lambda_e,\lambda_e),\qquad\lambda_e\in[0,w_e], \tag*{(28)} e ∈ E − : M e = ( λ e , w e − λ e ) , e ∈ E + : M e = ( λ e , λ e ) , λ e ∈ [ 0 , w e ] , (28) and the dual set becomes the box (23). Substituting (28) into (27) leaves the same saddle problem with one dual coordinate per edge:
min x ∈ [ 0 , 1 ] n max λ ∈ Y ⟨ d , x ⟩ − ∑ e = { i , j } ∈ E − ( λ e x i + ( w e − λ e ) x j ) + ∑ e = { i , j } ∈ E + λ e ( x i + x j − 1 ) . (29) \min_{x\in[0,1]^n}\ \max_{\lambda\in\mathcal{Y}}\ \langle d,x\rangle-\sum_{e=\{i,j\}\in E_-}\bigl(\lambda_ex_i+(w_e-\lambda_e)x_j\bigr)+\sum_{e=\{i,j\}\in E_+}\lambda_e(x_i+x_j-1). \tag*{(29)} x ∈ [ 0 , 1 ] n min λ ∈ Y max ⟨ d , x ⟩ − e = { i , j } ∈ E − ∑ ( λ e x i + ( w e − λ e ) x j ) + e = { i , j } ∈ E + ∑ λ e ( x i + x j − 1 ) . (29) On E − E_- E − − w e x j -w_ex_j − w e x j c c c q q q A ~ \widetilde{A} A
Step 3: the inner maximum is F F F Substitute (28) into (27); everything but ⟨ d , x ⟩ \langle d,x\rangle ⟨ d , x ⟩ e = { i , j } ∈ E − e=\{i,j\}\in E_- e = { i , j } ∈ E − − λ e x i − ( w e − λ e ) x j = − w e x j + λ e ( x j − x i ) -\lambda_ex_i-(w_e-\lambda_e)x_j=-w_ex_j+\lambda_e(x_j-x_i) − λ e x i − ( w e − λ e ) x j = − w e x j + λ e ( x j − x i ) c c c λ e ( q e − ( A ~ x ) e ) \lambda_e(q_e-(\widetilde{A}x)_e) λ e ( q e − ( A x ) e ) [ 0 , w e ] [0,w_e] [ 0 , w e ] − w e x j + w e max { 0 , x j − x i } -w_ex_j+w_e\max\{0,x_j-x_i\} − w e x j + w e max { 0 , x j − x i } − w e x i -w_ex_i − w e x i x i ≤ x j x_i\le x_j x i ≤ x j − w e x j -w_ex_j − w e x j w e max { − x i , − x j } w_e\max\{-x_i,-x_j\} w e max { − x i , − x j } e ∈ E + e\in E_+ e ∈ E + λ e ( x i + x j − 1 ) = λ e ( q e − ( A ~ x ) e ) \lambda_e(x_i+x_j-1)=\lambda_e(q_e-(\widetilde{A}x)_e) λ e ( x i + x j − 1 ) = λ e ( q e − ( A x ) e ) [ 0 , w e ] [0,w_e] [ 0 , w e ] w e max { 0 , x i + x j − 1 } w_e\max\{0,x_i+x_j-1\} w e max { 0 , x i + x j − 1 } max λ ∈ Y L ( x , λ ) = F ( x ) \max_{\lambda\in\mathcal{Y}}L(x,\lambda)=F(x) max λ ∈ Y L ( x , λ ) = F ( x ) [ 0 , 1 ] n [0,1]^n [ 0 , 1 ] n □ \square □
Diagonal scaling. Λ X \Lambda_X Λ X i i i deg ( i ) \operatorname{deg}(i) deg ( i ) Λ Y \Lambda_Y Λ Y [47 ] at exponent 1, since ∑ e ∣ A ~ e i ∣ = deg ( i ) \sum_e|\widetilde{A}_{ei}|=\operatorname{deg}(i) ∑ e ∣ A e i ∣ = deg ( i ) ∑ i ∣ A ~ e i ∣ = 2 \sum_i|\widetilde{A}_{ei}|=2 ∑ i ∣ A e i ∣ = 2 i i i 1 / deg ( i ) 1/\operatorname{deg}(i) 1/ deg ( i ) 1 / 2 1/2 1/2
Step size. Under these scalings the coupling operator satisfies
∥ A ~ ∥ Λ = ∥ Λ y − 1 / 2 A ~ Λ x − 1 / 2 ∥ 2 ≤ 1. (30) \left\lVert\widetilde{A} \right\rVert_{\Lambda} = \left\lVert\Lambda_{y}^{-1/2}\widetilde{A}\Lambda_{x}^{-1/2} \right\rVert_{2} \le1. \tag*{(30)} A Λ = Λ y − 1/2 A Λ x − 1/2 2 ≤ 1. (30) This is the Pock–Chambolle bound at exponent 1 [47 ] . Any 0 < η < 1 0 < \eta< 1 0 < η < 1 η = 0.999 \eta= 0.999 η = 0.999
Per-iteration cost. With the scalings specified above, the PDHG updates (2) take the form
x k + 1 : = proj [ 0 , 1 ] n ( x k − η ω Λ x − 1 ( c − A ~ T λ k ) ) , x^{k+1} := \operatorname{proj}_{[0,1]^{n}}\left(x^{k} - \frac{\eta}{\omega}\Lambda_{x}^{-1}\left(c - \widetilde{A}^{\mathsf{T}}\lambda^{k}\right)\right), x k + 1 := proj [ 0 , 1 ] n ( x k − ω η Λ x − 1 ( c − A T λ k ) ) , λ k + 1 : = proj Y ( λ k + η ω 2 ( q − A ~ ( 2 x k + 1 − x k ) ) ) . (31) \lambda^{k+1} := \operatorname{proj}_{\mathcal{Y}}\left(\lambda^{k} + \frac{\eta\omega}{2}\left(q - \widetilde{A}\left(2x^{k+1} - x^{k}\right)\right)\right). \tag*{(31)} λ k + 1 := proj Y ( λ k + 2 η ω ( q − A ( 2 x k + 1 − x k ) ) ) . (31) Each iteration requires two box projections and one application each of A ~ \widetilde{A} A A ~ T \widetilde{A}^{\mathsf{T}} A T O ( n + m ) O(n+m) O ( n + m ) A ~ \widetilde{A} A ± 1 \pm1 ± 1 W i j W_{ij} W ij
Termination. The QUBO reformulation allows termination to be checked using the primal–dual gap alone. At each scheduled check, we evaluate the projected PDHG point ( x , λ ) : = T ω ( v ) (x,\lambda) := T_{\omega}(v) ( x , λ ) := T ω ( v ) x ∈ [ 0 , 1 ] n x \in[0,1]^{n} x ∈ [ 0 , 1 ] n λ ∈ Y \lambda\in\mathcal{Y} λ ∈ Y L P L_P L P L D L_D L D g ( λ ) : = c − A ~ T λ g(\lambda) := c - \widetilde{A}^{\mathsf{T}}\lambda g ( λ ) := c − A T λ
L D ( λ ) = ∑ i = 1 n min { 0 , g i ( λ ) } − ∑ e ∈ E + λ e , L P ( x ) = F ( x ) . (32) L_D(\lambda) = \sum_{i=1}^{n}\min\{0,g_i(\lambda)\} - \sum_{e\in E_{+}}\lambda_e,\qquad L_P(x) = F(x). \tag*{(32)} L D ( λ ) = i = 1 ∑ n min { 0 , g i ( λ )} − e ∈ E + ∑ λ e , L P ( x ) = F ( x ) . (32) The lower bound follows by minimizing L ( x , λ ) = ⟨ g ( λ ) , x ⟩ + ⟨ λ , q ⟩ L(x,\lambda) = \langle g(\lambda),x\rangle+ \langle\lambda,q\rangle L ( x , λ ) = ⟨ g ( λ ) , x ⟩ + ⟨ λ , q ⟩
We terminate when the relative duality gap satisfies
L P ( x ) − L D ( λ ) 1 + ∣ L P ( x ) ∣ + ∣ L D ( λ ) ∣ ≤ ε , \frac{L_P(x)-L_D(\lambda)}{1+\lvert L_P(x)\rvert+\lvert L_D(\lambda)\rvert}\le\varepsilon, 1 + ∣ L P ( x )∣ + ∣ L D ( λ )∣ L P ( x ) − L D ( λ ) ≤ ε , where the objective values in this test are evaluated at the original data scale. The monitored point satisfies both box constraints and the eliminated variables are recovered exactly, so the gap alone certifies.
Theoretical improvements
Assumption 1 (Integral data). d ∈ Z n d \in\mathbb{Z}^{n} d ∈ Z n w ∈ Z E w \in\mathbb{Z}^{E} w ∈ Z E
Theorem 1 (Convergence of PRAQ). For the roof dual (21) of an instance W W W n ≥ 5 n \ge5 n ≥ 5 ω = 1 \omega= 1 ω = 1 0 < η < 1 0 < \eta< 1 0 < η < 1 ( d ; w ) / ∥ ( d ; w ) ∥ 2 (d;w)/\lVert(d;w)\rVert_{2} ( d ; w ) / ∥( d ; w ) ∥ 2 v k v^{k} v k v j , k v^{j,k} v j , k j j j ∥ ⋅ ∥ P \lVert\cdot\rVert_{P} ∥ ⋅ ∥ P
(a) Sublinear convergence without restarts. For every $k \ge 0,
∥ v k − T ω ( v k ) ∥ P ≤ 2 k + 1 2 m η . \left\lVert v^{k}-T_{\omega}(v^{k})\right\rVert_{P}\le\frac{2}{k+1}\sqrt{\frac{2m}{\eta}}. v k − T ω ( v k ) P ≤ k + 1 2 η 2 m . (b) Linear convergence with adaptive restarts. Restarting by sufficient decay (16), every epoch has length at most
⌈ e ( 1 + 16 η ∥ W ∥ F m n ) ⌉ \left\lceil e\left(1+\frac{16}{\eta}\lVert W\rVert_F m\sqrt{n}\right)\right\rceil ⌈ e ( 1 + η 16 ∥ W ∥ F m n ) ⌉ and, for every j ≥ 0 j\ge0 j ≥ 0
dist P ( v j + 1 , 0 , Z ⋆ ) ≤ e − j ( 1 + 16 η ∥ W ∥ F m n ) 2 m η . \operatorname{dist}_P(v^{j+1,0},\mathcal{Z}^{\star})\le e^{-j}\left(1+\frac{16}{\eta}\lVert W\rVert_F m\sqrt{n}\right)\sqrt{\frac{2m}{\eta}}. dist P ( v j + 1 , 0 , Z ⋆ ) ≤ e − j ( 1 + η 16 ∥ W ∥ F m n ) η 2 m . Proof. The reformulation has
∥ A ~ ∥ Λ ≤ 1 , diam Λ ( Z ) 2 = 2 m + 2 ∥ w ∥ 2 2 ∥ ( d ; w ) ∥ 2 2 ≤ 2 m + 2 , (33) \lVert\widetilde{A}\rVert_{\Lambda}\le1,\qquad\operatorname{diam}_{\Lambda}(\mathcal{Z})^2=2m+\frac{2\lVert w\rVert_2^2}{\lVert(d;w)\rVert_2^2}\le2m+2, \tag*{(33)} ∥ A ∥ Λ ≤ 1 , diam Λ ( Z ) 2 = 2 m + ∥( d ; w ) ∥ 2 2 2 ∥ w ∥ 2 2 ≤ 2 m + 2 , (33) by (30), and by Λ X \Lambda_X Λ X i i i deg ( i ) \deg(i) deg ( i ) Λ Y \Lambda_Y Λ Y ∑ i deg ( i ) = 2 m \sum_i\deg(i)=2m ∑ i deg ( i ) = 2 m ∥ w ∥ 2 ≤ ∥ ( d ; w ) ∥ 2 \lVert w\rVert_2\le\lVert(d;w)\rVert_2 ∥ w ∥ 2 ≤ ∥( d ; w ) ∥ 2
1 α ≤ 1 + 16 η ∥ W ∥ F m n (34) \frac{1}{\alpha}\le1+\frac{16}{\eta}\lVert W\rVert_Fm\sqrt{n} \tag*{(34)} α 1 ≤ 1 + η 16 ∥ W ∥ F m n (34) by Proposition 4, proved in Appendix § 8 . The initial point lies in Z \mathcal{Z} Z ω = 1 \omega=1 ω = 1 m ≥ 1 m\ge1 m ≥ 1
dist P ( v 0 , 0 , Z ⋆ ) ≤ 2 η diam Λ ( Z ) ≤ 2 2 m η . \operatorname{dist}_P(v^{0,0},\mathcal{Z}^{\star})\le\sqrt{\frac{2}{\eta}}\operatorname{diam}_{\Lambda}(\mathcal{Z})\le2\sqrt{\frac{2m}{\eta}}. dist P ( v 0 , 0 , Z ⋆ ) ≤ η 2 diam Λ ( Z ) ≤ 2 η 2 m . Part (a) is Proposition 1(a) with this bound, part (b) Proposition 1(b) with this α \alpha α τ 0 ≥ 1 \tau^0\ge1 τ 0 ≥ 1
Table 3 presents the effect of reformulation on the constants governing convergence. For the standard PDLP formulation, the sublinear bound involves the full coupling matrix through ∥ A ∥ Λ \lVert A\rVert_{\Lambda} ∥ A ∥ Λ R R R ∥ A ~ ∥ Λ ≤ 1 \lVert\widetilde{A}\rVert_{\Lambda}\le1 ∥ A ∥ Λ ≤ 1 m m m n n n m m m H K K T H_{\mathrm{KKT}} H KKT ρ \rho ρ m m m n n n ∥ W ∥ F \lVert W\rVert_F ∥ W ∥ F
Table 3. Applications of T ω T_\omega T ω ε \varepsilon ε n n n m m m dist Λ ( v , Z ⋆ ) \operatorname{dist}_{\Lambda}(v,\mathcal{Z}^{\star}) dist Λ ( v , Z ⋆ ) ∥ ⋅ ∥ P \lVert\cdot\rVert_P ∥ ⋅ ∥ P dist P ( v , Z ⋆ ) \operatorname{dist}_P(v,\mathcal{Z}^{\star}) dist P ( v , Z ⋆ )
Reformulation also changes the operations within each iteration. Although both formulations have O ( n + m ) O(n+m) O ( n + m )
Computational results
Implementation. The product A ~ x \tilde{A}x A ~ x A ~ ⊤ λ \tilde{A}^{\top}\lambda A ~ ⊤ λ [ 0 , 1 ] [0,1] [ 0 , 1 ] [ 0 , w e ] [0,w_e] [ 0 , w e ]
Instances. Three generated families, 2282 instances with varying edge counts. Given n n n 10 n 10n 10 n 1 4 \frac{1}{4} 4 1 − 1 , 1 -1,1 − 1 , 1
Protocol. Every solver that takes a tolerance is asked for a relative duality gap of 10 − 8 10^{-8} 1 0 − 8 10 − 9 10^{-9} 1 0 − 9 [35 ] run on the same NVIDIA H100 80 GB, one instance at a time. The third baseline computes the roof dual as a maximum flow with OR-Tools’ SimpleMaxFlow [46 ] on CPU. Those three each solved all 2282 instances. Gurobi [21 ] was run on the same instances on 16 CPU cores under the same limit, solving 1954 instances.
Results. Table 4 is the running time comparison. PRAQ gains on cuPDLPx and Gurobi in solve time at every instance size, and is faster than the maximum flow solver on all but the smallest instances.
Table 4. Solve time in seconds, shifted (by 1 ms) geometric mean by edge count. The paired column gives how many instances Gurobi finished within the limit; an unfinished bucket’s entry is a floor. PRAQ, cuPDLPx and max-flow finished all instances.
Tables 5 and 6 show that the running-time advantage of PRAQ comes from both cheaper iterations and fewer iterations. Across the five size groups, PRAQ reduces the median time per iteration by factors of 2.3–4.8 and the iteration count by factors of 1.7–4.3 relative to cuPDLPx. The lower iteration cost is consistent with evaluating the coupling through custom gather and scatter operations over the edge list, together with coordinatewise box projections. The reduction in iteration counts provides an additional benefit beyond the lower computational cost of each step, illustrating the two mechanisms discussed in § 4.1 .
Table 7 shows a corresponding reduction in peak device memory. For the two largest size groups, cuPDLPx uses 3.0–3.5 times as much device memory as PRAQ. In particular, for instances with 10 7 10^7 1 0 7 10 8 10^8 1 0 8
Table 7. Largest device memory held during a solve, sampled with NVML.
Table 5. Median microseconds per iteration.
Table 6. Iterations to converge, shifted geometric mean with shift 10.
The joint replenishment problem (JRP)
The joint replenishment problem (JRP) is a classical model in deterministic inventory theory, studied since the 1960s as the multi-item extension of the Wagner–Whitin lot-sizing problem [48 ] . A review of its literature and industrial applications is in [28 ] . Its stochastic variant remains an active computational target in high dimensions [6 ] . Over a horizon of T T T N N N i i i K i K_i K i some item is ordered additionally costs K 0 K_0 K 0 i i i t t t s ≤ t s \le t s ≤ t H i s t H_{ist} H i s t [4 ] and remains APX-hard in natural special cases [42 ] .
The strongest approximation algorithms for the JRP are LP-based. Levi, Roundy and Shmoys gave a primal–dual 2-approximation and showed that the natural LP relaxation has integrality gap at most 2 [29 , 30 ] . LP-rounding refinements improved the ratio to 1.8 [31 ] and then 1.791 [9 ] . Every one of these algorithms must first solve the relaxation, so the relaxation is the computational bottleneck at scale: it has Θ ( N T 2 ) \Theta(NT^2) Θ ( N T 2 )
The LP relaxation
The variable x i s t x_{ist} x i s t i i i t t t s ≤ t s \le t s ≤ t y i s y_{is} y i s i i i s s s z s z_s z s s s s
minimize x , y , z ∑ s K 0 z s + ∑ i , s K i y i s + ∑ i , t , s ≤ t H i s t x i s t subject to ∑ s ≤ t x i s t = 1 , i ∈ [ N ] , t ∈ [ T ] , x i s t − y i s ≤ 0 , i ∈ [ N ] , t ∈ [ T ] , s ≤ t , x i s t − z s ≤ 0 , i ∈ [ N ] , t ∈ [ T ] , s ≤ t , x i s t , y i s , z s ≥ 0. (35) \begin{aligned}
\underset{x,y,z}{\operatorname{minimize}}\quad& \sum_s K_0z_s + \sum_{i,s} K_i y_{is} + \sum_{i,t,s\le t} H_{ist}x_{ist}\\
\text{subject to}\quad& \sum_{s\le t} x_{ist}=1,\quad i\in[N],\ t\in[T],\\
& x_{ist}-y_{is}\le0,\quad i\in[N],\ t\in[T],\ s\le t,\\
& x_{ist}-z_s\le0,\quad i\in[N],\ t\in[T],\ s\le t,\\
& x_{ist},\ y_{is},\ z_s\ge0. \tag*{(35)}
\end{aligned} x , y , z minimize subject to s ∑ K 0 z s + i , s ∑ K i y i s + i , t , s ≤ t ∑ H i s t x i s t s ≤ t ∑ x i s t = 1 , i ∈ [ N ] , t ∈ [ T ] , x i s t − y i s ≤ 0 , i ∈ [ N ] , t ∈ [ T ] , s ≤ t , x i s t − z s ≤ 0 , i ∈ [ N ] , t ∈ [ T ] , s ≤ t , x i s t , y i s , z s ≥ 0. (35) The assignment rows force every demand to be served; the covering rows charge the item and joint order costs whenever an assignment uses period s s s K , K 0 ≥ 0 K,K_0\ge0 K , K 0 ≥ 0
y i s = max t ≥ s x i s t , z s = max i , t ≥ s x i s t , y_{is}=\max_{t\ge s}x_{ist},\qquad z_s=\max_{i,\ t\ge s}x_{ist}, y i s = t ≥ s max x i s t , z s = i , t ≥ s max x i s t , which is the max-selection structure of § 2.4 with two dual blocks, whose groups are indexed by ( i , s ) (i,s) ( i , s ) s s s ( vec ( H ) ; K ; K 0 ) / ∥ ( vec ( H ) ; K ; K 0 ) ∥ 2 (\operatorname{vec}(H);K;K_0)/\lVert(\operatorname{vec}(H);K;K_0)\rVert_2 ( vec ( H ) ; K ; K 0 ) / ∥( vec ( H ) ; K ; K 0 ) ∥ 2
The reformulation
Writing R u N × T × T \mathbb{R}_u^{N\times T\times T} R u N × T × T
Δ 1 , 3 1 : = { X ∈ R u N × T × T ∣ ∑ s = 1 t X i s t = 1 , i ∈ [ N ] , t ∈ [ T ] } , \Delta_{1,3}^{1}:=\left\{X\in\mathbb{R}_u^{N\times T\times T}\mathrel{\Big|}\sum_{s=1}^{t}X_{ist}=1,\ i\in[N],\ t\in[T]\right\}, Δ 1 , 3 1 := { X ∈ R u N × T × T s = 1 ∑ t X i s t = 1 , i ∈ [ N ] , t ∈ [ T ] } , Δ 1 , 2 K : = { Y ~ ∈ R u N × T × T ∣ ∑ t = s T Y ~ i s t = K i , i ∈ [ N ] , s ∈ [ T ] } , \Delta_{1,2}^{K}:=\left\{\widetilde{Y}\in\mathbb{R}_u^{N\times T\times T}\mathrel{\Big|}\sum_{t=s}^{T}\widetilde{Y}_{ist}=K_i,\ i\in[N],\ s\in[T]\right\}, Δ 1 , 2 K := { Y ∈ R u N × T × T t = s ∑ T Y i s t = K i , i ∈ [ N ] , s ∈ [ T ] } , Δ 2 K 0 : = { Z ~ ∈ R u N × T × T ∣ ∑ i = 1 N ∑ t = s T Z ~ i s t = K 0 , s ∈ [ T ] } . \Delta_{2}^{K_0}:=\left\{\widetilde{Z}\in\mathbb{R}_u^{N\times T\times T}\mathrel{\Big|}\sum_{i=1}^{N}\sum_{t=s}^{T}\widetilde{Z}_{ist}=K_0,\ s\in[T]\right\}. Δ 2 K 0 := { Z ∈ R u N × T × T i = 1 ∑ N t = s ∑ T Z i s t = K 0 , s ∈ [ T ] } . the superscript giving each simplex’s total and the subscript the indices held fixed in the summation.
Proposition 3 (The JRP reformulation). Problem (35) is the saddle problem
min X ∈ Δ 1 , 3 1 max Y ~ ∈ Δ 1 , 2 K , Z ~ ∈ Δ 2 K 0 L ( X , Y ~ , Z ~ ) : = ⟨ H + Y ~ + Z ~ , X ⟩ , (36) \min_{X\in\Delta_{1,3}^{1}}\ \max_{\widetilde{Y}\in\Delta_{1,2}^{K},\ \widetilde{Z}\in\Delta_{2}^{K_0}} L(X,\widetilde{Y},\widetilde{Z}):=\langle H+\widetilde{Y}+\widetilde{Z},X\rangle, \tag*{(36)} X ∈ Δ 1 , 3 1 min Y ∈ Δ 1 , 2 K , Z ∈ Δ 2 K 0 max L ( X , Y , Z ) := ⟨ H + Y + Z , X ⟩ , (36) which is (1) with X = Δ 1 , 3 1 \mathcal{X}=\Delta_{1,3}^{1} X = Δ 1 , 3 1 Y = Δ 1 , 2 K × Δ 2 K 0 \mathcal{Y}=\Delta_{1,2}^{K}\times\Delta_{2}^{K_0} Y = Δ 1 , 2 K × Δ 2 K 0 c = vec ( H ) c=\operatorname{vec}(H) c = vec ( H ) q = 0 q=0 q = 0 A : X ↦ ( − X , − X ) A:X\mapsto(-X,-X) A : X ↦ ( − X , − X ) (35) are recovered as Y = Y ~ ⊘ K Y=\widetilde{Y}\oslash K Y = Y ⊘ K Z = Z ~ / K 0 Z=\widetilde{Z}/K_0 Z = Z / K 0
Proof. We prove Proposition 3 in two steps: price the covering rows, and absorb the simplices that pricing produces together with the assignment rows.
Step 1: the dual of each group is a simplex. Price x i s t ≤ y i s x_{ist}\le y_{is} x i s t ≤ y i s Y ~ i s t ≥ 0 \widetilde{Y}_{ist}\ge0 Y i s t ≥ 0 x i s t ≤ z s x_{ist}\le z_s x i s t ≤ z s Z ~ i s t ≥ 0 \widetilde{Z}_{ist}\ge0 Z i s t ≥ 0 (35) reads
L ( X , y , z ; Y ~ , Z ~ ) = ⟨ H , X ⟩ + ∑ i , s y i s ( K i − ∑ t ≥ s Y ~ i s t ) + ∑ s z s ( K 0 − ∑ i ∑ t ≥ s Z ~ i s t ) + ⟨ Y ~ + Z ~ , X ⟩ . L(X,y,z;\widetilde{Y},\widetilde{Z})=\langle H,X\rangle+\sum_{i,s}y_{is}\left(K_i-\sum_{t\ge s}\widetilde{Y}_{ist}\right)+\sum_s z_s\left(K_0-\sum_i\sum_{t\ge s}\widetilde{Z}_{ist}\right)+\langle\widetilde{Y}+\widetilde{Z},X\rangle. L ( X , y , z ; Y , Z ) = ⟨ H , X ⟩ + i , s ∑ y i s ( K i − t ≥ s ∑ Y i s t ) + s ∑ z s ( K 0 − i ∑ t ≥ s ∑ Z i s t ) + ⟨ Y + Z , X ⟩ . Each y i s y_{is} y i s z s z_s z s x ≥ 0 x\ge0 x ≥ 0 ∑ t ≥ s Y ~ i s t = K i \sum_{t\ge s}\widetilde{Y}_{ist}=K_i ∑ t ≥ s Y i s t = K i ( i , s ) (i,s) ( i , s ) ∑ i ∑ t ≥ s Z ~ i s t = K 0 \sum_i\sum_{t\ge s}\widetilde{Z}_{ist}=K_0 ∑ i ∑ t ≥ s Z i s t = K 0 s s s
Step 2: absorb, on both sides. Definition 1 applied to those rows puts Y ~ \widetilde{Y} Y Δ 1 , 2 K \Delta_{1,2}^{K} Δ 1 , 2 K Z ~ \widetilde{Z} Z Δ 2 K 0 \Delta_{2}^{K_0} Δ 2 K 0 (35) puts X X X Δ 1 , 3 1 \Delta_{1,3}^{1} Δ 1 , 3 1 (35) has been absorbed, on both sides: all structure now resides in three products of scaled simplices, and what is left is (36).
Step size. The groups of both dual families partition the primal index set, so Y ~ \widetilde{Y} Y Z ~ \widetilde{Z} Z X X X A : X ↦ ( − X , − X ) A:X\mapsto(-X,-X) A : X ↦ ( − X , − X ) Λ = I \Lambda=I Λ = I
A T A = 2 I , ∥ A ∥ Λ = ∥ A ∥ 2 = 2 , (37) A^{\mathsf{T}}A=2I,\qquad\lVert A\rVert_{\Lambda}=\lVert A\rVert_2=\sqrt{2}, \tag*{(37)} A T A = 2 I , ∥ A ∥ Λ = ∥ A ∥ 2 = 2 , (37) so any η < 1 / 2 \eta<1/\sqrt{2} η < 1/ 2 X i s t X_{ist} X i s t ( i , s ) (i,s) ( i , s ) s s s 2 I 2I 2 I
Per-iteration cost. With Λ = I \Lambda=I Λ = I c = vec ( H ) c=\operatorname{vec}(H) c = vec ( H ) q = 0 q=0 q = 0 A T ( Y ~ , Z ~ ) = − Y ~ − Z ~ A^{\mathsf{T}}(\widetilde{Y},\widetilde{Z})=-\widetilde{Y}-\widetilde{Z} A T ( Y , Z ) = − Y − Z
X k + 1 : = proj Δ 1 , 3 1 ( X k − η ω ( H + Y ~ k + Z ~ k ) ) , X^{k+1}:=\operatorname{proj}_{\Delta_{1,3}^{1}}\left(X^k-\frac{\eta}{\omega}\left(H+\widetilde{Y}^k+\widetilde{Z}^k\right)\right), X k + 1 := proj Δ 1 , 3 1 ( X k − ω η ( H + Y k + Z k ) ) , Y ~ k + 1 : = proj Δ 1 , 2 K ( Y ~ k + η ω ( 2 X k + 1 − X k ) ) , Z ~ k + 1 : = proj Δ 2 K 0 ( Z ~ k + η ω ( 2 X k + 1 − X k ) ) . (38) \widetilde{Y}^{k+1}:=\operatorname{proj}_{\Delta_{1,2}^{K}}\left(\widetilde{Y}^k+\eta\omega\left(2X^{k+1}-X^k\right)\right),\qquad\widetilde{Z}^{k+1}:=\operatorname{proj}_{\Delta_{2}^{K_0}}\left(\widetilde{Z}^k+\eta\omega\left(2X^{k+1}-X^k\right)\right). \tag*{(38)} Y k + 1 := proj Δ 1 , 2 K ( Y k + η ω ( 2 X k + 1 − X k ) ) , Z k + 1 := proj Δ 2 K 0 ( Z k + η ω ( 2 X k + 1 − X k ) ) . (38) which is elementwise: forming H + Y ~ + Z ~ H + \widetilde{Y} + \widetilde{Z} H + Y + Z 2 X k + 1 − X k 2X^{k+1} - X^k 2 X k + 1 − X k O ( N T 2 ) O(NT^2) O ( N T 2 ) N T NT N T T T T Δ 1 , 3 K \Delta^K_{1,3} Δ 1 , 3 K Δ 1 , 2 K \Delta^K_{1,2} Δ 1 , 2 K T T T N T NT N T Δ 2 K 0 \Delta^{K_0}_2 Δ 2 K 0 [11 , 44 ] is linear in its length, so a step is O ( N T 2 ) O(NT^2) O ( N T 2 ) [25 ] costs a logarithmic factor more. The implementation uses active-set methods [40 , 45 ] , which trade the worst-case guarantee for data-dependent control flow that is faster in practice. PDHG applied directly to (35) multiplies twice per iteration by a constraint matrix with roughly 5 2 N T 2 \frac{5}{2}NT^2 2 5 N T 2
Termination. Both iterates are feasible at every step, so the gap is computable in closed form:
L D ( Y ~ , Z ~ ) = ∑ i , t min s ≤ t ( H + Y ~ + Z ~ ) i s t , L P ( X ) = ⟨ H , X ⟩ + ∑ i , s K i max t ≥ s X i s t + ∑ s K 0 max i , t ≥ s X i s t , (39) L_D(\widetilde{Y},\widetilde{Z}) = \sum_{i,t} \min_{s \leq t}(H+\widetilde{Y}+\widetilde{Z})_{ist}, \qquad L_P(X) = \langle H,X\rangle+ \sum_{i,s}K_i\max_{t\geq s}X_{ist} + \sum_s K_0\max_{i,t\geq s}X_{ist}, \tag*{(39)} L D ( Y , Z ) = i , t ∑ s ≤ t min ( H + Y + Z ) i s t , L P ( X ) = ⟨ H , X ⟩ + i , s ∑ K i t ≥ s max X i s t + s ∑ K 0 i , t ≥ s max X i s t , (39) the first placing all mass on the cheapest source period, the second recovering the covering costs of the implied solution of (35) . Both are exact objective values for (35) , so the relative gap the solver monitors is exact.
Theoretical improvements
Theorem 2 (Convergence of PRIA). For the joint replenishment relaxation (35) , run the method of Section 3.1 on (36) at a fixed ω = 1 \omega= 1 ω = 1 and 0 < η < 1 / 2 0 < \eta< 1/\sqrt{2} 0 < η < 1/ 2 on the normalized instance ( vec ( H ) ; K ; K 0 ) / ∥ ( vec ( H ) ; K ; K 0 ) ∥ 2 (\operatorname{vec}(H); K; K_0)/\|(\operatorname{vec}(H); K; K_0)\|_2 ( vec ( H ) ; K ; K 0 ) /∥ ( vec ( H ) ; K ; K 0 ) ∥ 2 v k v^k v k v j , k v^{j,k} v j , k j j j ∥ ⋅ ∥ P \|\cdot\|_P ∥ ⋅ ∥ P
(a) Sublinear convergence without restarts. For every k ≥ 0 k \geq0 k ≥ 0
∥ v k − T ω ( v k ) ∥ P ≤ 2 k + 1 T ( N + 1 ) η . \|v^k-T_\omega(v^k)\|_P \leq\frac{2}{k+1}\sqrt{\frac{T(N+1)}{\eta}}. ∥ v k − T ω ( v k ) ∥ P ≤ k + 1 2 η T ( N + 1 ) . (b) Linear convergence with adaptive restarts. Restarting by sufficient decay (16), every epoch has length at most ⌈ e / α ⌉ \lceil e/\alpha\rceil ⌈ e / α ⌉ α > 0 \alpha> 0 α > 0 j ≥ 0 j \geq0 j ≥ 0
dist P ( v j + 1 , 0 , Z ⋆ ) ≤ e − j α T ( N + 1 ) η . \operatorname{dist}_P(v^{j+1,0},\mathcal{Z}^\star) \leq\frac{e^{-j}}{\alpha}\sqrt{\frac{T(N+1)}{\eta}}. dist P ( v j + 1 , 0 , Z ⋆ ) ≤ α e − j η T ( N + 1 ) . Proof. Every block is a product of scaled simplices, so Λ = I \Lambda= I Λ = I
∥ A ∥ Λ = 2 , diam Λ ( Z ) 2 ≤ 2 N T + 2 T ∥ K ∥ 2 2 + 2 T K 0 2 ≤ 2 N T + 2 T (40) \|A\|_\Lambda= \sqrt{2}, \qquad\operatorname{diam}_\Lambda(\mathcal{Z})^2 \leq2NT + 2T\|K\|_2^2 + 2TK_0^2 \leq2NT + 2T \tag*{(40)} ∥ A ∥ Λ = 2 , diam Λ ( Z ) 2 ≤ 2 N T + 2 T ∥ K ∥ 2 2 + 2 T K 0 2 ≤ 2 N T + 2 T (40) by (37) with A T A = 2 I A^\mathsf{T}A = 2I A T A = 2 I diam Λ ( ∏ g Δ r g ) 2 ≤ 2 ∑ g r g 2 \operatorname{diam}_\Lambda\left(\prod_g \Delta^{r_g}\right)^2 \leq2\sum_g r_g^2 diam Λ ( ∏ g Δ r g ) 2 ≤ 2 ∑ g r g 2 N T NT N T T T T K i K_i K i K 0 K_0 K 0 α > 0 \alpha> 0 α > 0 Z \mathcal{Z} Z ω = 1 \omega= 1 ω = 1
dist P ( v 0 , 0 , Z ⋆ ) < 2 η 2 T ( N + 1 ) = 2 T ( N + 1 ) η . \operatorname{dist}_P(v^{0,0},\mathcal{Z}^\star) < \sqrt{\frac{2}{\eta}}\sqrt{2T(N+1)} = 2\sqrt{\frac{T(N+1)}{\eta}}. dist P ( v 0 , 0 , Z ⋆ ) < η 2 2 T ( N + 1 ) = 2 η T ( N + 1 ) . Part (a) is Proposition 1(a) with this bound, part (b) Proposition 1(b) with this α \alpha α τ 0 ≥ 1 \tau^0 \geq1 τ 0 ≥ 1
Table 8 presents the effect of reformulation on the constants governing convergence. For the standard PDLP formulation, the sublinear bound involves the full constraint matrix through ∥ A ∥ Λ \|A\|_\Lambda ∥ A ∥ Λ R R R ± 1 \pm1 ± 1 ∥ A ∥ Λ \|A\|_\Lambda ∥ A ∥ Λ N N N T T T R R R A : X ↦ ( − X , − X ) A : X \mapsto(-X,-X) A : X ↦ ( − X , − X ) ∥ A ∥ Λ = 2 \lVert A\rVert_{\Lambda} = \sqrt{2} ∥ A ∥ Λ = 2 N N N T T T N N N T T T K K K K 0 K_0 K 0 H H H H K K T H_{\mathrm{KKT}} H KKT ρ \rho ρ α \alpha α N N N T T T (35) . The different residuals and norms specified in the caption must be taken into account when comparing these bounds.
Table 8. Applications of T ω T_\omega T ω ε \varepsilon ε N N N T T T Table 3 . Here α \alpha α
Reformulation also changes the operations within each iteration. Although both formulations have O ( N T 2 ) O(NT^2) O ( N T 2 ) ≈ 5 2 N T 2 \approx\frac{5}{2}NT^2 ≈ 2 5 N T 2 (35) with the negated duplication A : X ↦ ( − X , − X ) A : X \mapsto(-X,-X) A : X ↦ ( − X , − X ) X \mathcal{X} X Y \mathcal{Y} Y Δ 1 , 3 1 \Delta^1_{1,3} Δ 1 , 3 1 Δ 1 , 2 K \Delta^{K}_{1,2} Δ 1 , 2 K Δ 2 K 0 \Delta^{K_0}_{2} Δ 2 K 0 Table 11 ). The two channels of § 4.1 therefore act in opposite directions on the JRP, and the running-time advantage comes through the iteration count.
Computational results
Implementation. The product A X AX A X X X X A T ( Y ~ , Z ~ ) A^{\mathsf{T}}(\widetilde{Y},\widetilde{Z}) A T ( Y , Z ) Δ 1 , 3 1 \Delta^1_{1,3} Δ 1 , 3 1 Δ 1 , 2 K \Delta^{K}_{1,2} Δ 1 , 2 K Δ 2 K 0 \Delta^{K_0}_{2} Δ 2 K 0 [40 ] with warm-started thresholds and a compacted active set for the large per-period segments of Z ~ \widetilde{Z} Z T T T t t t T − 1 − t T-1-t T − 1 − t
Instances. Given N N N T T T N + 1 N+1 N + 1 u ∼ Uniform ( 1 , T ) u \sim\operatorname{Uniform}(1,T) u ∼ Uniform ( 1 , T ) N N N K i K_i K i K 0 K_0 K 0 K 0 K_0 K 0 a i ∼ Uniform ( 0 , 1 ) a_i \sim\operatorname{Uniform}(0,1) a i ∼ Uniform ( 0 , 1 ) ε i s t ∼ Uniform ( 0 , a i ) \varepsilon_{\mathrm{ist}} \sim\operatorname{Uniform}(0,a_i) ε ist ∼ Uniform ( 0 , a i ) s < t s < t s < t H i s t = round ( a i ( t − s ) + ε i s t ) H_{\mathrm{ist}} = \operatorname{round}(a_i(t-s)+\varepsilon_{\mathrm{ist}}) H ist = round ( a i ( t − s ) + ε ist ) H i s t = 0 H_{\mathrm{ist}}=0 H ist = 0 s ≥ t s \ge t s ≥ t
Protocol. We test the eleven shapes of Table 9 , with more repetitions at the cheaper sizes, at a relative optimality tolerance 10 − 6 10^{-6} 1 0 − 6 10 − 9 10^{-9} 1 0 − 9 [35 ] and Gurobi [21 ] , each solving (35) . PRIA and cuPDLPx each ran on one NVIDIA H100 80 GB. Gurobi ran on one 64-core CPU node.
Table 9. The benchmark ladder: eleven shapes over four decades of x x x ≈ 5 2 N T 2 \approx\frac{5}{2}NT^2 ≈ 2 5 N T 2
Results. Table 10 is the running time comparison. PRIA gains on cuPDLPx and Gurobi in solve time at every instance size, by factors of 6.1–28.9 over cuPDLPx. PRIA solved all 363 instances, where cuPDLPx solved 347 and Gurobi 308, so the entries above 10 5 10^5 1 0 5
Table 10. JRP relaxation: solve time in seconds, shifted (by 10 − 2 10^{-2} 1 0 − 2 x x x
Table 11 and Table 12 show that the running-time advantage of PRIA comes from the iteration count alone. Across the five size groups, PRIA takes 1.3–2.9 times as long per iteration as cuPDLPx and converges in 13.6–52.3 times fewer iterations. The higher iteration cost is consistent with projecting onto products of scaled simplices where the standard form clips coordinates, and in particular with the Δ K 0 2 \Delta K_0^2 Δ K 0 2 T T T N T NT N T § 4.1 therefore act against each other here, and the iteration count decides the comparison.
Table 12. Iterations to converge, shifted geometric mean with shift 10. Instances where cuPDLPx failed to converge were not included.
Table 13 shows a corresponding reduction in peak device memory. For the two largest size groups, cuPDLPx uses 3.5–3.8 times as much device memory as PRIA. In particular, for instances with 10 7 10^7 1 0 7 10 8 10^8 1 0 8 ≈ 5 2 N T 2 \approx\frac{5}{2}NT^2 ≈ 2 5 N T 2
Table 13. Largest device memory held during a solve, sampled with NVML.
Conclusions
Reformulation, long considered of limited value for interior-point methods, is a substantial lever for PDHG-type LP solvers, and it acts through two independent channels. On both of our instances, absorbing every constraint into products of projection-friendly polyhedra removed the sparse constraint matrix from the iteration and changed the geometry of the saddle problem. Both reformulated solvers run faster than cuPDLPx on the original formulations at every rung of both ladders (Table 4 and Table 10 ).
We study what to expect from a new instantiation and what to check. The coupling norm is a property of the split, made data-free by the group structure it absorbs. The diameter is computed from the absorbed sets once they are known to be bounded. Sharpness is finite for polyhedral blocks, and bounded when the reformulated polytope carries a lattice of vertex values.
Several directions follow. The recipe applies wherever auxiliary variables select maxima over groups, beyond inventory and binary quadratic models: assignment blocks, ℓ 1 \ell_{1} ℓ 1 Table 1 . A presolve pass that detects embedded projection-friendly blocks and hands a general-purpose PDHG solver the reformulated problem would bring these gains without modelling effort. On the theory side, the restriction to polyhedral blocks and the restriction to bounded ones are likely both removable, at the cost of a normalized progress measure and a Hölderian error bound.
Acknowledgments
We are grateful for the very helpful discussion with Jinwen Yang on the JRP reformulation scheme that initiated this line of research. We thank Edward Rothberg for suggesting, during a visit to MIT, that QUBO would be a good candidate for reformulation. We thank Zedong Peng for the discussions about cuPDLPx. Haihao Lu thanks Oliver Hinder, Warren Schudy, Miles Lubin and David Applegate for their early discussions about reformulation for PDLP.
Jacob Dentes is supported by the National Science Foundation Graduate Research Fellowship Program under Grant No. 2141064. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation. Haihao Lu is partially supported by AFOSR Grant No. FA9550-24-1-0051, ONR Grant No. N000142412735 and a Sloan Research Fellowship.
Sharpness of the QUBO reformulation
This appendix proves (34) , the sharpness constant of Theorem 1. The rates of that theorem are Proposition 1 applied to (33) and (34) .
Proposition 4 (Sharpness of the QUBO reformulation). Let W W W n ≥ 5 n \ge5 n ≥ 5 ω = 1 \omega= 1 ω = 1 0 < η < 1 0 < \eta< 1 0 < η < 1 ( d ; w ) / ∥ ( d ; w ) ∥ 2 (d;w)/\lVert(d;w)\rVert_2 ( d ; w ) / ∥( d ; w ) ∥ 2 ∥ ⋅ ∥ P \lVert\cdot\rVert_P ∥ ⋅ ∥ P (9) at that primal weight. Then
dist P ( v , Z ⋆ ) ≤ ( 1 + 16 η ∥ W ∥ F m n ) ∥ v − T ω ( v ) ∥ P for every v ∈ R n + m , (41) \operatorname{dist}_P(v,Z^\star) \le\left(1+\frac{16}{\eta}\lVert W\rVert_Fm\sqrt{n}\right)\lVert v-T_\omega(v)\rVert_P \qquad\text{for every } v \in\mathbb{R}^{n+m}, \tag*{(41)} dist P ( v , Z ⋆ ) ≤ ( 1 + η 16 ∥ W ∥ F m n ) ∥ v − T ω ( v ) ∥ P for every v ∈ R n + m , (41) where ∥ W ∥ F \lVert W\rVert_F ∥ W ∥ F W W W
Every node must carry an edge, or A λ A_\lambda A λ P P P T ω T_\omega T ω γ \gamma γ (41) .
Vertex margins
Write γ \gamma γ dist Λ ( v , Z ∗ ) ≤ γ gap ( v ) \operatorname{dist}_{\Lambda}(v,\mathcal{Z}^{*}) \le\gamma\operatorname{gap}(v) dist Λ ( v , Z ∗ ) ≤ γ gap ( v ) Z \mathcal{Z} Z γ \gamma γ 1 / α 1/\alpha 1/ α γ \gamma γ
Lemma 1 (Vertex margins give sharpness). Let Q ⊆ R N Q \subseteq\mathbb{R}^{N} Q ⊆ R N g ∈ R N g \in\mathbb{R}^{N} g ∈ R N OPT : = min u ∈ Q ⟨ g , u ⟩ \operatorname{OPT} := \min_{u\in Q}\langle g,u\rangle OPT := min u ∈ Q ⟨ g , u ⟩ F \mathcal{F} F
δ : = min { ⟨ g , u ˉ ⟩ − OPT : u ˉ a vertex of Q with ⟨ g , u ˉ ⟩ > OPT } (42) \delta:= \min\{\langle g,\bar{u}\rangle-\operatorname{OPT}:\bar{u}\text{ a vertex of }Q\text{ with }\langle g,\bar{u}\rangle>\operatorname{OPT}\} \tag*{(42)} δ := min {⟨ g , u ˉ ⟩ − OPT : u ˉ a vertex of Q with ⟨ g , u ˉ ⟩ > OPT } (42) (with δ : = + ∞ \delta:= +\infty δ := + ∞ diam 2 ( Q ) > 0 \operatorname{diam}_{2}(Q)>0 diam 2 ( Q ) > 0
⟨ g , u ⟩ − OPT ≥ δ diam 2 ( Q ) dist 2 ( u , F ) \langle g,u\rangle-\operatorname{OPT}\ge\frac{\delta}{\operatorname{diam}_{2}(Q)}\operatorname{dist}_{2}(u,\mathcal{F}) ⟨ g , u ⟩ − OPT ≥ diam 2 ( Q ) δ dist 2 ( u , F ) for every u ∈ Q u\in Q u ∈ Q
Proof. Write u = ∑ u ˉ μ u ˉ u ˉ u=\sum_{\bar{u}}\mu_{\bar{u}}\bar{u} u = ∑ u ˉ μ u ˉ u ˉ Q Q Q p : = ∑ u ˉ ∈ F μ u ˉ p:=\sum_{\bar{u}\in\mathcal{F}}\mu_{\bar{u}} p := ∑ u ˉ ∈ F μ u ˉ p = 1 p=1 p = 1 u ∈ F u\in\mathcal{F} u ∈ F p < 1 p<1 p < 1 b : = ( 1 − p ) − 1 ∑ u ˉ ∉ F μ u ˉ u ˉ ∈ Q b:=(1-p)^{-1}\sum_{\bar{u}\notin\mathcal{F}}\mu_{\bar{u}}\bar{u}\in Q b := ( 1 − p ) − 1 ∑ u ˉ ∈ / F μ u ˉ u ˉ ∈ Q a : = p − 1 ∑ u ˉ ∈ F μ u ˉ u ˉ ∈ F a:=p^{-1}\sum_{\bar{u}\in\mathcal{F}}\mu_{\bar{u}}\bar{u}\in\mathcal{F} a := p − 1 ∑ u ˉ ∈ F μ u ˉ u ˉ ∈ F p > 0 p>0 p > 0 a a a F \mathcal{F} F p = 0 p=0 p = 0 u = p a + ( 1 − p ) b u=pa+(1-p)b u = p a + ( 1 − p ) b
dist 2 ( u , F ) ≤ ∥ u − a ∥ = ( 1 − p ) ∥ b − a ∥ ≤ ( 1 − p ) diam 2 ( Q ) , \operatorname{dist}_{2}(u,\mathcal{F})\le\lVert u-a\rVert=(1-p)\lVert b-a\rVert\le(1-p)\operatorname{diam}_{2}(Q), dist 2 ( u , F ) ≤ ∥ u − a ∥ = ( 1 − p ) ∥ b − a ∥ ≤ ( 1 − p ) diam 2 ( Q ) , ⟨ g , u ⟩ − OPT = ∑ u ˉ μ u ˉ ( ⟨ g , u ˉ ⟩ − OPT ) ≥ ( 1 − p ) δ . \langle g,u\rangle-\operatorname{OPT}=\sum_{\bar{u}}\mu_{\bar{u}}\bigl(\langle g,\bar{u}\rangle-\operatorname{OPT}\bigr)\ge(1-p)\delta. ⟨ g , u ⟩ − OPT = u ˉ ∑ μ u ˉ ( ⟨ g , u ˉ ⟩ − OPT ) ≥ ( 1 − p ) δ . The first display uses a ∈ F a\in\mathcal{F} a ∈ F a , b ∈ Q a,b\in Q a , b ∈ Q ⟨ g , ⋅ ⟩ \langle g,\cdot\rangle ⟨ g , ⋅ ⟩ F \mathcal{F} F δ \delta δ
L P L_{P} L P L D L_{D} L D L P L_{P} L P Q P Q_{P} Q P x ∈ [ 0 , 1 ] n x\in[0,1]^{n} x ∈ [ 0 , 1 ] n z ∈ [ 0 , 1 ] E z\in[0,1]^{E} z ∈ [ 0 , 1 ] E z z z L P L_{P} L P L P L_{P} L P L D L_{D} L D Y \mathcal{Y} Y g P , g D g_{P},g_{D} g P , g D D P , D D D_{P},D_{D} D P , D D Q P , Q D Q_{P},Q_{D} Q P , Q D δ P , δ D \delta_{P},\delta_{D} δ P , δ D − g D -g_{D} − g D
Lemma 2 (Value polytopes give sharpness). If Q P Q_{P} Q P Q D Q_{D} Q D
γ ≤ max { ∥ Λ X ∥ 2 1 / 2 D P δ P , ∥ Λ Y ∥ 2 1 / 2 D D δ D } . \gamma\le\max\left\{\lVert\Lambda_{X}\rVert_{2}^{1/2}\frac{D_{P}}{\delta_{P}},\lVert\Lambda_{Y}\rVert_{2}^{1/2}\frac{D_{D}}{\delta_{D}}\right\}. γ ≤ max { ∥ Λ X ∥ 2 1/2 δ P D P , ∥ Λ Y ∥ 2 1/2 δ D D D } . with all four quantities measured on the same instance.
Proof. Fix x ∈ X x\in\mathcal{X} x ∈ X u ∈ Q P u\in Q_{P} u ∈ Q P ⟨ g P , u ⟩ = L P ( x ) \langle g_{P},u\rangle=L_{P}(x) ⟨ g P , u ⟩ = L P ( x ) π \pi π x x x π \pi π F P \mathcal{F}_{P} F P X ∗ \mathcal{X}^{*} X ∗ ( x ′ , z ′ ) ∈ F P (x',z')\in\mathcal{F}_{P} ( x ′ , z ′ ) ∈ F P L P ( x ′ ) ≤ ⟨ g P , ( x ′ , z ′ ) ⟩ = OPT ≤ L P ( x ′ ) L_{P}(x')\le\langle g_{P},(x',z')\rangle=\operatorname{OPT}\le L_{P}(x') L P ( x ′ ) ≤ ⟨ g P , ( x ′ , z ′ )⟩ = OPT ≤ L P ( x ′ ) f ∈ F P f\in\mathcal{F}_{P} f ∈ F P dist 2 ( x , X ∗ ) ≤ ∥ π u − π f ∥ ≤ ∥ u − f ∥ \operatorname{dist}_{2}(x,\mathcal{X}^{*})\le\lVert\pi u-\pi f\rVert\le\lVert u-f\rVert dist 2 ( x , X ∗ ) ≤ ∥ π u − π f ∥ ≤ ∥ u − f ∥ f f f dist 2 ( x , X ∗ ) ≤ dist 2 ( u , F P ) \operatorname{dist}_{2}(x,\mathcal{X}^{*})\le\operatorname{dist}_{2}(u,\mathcal{F}_{P}) dist 2 ( x , X ∗ ) ≤ dist 2 ( u , F P ) L P ( x ) − OPT ≥ ( δ P / D P ) dist 2 ( x , X ∗ ) L_{P}(x)-\operatorname{OPT}\ge(\delta_{P}/D_{P})\operatorname{dist}_{2}(x,\mathcal{X}^{*}) L P ( x ) − OPT ≥ ( δ P / D P ) dist 2 ( x , X ∗ ) − g D -g_{D} − g D OPT − L D ( λ ) ≥ ( δ D / D D ) dist 2 ( λ , Y ∗ ) \operatorname{OPT}-L_{D}(\lambda)\ge(\delta_{D}/D_{D})\operatorname{dist}_{2}(\lambda,\mathcal{Y}^{*}) OPT − L D ( λ ) ≥ ( δ D / D D ) dist 2 ( λ , Y ∗ )
Write s : = ( δ P / D P ) dist 2 ( x , X ∗ ) s:=(\delta_{P}/D_{P})\operatorname{dist}_{2}(x,\mathcal{X}^{*}) s := ( δ P / D P ) dist 2 ( x , X ∗ ) t : = ( δ D / D D ) dist 2 ( λ , Y ∗ ) t:=(\delta_{D}/D_{D})\operatorname{dist}_{2}(\lambda,\mathcal{Y}^{*}) t := ( δ D / D D ) dist 2 ( λ , Y ∗ ) C P , C D C_{P},C_{D} C P , C D
gap ( v ) = ( L P ( x ) − OPT ) + ( OPT − L D ( λ ) ) ≥ s + t , \operatorname{gap}(v)=(L_{P}(x)-\operatorname{OPT})+(\operatorname{OPT}-L_{D}(\lambda))\ge s+t, gap ( v ) = ( L P ( x ) − OPT ) + ( OPT − L D ( λ )) ≥ s + t , dist Λ ( v , Z ∗ ) 2 ≤ C P s 2 + C D t 2 ≤ max { C P , C D } ( s + t ) 2 . \operatorname{dist}_{\Lambda}(v,\mathcal{Z}^{*})^{2}\le C_{P}s^{2}+C_{D}t^{2}\le\max\{C_{P},C_{D}\}(s+t)^{2}. dist Λ ( v , Z ∗ ) 2 ≤ C P s 2 + C D t 2 ≤ max { C P , C D } ( s + t ) 2 . The first line adds the two bounds. The second uses that Z ⋆ = X ⋆ × Y ⋆ Z^\star= X^\star\times Y^\star Z ⋆ = X ⋆ × Y ⋆ Λ \Lambda Λ
dist Λ ( v , Z ⋆ ) 2 = dist Λ X ( x , X ⋆ ) 2 + dist Λ Y ( λ , Y ⋆ ) 2 , dist Λ X ( x , X ⋆ ) 2 ≤ ∥ Λ X ∥ 2 dist 2 ( x , X ⋆ ) 2 = C P s 2 , dist Λ Y ( λ , Y ⋆ ) 2 ≤ ∥ Λ Y ∥ 2 dist 2 ( λ , Y ⋆ ) 2 = C D t 2 . \begin{aligned}
\operatorname{dist}_{\Lambda}(v,Z^\star)^2 &= \operatorname{dist}_{\Lambda_X}(x,X^\star)^2 + \operatorname{dist}_{\Lambda_Y}(\lambda,Y^\star)^2,\\
\operatorname{dist}_{\Lambda_X}(x,X^\star)^2 &\le\lVert\Lambda_X\rVert_2\operatorname{dist}_2(x,X^\star)^2 = C_Ps^2,\\
\operatorname{dist}_{\Lambda_Y}(\lambda,Y^\star)^2 &\le\lVert\Lambda_Y\rVert_2\operatorname{dist}_2(\lambda,Y^\star)^2 = C_Dt^2.
\end{aligned} dist Λ ( v , Z ⋆ ) 2 dist Λ X ( x , X ⋆ ) 2 dist Λ Y ( λ , Y ⋆ ) 2 = dist Λ X ( x , X ⋆ ) 2 + dist Λ Y ( λ , Y ⋆ ) 2 , ≤ ∥ Λ X ∥ 2 dist 2 ( x , X ⋆ ) 2 = C P s 2 , ≤ ∥ Λ Y ∥ 2 dist 2 ( λ , Y ⋆ ) 2 = C D t 2 . and then adds the nonnegative cross term 2 max { C P , C D } s t 2\max\{C_P,C_D\}st 2 max { C P , C D } s t dist Λ ( v , Z ⋆ ) ≤ max { C P , C D } gap ( v ) \operatorname{dist}_{\Lambda}(v,Z^\star) \le\sqrt{\max\{C_P,C_D\}}\operatorname{gap}(v) dist Λ ( v , Z ⋆ ) ≤ max { C P , C D } gap ( v )
One hypothesis supplies both margins.
Corollary 1 (A lattice of vertices makes γ \gamma γ Suppose the objective data of both value LPs is integral and every vertex of Q P Q_P Q P Q D Q_D Q D 1 q Z \frac{1}{q}\mathbb{Z} q 1 Z q ≥ 1 q \ge1 q ≥ 1 δ P , δ D ≥ 1 / q \delta_P,\delta_D \ge1/q δ P , δ D ≥ 1/ q
Proof. An integral linear functional evaluated at a point of 1 q Z N \frac{1}{q}\mathbb{Z}^N q 1 Z N 1 q Z \frac{1}{q}\mathbb{Z} q 1 Z 1 q Z \frac{1}{q}\mathbb{Z} q 1 Z 1 / q 1/q 1/ q
Both value polytopes are half-integral
Lemma 3 (Every vertex is half-integral). Every vertex of Q P Q_P Q P ( x , z ) ∈ { 0 , 1 2 , 1 } n + m (x,z) \in\{0,\frac{1}{2},1\}^{n+m} ( x , z ) ∈ { 0 , 2 1 , 1 } n + m
Proof. Let ( x , z ) (x,z) ( x , z ) Q P Q_P Q P B \mathcal{B} B n + m n+m n + m E ⊆ { i < j } E \subseteq\{i<j\} E ⊆ { i < j } 2 2 2
The variable z e z_e z e B \mathcal{B} B z e z_e z e B \mathcal{B} B t e ≥ 1 t_e \ge1 t e ≥ 1 p p p B \mathcal{B} B x x x
∑ e t e + p = n + m . \sum_e t_e + p = n+m. e ∑ t e + p = n + m . Choose one row per edge. Its z e z_e z e ± 1 \pm1 ± 1 z e ∈ { 0 , 1 , x i , x j , x i + x j − 1 } z_e \in\{0,1,x_i,x_j,x_i+x_j-1\} z e ∈ { 0 , 1 , x i , x j , x i + x j − 1 } ( z , x ) (z,x) ( z , x ) ± 1 \pm1 ± 1 B \mathcal{B} B det B = ± det S \det\mathcal{B}=\pm\det S det B = ± det S S S S x x x S S S ∑ e ( t e − 1 ) + p = n \sum_e(t_e-1)+p=n ∑ e ( t e − 1 ) + p = n n n n
Each remaining row of edge e e e z e z_e z e x x x
E + : { cover , z e ≥ 0 } ⇒ x i + x j = 1 , { cover , z e ≤ 1 } ⇒ x i + x j = 2 , E − : { z e ≤ x i , z e ≤ x j } ⇒ x i = x j , { z e ≤ x i , z e ≥ 0 } ⇒ x i = 0 , { z e ≤ x i , z e ≤ 1 } ⇒ x i = 1. \begin{aligned}
E_+:\ &\{\text{cover},\ z_e\ge0\}\Rightarrow x_i+x_j=1,\qquad\{\text{cover},\ z_e\le1\}\Rightarrow x_i+x_j=2,\\
E_-:\ &\{z_e\le x_i,\ z_e\le x_j\}\Rightarrow x_i=x_j,\qquad\{z_e\le x_i,\ z_e\ge0\}\Rightarrow x_i=0,\\
&\{z_e\le x_i,\ z_e\le1\}\Rightarrow x_i=1.
\end{aligned} E + : E − : { cover , z e ≥ 0 } ⇒ x i + x j = 1 , { cover , z e ≤ 1 } ⇒ x i + x j = 2 , { z e ≤ x i , z e ≤ x j } ⇒ x i = x j , { z e ≤ x i , z e ≥ 0 } ⇒ x i = 0 , { z e ≤ x i , z e ≤ 1 } ⇒ x i = 1. The pair { z e ≥ 0 , z e ≤ 1 } \{z_e\ge0,z_e\le1\} { z e ≥ 0 , z e ≤ 1 } x j x_j x j p p p S S S
x j = ϵ x i + c ( ϵ ∈ { ± 1 } , c ∈ Z ) , or x i = c ( c ∈ { 0 , 1 } ) . (43) x_j=\epsilon x_i+c\quad(\epsilon\in\{\pm1\},\ c\in\mathbb{Z}),\qquad\text{or}\qquad x_i=c\quad(c\in\{0,1\}). \tag*{(43)} x j = ϵ x i + c ( ϵ ∈ { ± 1 } , c ∈ Z ) , or x i = c ( c ∈ { 0 , 1 }) . (43) Each edge gives at most one equation of the first kind. Distinct edges join distinct pairs. Those equations therefore form a simple graph on [ n ] [n] [ n ]
Fix a component C C C r r r x r = t x_r=t x r = t x v = ϵ v t + c v x_v=\epsilon_vt+c_v x v = ϵ v t + c v ϵ v ∈ { ± 1 } \epsilon_v\in\{\pm1\} ϵ v ∈ { ± 1 } c v ∈ Z c_v\in\mathbb{Z} c v ∈ Z S S S S S S ∣ C ∣ |C| ∣ C ∣ n = ∑ C ∣ C ∣ n=\sum_C|C| n = ∑ C ∣ C ∣ ∣ C ∣ |C| ∣ C ∣ x v = c x_v=c x v = c t = ϵ v ( c − c v ) ∈ Z t=\epsilon_v(c-c_v)\in\mathbb{Z} t = ϵ v ( c − c v ) ∈ Z ( ϵ v ) v ∈ C (\epsilon_v)_{v\in C} ( ϵ v ) v ∈ C ϵ t + c = − ϵ t + c ′ \epsilon_t+c=-\epsilon_t+c' ϵ t + c = − ϵ t + c ′ 2 t ∈ Z 2t\in\mathbb{Z} 2 t ∈ Z t ∈ 1 2 Z t\in\frac{1}{2}\mathbb{Z} t ∈ 2 1 Z x ∈ 1 2 Z n x\in\frac{1}{2}\mathbb{Z}^n x ∈ 2 1 Z n 0 ≤ x ≤ 1 0\leq x\leq1 0 ≤ x ≤ 1 x ∈ { 0 , 1 2 , 1 } n x\in\{0,\frac{1}{2},1\}^n x ∈ { 0 , 2 1 , 1 } n z e z_e z e 0 , 1 , x i , x j , x i + x j − 1 0,1,x_i,x_j,x_i+x_j-1 0 , 1 , x i , x j , x i + x j − 1 z ∈ [ 0 , 1 ] E z\in[0,1]^E z ∈ [ 0 , 1 ] E z e ∈ 1 2 Z ∩ [ 0 , 1 ] z_e\in\frac{1}{2}\mathbb{Z}\cap[0,1] z e ∈ 2 1 Z ∩ [ 0 , 1 ]
The classical statement [23 ] is about optima. Lemma 3 extends it to all vertices, which is what Corollary 1 needs. The dual side needs a lemma on systems in which every variable touches at most two equations. Such a system can have an exponentially large determinant. The lemma bounds the solution instead, because the blocks are solved independently.
Lemma 4 (Two-nonzero systems). Let A A A nonsingular square matrix with entries in { 0 , ± 1 } \{0,\pm1\} { 0 , ± 1 } A − 1 r ∈ 1 2 Z N A^{-1}r\in\frac{1}{2}\mathbb{Z}^N A − 1 r ∈ 2 1 Z N r r r
Proof. Form the multigraph G G G A A A ∞ \infty ∞ ∞ \infty ∞ A A A C 1 , … , C k C_1,\ldots,C_k C 1 , … , C k G − ∞ G-\infty G − ∞ A = diag ( A 1 , … , A k ) A=\operatorname{diag}(A_1,\ldots,A_k) A = diag ( A 1 , … , A k ) A t A_t A t ( A − 1 r ) ∣ C t = A t − 1 r C t (A^{-1}r)|_{C_t}=A_t^{-1}r_{C_t} ( A − 1 r ) ∣ C t = A t − 1 r C t det A t \det A_t det A t ∣ det A t ∣ ∈ { 1 , 2 } |\det A_t|\in\{1,2\} ∣ det A t ∣ ∈ { 1 , 2 } A t A_t A t + 1 +1 + 1 − 1 -1 − 1 1 \mathbf{1} 1
Fix a block with n t n_t n t H H H G G G C t ∪ { ∞ } C_t\cup\{\infty\} C t ∪ { ∞ } n t + 1 n_t+1 n t + 1 n t n_t n t H H H A t A_t A t ∞ \infty ∞ ∞ \infty ∞ ∣ det A t ∣ = 1 |\det A_t|=1 ∣ det A t ∣ = 1
Now suppose H H H C t C_t C t C t C_t C t ∞ \infty ∞ C t C_t C t n t n_t n t n t n_t n t # V = # E \#V=\#E # V = # E # V = # E \#V=\#E # V = # E 2 # E ≥ 2 # V = 2 # E 2\#E\geq2\#V=2\#E 2# E ≥ 2# V = 2# E
det A t = ± ∏ l s l ± ∏ l t l , s l , t l ∈ { ± 1 } the two entries of column l . \det A_t=\pm\prod_l s_l\pm\prod_l t_l,\qquad s_l,t_l\in\{\pm1\}\text{ the two entries of column }l. det A t = ± l ∏ s l ± l ∏ t l , s l , t l ∈ { ± 1 } the two entries of column l . Hence ∣ det A t ∣ ∈ { 0 , 2 } |\det A_t|\in\{0,2\} ∣ det A t ∣ ∈ { 0 , 2 } 2 2 2
Lemma 5 (Dual margin). Under Assumption 1 , every vertex of Q D Q_D Q D 1 2 Z \frac{1}{2}\mathbb{Z} 2 1 Z
Proof. The absorbed dual on W W W λ e ∈ [ 0 , w e ] \lambda_e\in[0,w_e] λ e ∈ [ 0 , w e ] a e : = λ e a_e:=\lambda_e a e := λ e E − E_- E − b e : = λ e b_e:=\lambda_e b e := λ e E + E_+ E + ( a , b ) (a,b) ( a , b ) m m m
g i ( a , b ) = d i − ∑ e ∈ E − , i = min e a e − ∑ e ∈ E − , i = max e ( w e − a e ) + ∑ e ∈ E + , e ∋ i b e , g_i(a,b)=d_i-\sum_{e\in E_-,\,i=\min e}a_e-\sum_{e\in E_-,\,i=\max e}(w_e-a_e)+\sum_{e\in E_+,\,e\ni i}b_e, g i ( a , b ) = d i − e ∈ E − , i = m i n e ∑ a e − e ∈ E − , i = m a x e ∑ ( w e − a e ) + e ∈ E + , e ∋ i ∑ b e , the dual value LP is max { ∑ i φ i − ∑ e ∈ E + b e } \max\{\sum_i\varphi_i-\sum_{e\in E_+}b_e\} max { ∑ i φ i − ∑ e ∈ E + b e }
Q D = { ( a , b , φ ) : 0 ≤ a ≤ w , 0 ≤ b ≤ w , − B ≤ φ ≤ 0 , φ i ≤ g i ( a , b ) } , B i : = ∣ d i ∣ + ∑ e ∋ i w e . Q_D=\{(a,b,\varphi):0\leq a\leq w,\ 0\leq b\leq w,\ -B\leq\varphi\leq0,\ \varphi_i\leq g_i(a,b)\},\qquad B_i:=|d_i|+\sum_{e\ni i}w_e. Q D = {( a , b , φ ) : 0 ≤ a ≤ w , 0 ≤ b ≤ w , − B ≤ φ ≤ 0 , φ i ≤ g i ( a , b )} , B i := ∣ d i ∣ + e ∋ i ∑ w e . The bound φ ≥ − B \varphi\ge-B φ ≥ − B Q D Q_D Q D φ i = min { 0 , g i ( a , b ) } \varphi_i = \min\{0,g_i(a,b)\} φ i = min { 0 , g i ( a , b )} φ i \varphi_i φ i E − E_- E − g i g_i g i [ − w e , 0 ] [-w_e,0] [ − w e , 0 ] E + E_+ E + [ 0 , w e ] [0,w_e] [ 0 , w e ] ∣ g i ( a , b ) ∣ ≤ B i |g_i(a,b)| \le B_i ∣ g i ( a , b ) ∣ ≤ B i φ i ≥ − B i \varphi_i \ge-B_i φ i ≥ − B i Q D Q_D Q D B i B_i B i g i ( a , b ) = − B i g_i(a,b)=-B_i g i ( a , b ) = − B i
Let ( a , b , φ ) (a,b,\varphi) ( a , b , φ ) Q D Q_D Q D B \mathcal{B} B m + n m+n m + n φ i \varphi_i φ i φ i ≤ 0 \varphi_i \le0 φ i ≤ 0 φ i ≤ g i \varphi_i \le g_i φ i ≤ g i φ i ≥ − B i \varphi_i \ge-B_i φ i ≥ − B i φ i \varphi_i φ i n i n_i n i B \mathcal{B} B φ i \varphi_i φ i p p p ( a , b ) (a,b) ( a , b )
∑ i n i + p = m + n . \sum_i n_i+p=m+n. i ∑ n i + p = m + n . Here n i ≥ 1 n_i\ge1 n i ≥ 1 φ i \varphi_i φ i n i ≤ 2 n_i\le2 n i ≤ 2 0 = φ i = − B i 0=\varphi_i=-B_i 0 = φ i = − B i B i = 0 B_i=0 B i = 0 φ i ≤ 0 \varphi_i\le0 φ i ≤ 0 − φ i ≤ 0 -\varphi_i\le0 − φ i ≤ 0 B \mathcal{B} B i i i φ i \varphi_i φ i ± 1 \pm1 ± 1 ∑ i ( n i − 1 ) + p = m \sum_i(n_i-1)+p=m ∑ i ( n i − 1 ) + p = m ( a , b ) (a,b) ( a , b ) m m m n i = 2 n_i=2 n i = 2
{ φ i ≤ 0 , φ i ≤ g i } ⇒ g i ( a , b ) = 0 , { φ i ≤ g i , φ i ≥ − B i } ⇒ g i ( a , b ) = − B i . \{\varphi_i\le0,\ \varphi_i\le g_i\}\Rightarrow g_i(a,b)=0,\qquad
\{\varphi_i\le g_i,\ \varphi_i\ge-B_i\}\Rightarrow g_i(a,b)=-B_i. { φ i ≤ 0 , φ i ≤ g i } ⇒ g i ( a , b ) = 0 , { φ i ≤ g i , φ i ≥ − B i } ⇒ g i ( a , b ) = − B i . The third pair is the dependent one just excluded. Both right-hand sides are integral. Writing a g g g ( a , b ) (a,b) ( a , b ) − d i + ∑ e ∈ E − , i = max e w e -d_i+\sum_{e\in E_-,\,i=\max e}w_e − d i + ∑ e ∈ E − , i = m a x e w e p p p 0 0 0 w e w_e w e p p p w e = 2 ∣ W i j ∣ > 0 w_e=2|W_{ij}|>0 w e = 2∣ W ij ∣ > 0 a e ≥ 0 a_e\ge0 a e ≥ 0 a e ≤ w e a_e\le w_e a e ≤ w e
Substituting the pins leaves m − p m-p m − p m − p m-p m − p g g g
( ± I p 0 ∗ S ) , \begin{pmatrix}
\pm I_p & 0\\
* & S
\end{pmatrix}, ( ± I p ∗ 0 S ) , so det S = ± det ( residual ) ≠ 0 \det S=\pm\det(\text{residual})\ne0 det S = ± det ( residual ) = 0 { 0 , ± 1 } \{0,\pm1\} { 0 , ± 1 } a e a_e a e b e b_e b e g i g_i g i e e e a e a_e a e − 1 -1 − 1 min e \min e min e + 1 +1 + 1 max e \max e max e b e b_e b e + 1 +1 + 1 ( a , b ) ∈ 1 2 Z m (a,b)\in\frac{1}{2}\mathbb{Z}^m ( a , b ) ∈ 2 1 Z m g i ( a , b ) ∈ 1 2 Z g_i(a,b)\in\frac{1}{2}\mathbb{Z} g i ( a , b ) ∈ 2 1 Z φ i ∈ { 0 , g i ( a , b ) , − B i } ⊆ 1 2 Z \varphi_i\in\{0,g_i(a,b),-B_i\}\subseteq\frac{1}{2}\mathbb{Z} φ i ∈ { 0 , g i ( a , b ) , − B i } ⊆ 2 1 Z
Both propositions are tight. On the primal side a gap of exactly 1 2 \frac{1}{2} 2 1 n = 3 n=3 n = 3 W = ( − 2 3 3 3 − 2 3 3 3 − 3 ) W=\begin{pmatrix}-2&3&3\\3&-2&3\\3&3&-3\end{pmatrix} W = − 2 3 3 3 − 2 3 3 3 − 3 − 2 -2 − 2 − 7 2 -\frac{7}{2} − 2 7 b = ( 1 2 , 3 2 , 3 2 ) b=\left(\frac{1}{2},\frac{3}{2},\frac{3}{2}\right) b = ( 2 1 , 2 3 , 2 3 ) φ = 0 \varphi=0 φ = 0 − 4 -4 − 4
From the gap to the fixed-point residual
The two lemmas of this subsection hold for any complete split; Section 5.2.3 applies the second of them to the JRP. Distances are measured in ∥ ⋅ ∥ Λ \|\cdot\|_{\Lambda} ∥ ⋅ ∥ Λ ∥ ξ ∥ Λ − 1 : = ξ T Λ − 1 ξ \|\xi\|_{\Lambda^{-1}}:=\sqrt{\xi^{\mathsf T}\Lambda^{-1}\xi} ∥ ξ ∥ Λ − 1 := ξ T Λ − 1 ξ λ max \lambda_{\max} λ m a x P P P Λ \Lambda Λ P ⪯ λ max Λ P\preceq\lambda_{\max}\Lambda P ⪯ λ m a x Λ ω = 1 \omega=1 ω = 1 λ max < 2 / η \lambda_{\max}<2/\eta λ m a x < 2/ η
Lemma 6 (Gap bounded by the residual). For every v ∈ Z v\in\mathcal{Z} v ∈ Z gap ( v ) ≤ diam Λ ( Z ) ⋅ dist Λ − 1 ( 0 , F ( v ) ) \operatorname{gap}(v)\le\operatorname{diam}_{\Lambda}(\mathcal{Z})\cdot\operatorname{dist}_{\Lambda^{-1}}(0,F(v)) gap ( v ) ≤ diam Λ ( Z ) ⋅ dist Λ − 1 ( 0 , F ( v )) v v v v + = T ω ( v ) v^+=T_\omega(v) v + = T ω ( v )
dist Λ − 1 ( 0 , F ( v + ) ) ≤ ∥ P ( v + − v ) ∥ Λ − 1 ≤ λ max ∥ v − T ω ( v ) ∥ P . \operatorname{dist}_{\Lambda^{-1}}(0,F(v^+))\le\|P(v^+-v)\|_{\Lambda^{-1}}\le\sqrt{\lambda_{\max}}\|v-T_\omega(v)\|_P. dist Λ − 1 ( 0 , F ( v + )) ≤ ∥ P ( v + − v ) ∥ Λ − 1 ≤ λ m a x ∥ v − T ω ( v ) ∥ P . Consequently gap ( T ω ( v ) ) ≤ λ max diam Λ ( Z ) ∥ v − T ω ( v ) ∥ P \operatorname{gap}(T_\omega(v))\le\sqrt{\lambda_{\max}}\operatorname{diam}_{\Lambda}(\mathcal{Z})\|v-T_\omega(v)\|_P gap ( T ω ( v )) ≤ λ m a x diam Λ ( Z ) ∥ v − T ω ( v ) ∥ P
Proof. Write S v + ϱ : = ( c − A T y , A x − q ) Sv+\varrho:=(c-A^{\mathsf T}y,Ax-q) S v + ϱ := ( c − A T y , A x − q ) F ( v ) = S v + ϱ + N Z ( v ) F(v)=Sv+\varrho+N_{\mathcal{Z}}(v) F ( v ) = S v + ϱ + N Z ( v ) u = ( x ′ , y ′ ) ∈ Z u=(x',y')\in\mathcal{Z} u = ( x ′ , y ′ ) ∈ Z
⟨ S v + ϱ , v − u ⟩ = ⟨ c − A T y , x − x ′ ⟩ + ⟨ A x − q , y − y ′ ⟩ = L ( x , y ) − L ( x ′ , y ) , (44) \langle Sv+\varrho,v-u\rangle=\langle c-A^{\mathsf T}y,x-x'\rangle+\langle Ax-q,y-y'\rangle=L(x,y)-L(x',y), \tag*{(44)} ⟨ S v + ϱ , v − u ⟩ = ⟨ c − A T y , x − x ′ ⟩ + ⟨ A x − q , y − y ′ ⟩ = L ( x , y ) − L ( x ′ , y ) , (44) after cancelling the two copies of ⟨ y , A x ⟩ \langle y, Ax\rangle ⟨ y , A x ⟩ u ∈ X × Y u \in\mathcal{X} \times\mathcal{Y} u ∈ X × Y max y L ( x , ⋅ ) \max_{y} L(x,\cdot) max y L ( x , ⋅ ) − min x L ( ⋅ , y ) -\min_{x} L(\cdot,y) − min x L ( ⋅ , y ) max u ∈ Z ⟨ S v + ϱ , v − u ⟩ = gap ( v ) \max_{u\in\mathcal{Z}} \langle Sv+\varrho,v-u\rangle=\operatorname{gap}(v) max u ∈ Z ⟨ S v + ϱ , v − u ⟩ = gap ( v )
Let ξ ∈ F ( v ) \xi\in F(v) ξ ∈ F ( v ) ∥ ⋅ ∥ Λ − 1 \|\cdot\|_{\Lambda^{-1}} ∥ ⋅ ∥ Λ − 1 ξ = S v + ϱ + ζ \xi=Sv+\varrho+\zeta ξ = S v + ϱ + ζ ζ ∈ N Z ( v ) \zeta\in N_{\mathcal{Z}}(v) ζ ∈ N Z ( v ) u ∈ Z u\in\mathcal{Z} u ∈ Z
⟨ S v + ϱ , v − u ⟩ ≤ ⟨ ξ , v − u ⟩ ≤ ∥ ξ ∥ Λ − 1 ∥ v − u ∥ Λ ≤ ∥ ξ ∥ Λ − 1 diam Λ ( Z ) . \langle Sv+\varrho,v-u\rangle\leq\langle\xi,v-u\rangle\leq\|\xi\|_{\Lambda^{-1}}\|v-u\|_{\Lambda}\leq\|\xi\|_{\Lambda^{-1}}\operatorname{diam}_{\Lambda}(\mathcal{Z}). ⟨ S v + ϱ , v − u ⟩ ≤ ⟨ ξ , v − u ⟩ ≤ ∥ ξ ∥ Λ − 1 ∥ v − u ∥ Λ ≤ ∥ ξ ∥ Λ − 1 diam Λ ( Z ) . The first inequality uses ⟨ ζ , u − v ⟩ ≤ 0 \langle\zeta,u-v\rangle\leq0 ⟨ ζ , u − v ⟩ ≤ 0 ( Λ , Λ − 1 ) (\Lambda,\Lambda^{-1}) ( Λ , Λ − 1 ) u , v ∈ Z u,v\in\mathcal{Z} u , v ∈ Z u u u
For the second, − P ( v + − v ) ∈ F ( v + ) -P(v^{+}-v)\in F(v^{+}) − P ( v + − v ) ∈ F ( v + ) v + = T ω ( v ) v^{+}=T_{\omega}(v) v + = T ω ( v ) P ⪯ λ max Λ P\preceq\lambda_{\max}\Lambda P ⪯ λ m a x Λ P Λ − 1 P ⪯ λ max P P\Lambda^{-1}P\preceq\lambda_{\max}P P Λ − 1 P ⪯ λ m a x P Λ − 1 / 2 \Lambda^{-1/2} Λ − 1/2 M 2 ⪯ λ max M M^{2}\preceq\lambda_{\max}M M 2 ⪯ λ m a x M M : = Λ − 1 / 2 P Λ − 1 / 2 M:=\Lambda^{-1/2}P\Lambda^{-1/2} M := Λ − 1/2 P Λ − 1/2 M ⪯ λ max I M\preceq\lambda_{\max}I M ⪯ λ m a x I ∥ P u ∥ Λ − 1 2 = u T P Λ − 1 P u ≤ λ max u T P u = λ max ∥ u ∥ P 2 \|Pu\|_{\Lambda^{-1}}^{2}=u^{\mathsf{T}}P\Lambda^{-1}Pu\leq\lambda_{\max}u^{\mathsf{T}}Pu=\lambda_{\max}\|u\|_{P}^{2} ∥ P u ∥ Λ − 1 2 = u T P Λ − 1 P u ≤ λ m a x u T P u = λ m a x ∥ u ∥ P 2 T ω ( v ) ∈ Z T_{\omega}(v)\in\mathcal{Z} T ω ( v ) ∈ Z
Lemma 7 (Sharpness of the fixed-point residual). For every v ∈ R n + m v\in\mathbb{R}^{n+m} v ∈ R n + m
dist P ( v , Z ⋆ ) ≤ ( 1 + γ diam Λ ( Z ) λ max ) ∥ v − T ω ( v ) ∥ P . \operatorname{dist}_{P}(v,\mathcal{Z}^{\star})\leq\left(1+\gamma\operatorname{diam}_{\Lambda}(\mathcal{Z})\lambda_{\max}\right)\|v-T_{\omega}(v)\|_{P}. dist P ( v , Z ⋆ ) ≤ ( 1 + γ diam Λ ( Z ) λ m a x ) ∥ v − T ω ( v ) ∥ P . Proof. Write v + : = T ω ( v ) ∈ Z v^{+}:=T_{\omega}(v)\in\mathcal{Z} v + := T ω ( v ) ∈ Z z ⋆ z^{\star} z ⋆ Λ \Lambda Λ Z ⋆ \mathcal{Z}^{\star} Z ⋆ v + v^{+} v +
dist P ( v , Z ⋆ ) ≤ ∥ v − v + ∥ P + ∥ v + − z ⋆ ∥ P ≤ ∥ v − v + ∥ P + λ max dist Λ ( v + , Z ⋆ ) ≤ ∥ v − v + ∥ P + λ max γ gap ( v + ) ≤ ( 1 + γ diam Λ ( Z ) λ max ) ∥ v − v + ∥ P , \begin{aligned}
\operatorname{dist}_{P}(v,\mathcal{Z}^{\star})
&\leq\|v-v^{+}\|_{P}+\|v^{+}-z^{\star}\|_{P}\\
&\leq\|v-v^{+}\|_{P}+\sqrt{\lambda_{\max}}\operatorname{dist}_{\Lambda}(v^{+},\mathcal{Z}^{\star})\\
&\leq\|v-v^{+}\|_{P}+\sqrt{\lambda_{\max}}\gamma\operatorname{gap}(v^{+})\\
&\leq\left(1+\gamma\operatorname{diam}_{\Lambda}(\mathcal{Z})\lambda_{\max}\right)\|v-v^{+}\|_{P},
\end{aligned} dist P ( v , Z ⋆ ) ≤ ∥ v − v + ∥ P + ∥ v + − z ⋆ ∥ P ≤ ∥ v − v + ∥ P + λ m a x dist Λ ( v + , Z ⋆ ) ≤ ∥ v − v + ∥ P + λ m a x γ gap ( v + ) ≤ ( 1 + γ diam Λ ( Z ) λ m a x ) ∥ v − v + ∥ P , where the first inequality is the triangle inequality in ∥ ⋅ ∥ P \|\cdot\|_{P} ∥ ⋅ ∥ P z ⋆ ∈ Z ⋆ z^{\star}\in\mathcal{Z}^{\star} z ⋆ ∈ Z ⋆ P ⪯ λ max Λ P\preceq\lambda_{\max}\Lambda P ⪯ λ m a x Λ v + − z ⋆ v^{+}-z^{\star} v + − z ⋆ γ \gamma γ v + ∈ Z v^{+}\in\mathcal{Z} v + ∈ Z
Proof of Proposition 4
Proof of Proposition 4. The two value polytopes are those of W W W W W W
The rows cutting out Q P Q_{P} Q P x ∈ [ 0 , 1 ] n x\in[0,1]^{n} x ∈ [ 0 , 1 ] n z ∈ [ 0 , 1 ] E z\in[0,1]^{E} z ∈ [ 0 , 1 ] E Q P Q_{P} Q P ∥ ( d ; w ) ∥ 2 \|(d;w)\|_{2} ∥ ( d ; w ) ∥ 2 ∑ i φ i − ∑ e ∈ E + b e \sum_i\varphi_i-\sum_{e\in E_+}b_e ∑ i φ i − ∑ e ∈ E + b e ± 1 \pm1 ± 1 Q D Q_{D} Q D ( d , w , a , b , φ ) (d,w,a,b,\varphi) ( d , w , a , b , φ ) 0 ≤ a , b ≤ w 0\leq a,b\leq w 0 ≤ a , b ≤ w B i B_i B i g i g_i g i Q D Q_{D} Q D Q D / ∥ ( d ; w ) ∥ 2 Q_{D}/\|(d;w)\|_{2} Q D /∥ ( d ; w ) ∥ 2 W W W ∥ ( d ; w ) ∥ 2 \|(d;w)\|_{2} ∥ ( d ; w ) ∥ 2 ∥ ( d ; w ) ∥ 2 \|(d;w)\|_{2} ∥ ( d ; w ) ∥ 2 q = 2 q=2 q = 2 δ P , δ D ≥ 1 / ( 2 ∥ ( d ; w ) ∥ 2 ) \delta_P,\delta_D\geq1/(2\|(d;w)\|_{2}) δ P , δ D ≥ 1/ ( 2∥ ( d ; w ) ∥ 2 )
With ∥ Λ X ∥ 2 = max i deg ( i ) \|\Lambda_X\|_{2}=\max_i\deg(i) ∥ Λ X ∥ 2 = max i deg ( i ) ∥ Λ Y ∥ 2 = 2 \|\Lambda_Y\|_{2}=2 ∥ Λ Y ∥ 2 = 2
γ ≤ 2 ∥ ( d ; w ) ∥ 2 max { max i deg ( i ) D P , 2 D D } ≤ 2 ∥ ( d ; w ) ∥ 2 max { 2 n m , 10 m } = 2 ∥ ( d ; w ) ∥ 2 2 n m ≤ 4 ∥ W ∥ F n m . (45) \begin{aligned}
\gamma
&\leq2\|(d;w)\|_{2}\max\left\{\sqrt{\max_i\deg(i)}D_P,\sqrt{2}D_D\right\}\\
&\leq2\|(d;w)\|_{2}\max\left\{\sqrt{2nm},\sqrt{10m}\right\}\\
&=2\|(d;w)\|_{2}\sqrt{2nm}\leq4\|W\|_{F}\sqrt{nm}. \tag*{(45)}
\end{aligned} γ ≤ 2∥ ( d ; w ) ∥ 2 max { i max deg ( i ) D P , 2 D D } ≤ 2∥ ( d ; w ) ∥ 2 max { 2 nm , 10 m } = 2∥ ( d ; w ) ∥ 2 2 nm ≤ 4∥ W ∥ F nm . (45) The second inequality bounds the two polytopes. The equality uses n ≥ 5 n\geq5 n ≥ 5 ∥ ( d ; w ) ∥ 2 \|(d;w)\|_{2} ∥ ( d ; w ) ∥ 2 ∥ W ∥ F \|W\|_{F} ∥ W ∥ F
Take the primal polytope first. Q P ⊆ [ 0 , 1 ] n + m Q_{P}\subseteq[0,1]^{n+m} Q P ⊆ [ 0 , 1 ] n + m D P 2 ≤ n + m D_{P}^{2}\leq n+m D P 2 ≤ n + m n n n m m m max i deg ( i ) ≤ min { n − 1 , m } \max_i\deg(i)\leq\min\{n-1,m\} max i deg ( i ) ≤ min { n − 1 , m } m ≥ n m\geq n m ≥ n max i deg ( i ) ( n + m ) ≤ ( n − 1 ) ⋅ 2 m ≤ 2 n m \max_i\deg(i)(n+m)\leq(n-1)\cdot2m\leq2nm max i deg ( i ) ( n + m ) ≤ ( n − 1 ) ⋅ 2 m ≤ 2 nm m < n m<n m < n max i deg ( i ) ( n + m ) ≤ m ⋅ 2 n = 2 n m \max_i\deg(i)(n+m)\leq m\cdot2n=2nm max i deg ( i ) ( n + m ) ≤ m ⋅ 2 n = 2 nm max i deg ( i ) D P 2 ≤ 2 n m \max_i\deg(i)D_{P}^{2}\leq2nm max i deg ( i ) D P 2 ≤ 2 nm E − E_{-} E − E + E_{+} E + E E E Q D Q_{D} Q D ∏ e ∈ E [ 0 , w e ] × ∏ i [ − B i , 0 ] \prod_{e\in E}[0,w_{e}]\times\prod_{i}[-B_{i},0] ∏ e ∈ E [ 0 , w e ] × ∏ i [ − B i , 0 ] ∥ ( d ; w ) ∥ 2 \lVert(d;w)\rVert_{2} ∥( d ; w ) ∥ 2
D D 2 ≤ ∥ w ∥ 2 2 + ∑ i B i 2 ∥ ( d ; w ) ∥ 2 2 ≤ 2 ∥ d ∥ 2 2 + ( 1 + 4 max i deg ( i ) ) ∥ w ∥ 2 2 ∥ ( d ; w ) ∥ 2 2 ≤ 1 + 4 max i deg ( i ) ≤ 5 m . D_{D}^{2}\le\frac{\lVert w\rVert_{2}^{2}+\sum_{i}B_{i}^{2}}{\lVert(d;w)\rVert_{2}^{2}}\le\frac{2\lVert d\rVert_{2}^{2}+(1+4\max_{i}\deg(i))\lVert w\rVert_{2}^{2}}{\lVert(d;w)\rVert_{2}^{2}}\le1+4\max_{i}\deg(i)\le5m. D D 2 ≤ ∥( d ; w ) ∥ 2 2 ∥ w ∥ 2 2 + ∑ i B i 2 ≤ ∥( d ; w ) ∥ 2 2 2 ∥ d ∥ 2 2 + ( 1 + 4 max i deg ( i )) ∥ w ∥ 2 2 ≤ 1 + 4 i max deg ( i ) ≤ 5 m . and 2 D D 2 ≤ 10 m 2D_{D}^{2}\le10m 2 D D 2 ≤ 10 m B i 2 ≤ 2 d i 2 + 2 ( ∑ e ∋ i w e ) 2 B_{i}^{2}\le2d_{i}^{2}+2\left(\sum_{e\ni i}w_{e}\right)^{2} B i 2 ≤ 2 d i 2 + 2 ( ∑ e ∋ i w e ) 2 ∑ i ( ∑ e ∋ i w e ) 2 ≤ ∑ i deg ( i ) ∑ e ∋ i w e 2 ≤ 2 max i deg ( i ) ∥ w ∥ 2 2 \sum_{i}\left(\sum_{e\ni i}w_{e}\right)^{2}\le\sum_{i}\deg(i)\sum_{e\ni i}w_{e}^{2}\le2\max_{i}\deg(i)\lVert w\rVert_{2}^{2} ∑ i ( ∑ e ∋ i w e ) 2 ≤ ∑ i deg ( i ) ∑ e ∋ i w e 2 ≤ 2 max i deg ( i ) ∥ w ∥ 2 2 max i deg ( i ) ≥ 1 \max_{i}\deg(i)\ge1 max i deg ( i ) ≥ 1 2 ≤ 1 + 4 max i deg ( i ) 2\le1+4\max_{i}\deg(i) 2 ≤ 1 + 4 max i deg ( i ) ∥ ( d ; w ) ∥ 2 2 = ∥ d ∥ 2 2 + ∥ w ∥ 2 2 \lVert(d;w)\rVert_{2}^{2}=\lVert d\rVert_{2}^{2}+\lVert w\rVert_{2}^{2} ∥( d ; w ) ∥ 2 2 = ∥ d ∥ 2 2 + ∥ w ∥ 2 2 max i deg ( i ) ≤ m \max_{i}\deg(i)\le m max i deg ( i ) ≤ m m ≥ 1 m\ge1 m ≥ 1 1 + 4 max i deg ( i ) ≤ 5 m 1+4\max_{i}\deg(i)\le5m 1 + 4 max i deg ( i ) ≤ 5 m 2 n m ≥ 10 m 2nm\ge10m 2 nm ≥ 10 m n ≥ 5 n\ge5 n ≥ 5
For the conversion, d i = W i i d_{i}=W_{ii} d i = W ii w e = 2 ∣ W i j ∣ w_{e}=2\lvert W_{ij}\rvert w e = 2 ∣ W ij ∣ ∥ W ∥ F 2 = ∥ d ∥ 2 2 + 1 2 ∥ w ∥ 2 2 \lVert W\rVert_{F}^{2}=\lVert d\rVert_{2}^{2}+\frac{1}{2}\lVert w\rVert_{2}^{2} ∥ W ∥ F 2 = ∥ d ∥ 2 2 + 2 1 ∥ w ∥ 2 2 ∥ ( d ; w ) ∥ 2 2 = ∥ d ∥ 2 2 + ∥ w ∥ 2 2 ≤ 2 ∥ W ∥ F 2 \lVert(d;w)\rVert_{2}^{2}=\lVert d\rVert_{2}^{2}+\lVert w\rVert_{2}^{2}\le2\lVert W\rVert_{F}^{2} ∥( d ; w ) ∥ 2 2 = ∥ d ∥ 2 2 + ∥ w ∥ 2 2 ≤ 2 ∥ W ∥ F 2 ∥ ( d ; w ) ∥ 2 ≤ 2 ∥ W ∥ F \lVert(d;w)\rVert_{2}\le\sqrt{2}\lVert W\rVert_{F} ∥( d ; w ) ∥ 2 ≤ 2 ∥ W ∥ F
Every node carries an edge, so deg ( i ) ≥ 1 \deg(i)\ge1 deg ( i ) ≥ 1
diam Λ ( Z ) 2 = ∑ i deg ( i ) + 2 ∥ w ∥ 2 2 ∥ ( d ; w ) ∥ 2 2 = 2 m + 2 ∥ w ∥ 2 2 ∥ ( d ; w ) ∥ 2 2 ≤ 2 m + 2 ≤ 4 m . (46) \operatorname{diam}_{\Lambda}(Z)^{2}=\sum_{i}\deg(i)+\frac{2\lVert w\rVert_{2}^{2}}{\lVert(d;w)\rVert_{2}^{2}}=2m+\frac{2\lVert w\rVert_{2}^{2}}{\lVert(d;w)\rVert_{2}^{2}}\le2m+2\le4m. \tag*{(46)} diam Λ ( Z ) 2 = i ∑ deg ( i ) + ∥( d ; w ) ∥ 2 2 2 ∥ w ∥ 2 2 = 2 m + ∥( d ; w ) ∥ 2 2 2 ∥ w ∥ 2 2 ≤ 2 m + 2 ≤ 4 m . (46) The first equality is the diameter of [ 0 , 1 ] n × ∏ e [ 0 , w e / ∥ ( d ; w ) ∥ 2 ] [0,1]^{n}\times\prod_{e}[0,w_{e}/\lVert(d;w)\rVert_{2}] [ 0 , 1 ] n × ∏ e [ 0 , w e / ∥( d ; w ) ∥ 2 ] ∥ w ∥ 2 ≤ ∥ ( d ; w ) ∥ 2 \lVert w\rVert_{2}\le\lVert(d;w)\rVert_{2} ∥ w ∥ 2 ≤ ∥( d ; w ) ∥ 2 m ≥ 1 m\ge1 m ≥ 1
Combining,
1 + γ diam Λ ( Z ) λ max ≤ 1 + 4 ∥ W ∥ F n m ⋅ 2 m ⋅ 2 η = 1 + 16 η ∥ W ∥ F m n . 1+\gamma\operatorname{diam}_{\Lambda}(Z)\lambda_{\max}\le1+4\lVert W\rVert_{F}\sqrt{nm}\cdot2\sqrt{m}\cdot\frac{2}{\eta}=1+\frac{16}{\eta}\lVert W\rVert_{F}m\sqrt{n}. 1 + γ diam Λ ( Z ) λ m a x ≤ 1 + 4 ∥ W ∥ F nm ⋅ 2 m ⋅ η 2 = 1 + η 16 ∥ W ∥ F m n . The inequality substitutes (44), (45) and λ max < 2 / η \lambda_{\max}<2/\eta λ m a x < 2/ η
Proof of Proposition 1(b)
[38 Theorem 2] carrying the reflection of [38 Section 6] , stated for T ω T_{\omega} T ω T ω T_{\omega} T ω ∥ ⋅ ∥ P \lVert\cdot\rVert_{P} ∥ ⋅ ∥ P Z ⋆ Z^{\star} Z ⋆ P ≻ 0 P\succ0 P ≻ 0 § 3.1 , and the second follows from it by applying [33 ] to 2 T ω − I 2T_{\omega}-I 2 T ω − I [38 ] uses the standard form: [38 Proposition 2] supplies, for each ρ > 0 \rho>0 ρ > 0 { v : ∥ T ω ( v ) ∥ 2 ≤ ρ } \{v:\lVert T_{\omega}(v)\rVert_{2}\le\rho\} { v : ∥ T ω ( v ) ∥ 2 ≤ ρ } v v v ρ \rho ρ [38 ] sets the length of the first one in advance. The argument below is that of [38 ] , recorded for completeness.
Proof of Proposition 1(b). Write τ j \tau^{j} τ j j j j k = τ j k=\tau^{j} k = τ j k = τ j − 1 k=\tau^{j}-1 k = τ j − 1 v j + 1 , 0 = T ω ( v j , τ j ) v^{j+1,0}=T_{\omega}(v^{j,\tau^{j}}) v j + 1 , 0 = T ω ( v j , τ j ) j j j v j , 0 ∈ Z ⋆ v^{j,0}\in Z^{\star} v j , 0 ∈ Z ⋆ τ j = 0 \tau^{j}=0 τ j = 0 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P > 0 \lVert v^{j,0}-T_{\omega}(v^{j,0})\rVert_{P}>0 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P > 0 τ j ≥ 1 \tau^{j}\ge1 τ j ≥ 1 ∥ v j , k − T ω ( v j , k ) ∥ P \lVert v^{j,k}-T_{\omega}(v^{j,k})\rVert_{P} ∥ v j , k − T ω ( v j , k ) ∥ P k k k τ j \tau^{j} τ j
For the epoch length,
1 ϵ ∥ v j , 0 − T ω ( v j , 0 ) ∥ P < ∥ v j , τ j − 1 − T ω ( v j , τ j − 1 ) ∥ P ≤ dist P ( v j , 0 , Z ⋆ ) τ j ≤ 1 α τ j ∥ v j , 0 − T ω ( v j , 0 ) ∥ P , \frac{1}{\epsilon}\lVert v^{j,0}-T_{\omega}(v^{j,0})\rVert_{P}<\lVert v^{j,\tau^{j}-1}-T_{\omega}(v^{j,\tau^{j}-1})\rVert_{P}\le\frac{\operatorname{dist}_{P}(v^{j,0},Z^{\star})}{\tau^{j}}\le\frac{1}{\alpha\tau^{j}}\lVert v^{j,0}-T_{\omega}(v^{j,0})\rVert_{P}, ϵ 1 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P < ∥ v j , τ j − 1 − T ω ( v j , τ j − 1 ) ∥ P ≤ τ j dist P ( v j , 0 , Z ⋆ ) ≤ α τ j 1 ∥ v j , 0 − T ω ( v j , 0 ) ∥ P , where the first inequality is the failure of (16) at k = τ j − 1 k=\tau^{j}-1 k = τ j − 1 k k k v j , 0 v^{j,0} v j , 0 v j , 0 v^{j,0} v j , 0 τ j < ϵ / α \tau^{j}<\epsilon/\alpha τ j < ϵ / α τ j ≤ ⌈ ϵ / α ⌉ \tau^{j}\le\lceil\epsilon/\alpha\rceil τ j ≤ ⌈ ϵ / α ⌉ T ω T_\omega T ω (16) at k = τ i k=\tau^i k = τ i i ≥ 0 i\ge0 i ≥ 0
∥ v i + 1 , 0 − T ω ( v i + 1 , 0 ) ∥ P = ∥ T ω ( v i , τ i ) − T ω ( T ω ( v i , τ i ) ) ∥ P ≤ ∥ v i , τ i − T ω ( v i , τ i ) ∥ P ≤ 1 e ∥ v i , τ i − T ω ( v i , 0 ) ∥ P . (47) \left\|v^{i+1,0}-T_\omega\left(v^{i+1,0}\right)\right\|_P=\left\|T_\omega\left(v^{i,\tau^i}\right)-T_\omega\left(T_\omega\left(v^{i,\tau^i}\right)\right)\right\|_P\le\left\|v^{i,\tau^i}-T_\omega\left(v^{i,\tau^i}\right)\right\|_P\le\frac{1}{e}\left\|v^{i,\tau^i}-T_\omega\left(v^{i,0}\right)\right\|_P. \tag*{(47)} v i + 1 , 0 − T ω ( v i + 1 , 0 ) P = T ω ( v i , τ i ) − T ω ( T ω ( v i , τ i ) ) P ≤ v i , τ i − T ω ( v i , τ i ) P ≤ e 1 v i , τ i − T ω ( v i , 0 ) P . (47) For the distance,
dist P ( v j + 1 , 0 , Z ∗ ) ≤ 1 α ∥ v j + 1 , 0 − T ω ( v j + 1 , 0 ) ∥ P ≤ e − j α ∥ v 1 , 0 − T ω ( v 1 , 0 ) ∥ P ≤ e − j α ∥ v 0 , τ 0 − T ω ( v 0 , τ 0 ) ∥ P ≤ e − j α ( τ 0 + 1 ) dist P ( v 0 , 0 , Z ∗ ) , \begin{aligned}
\operatorname{dist}_P\left(v^{j+1,0},Z^*\right)&\le\frac{1}{\alpha}\left\|v^{j+1,0}-T_\omega\left(v^{j+1,0}\right)\right\|_P\le\frac{e^{-j}}{\alpha}\left\|v^{1,0}-T_\omega\left(v^{1,0}\right)\right\|_P\\
&\le\frac{e^{-j}}{\alpha}\left\|v^{0,\tau^0}-T_\omega\left(v^{0,\tau^0}\right)\right\|_P\le\frac{e^{-j}}{\alpha(\tau^0+1)}\operatorname{dist}_P\left(v^{0,0},Z^*\right),
\end{aligned} dist P ( v j + 1 , 0 , Z ∗ ) ≤ α 1 v j + 1 , 0 − T ω ( v j + 1 , 0 ) P ≤ α e − j v 1 , 0 − T ω ( v 1 , 0 ) P ≤ α e − j v 0 , τ 0 − T ω ( v 0 , τ 0 ) P ≤ α ( τ 0 + 1 ) e − j dist P ( v 0 , 0 , Z ∗ ) , where the first inequality is Definition 2 at v j + 1 , 0 v^{j+1,0} v j + 1 , 0 i = j , … , 1 i=j,\ldots,1 i = j , … , 1 i = 0 i=0 i = 0 k = τ 0 k=\tau^0 k = τ 0 v 0 , 0 v^{0,0} v 0 , 0
Solver parameters 