Optimization

Motivation

Interior-point methods are, together with the simplex method, the standard algorithms for linear programming, and the primal–dual path-following method is the form in which they are implemented in most solvers. Unlike the simplex method, it is a one-phase method: it can start from any point whose primal and dual variables are strictly positive, feasible or not, and drives infeasibility and complementarity to zero simultaneously. The question every user of such a method eventually asks is how fast these three measures of non-optimality decrease.

Chapter 18 of R. J. Vanderbei, Linear Programming: Foundations and Extensions (4th ed., Springer 2014, DOI 10.1007/978-1-4614-7630-6) defines the method from scratch (Fig. 18.1, p. 273) and proves a rate statement, Theorem 18.1 (pp. 277–279): as long as the step lengths stay bounded below and the iterates stay bounded, the primal and dual infeasibilities decay geometrically, and so does the complementarity, at a slower rate. This mission formalizes that theorem and the one-step identities and estimates it is built from. It is the fifth mission of a series on the book; the missions are independent of each other.

Setting

Let $A$ be a real $m \times n$ matrix, $b \in \mathbb{R}^m$, $c \in \mathbb{R}^n$. The primal problem is to maximize $c^Tx$ subject to $Ax + w = b$, $x, w \ge 0$; the dual is to minimize $b^Ty$ subject to $A^Ty - z = c$, $y, z \ge 0$. A primal–dual point is a quadruple $(x, w, y, z)$ with $x, z \in \mathbb{R}^n$, $w, y \in \mathbb{R}^m$; it is strictly positive, $(x, w, y, z) > 0$, if every component is. Write $X, W, Y, Z$ for the diagonal matrices of $x, w, y, z$ and $e$ for the all-ones vector. The norms are $\|v\|_1 = \sum_j |v_j|$ and $\|v\|_\infty = \max_j |v_j|$.

At a point $(x, w, y, z)$ the three measures of progress are the primal infeasibility $\rho = b - Ax - w$, the dual infeasibility $\sigma = c - A^Ty + z$, and the complementarity $\gamma = z^Tx + y^Tw$. Fix parameters $0 < \delta < 1$ and $0 < r < 1$. One iteration of the method, from a strictly positive point, sets $\mu = \delta\gamma/(n+m)$, takes any solution $(\Delta x, \Delta w, \Delta y, \Delta z)$ of the Newton system

$$A\Delta x + \Delta w = \rho, \quad A^T\Delta y - \Delta z = \sigma, \quad Z\Delta x + X\Delta z = \mu e - XZe, \quad W\Delta y + Y\Delta w = \mu e - YWe,$$

computes the step length

$$\theta = r\left(\max_{i,j}\left\{\left|\tfrac{\Delta x_j}{x_j}\right|, \left|\tfrac{\Delta w_i}{w_i}\right|, \left|\tfrac{\Delta y_i}{y_i}\right|, \left|\tfrac{\Delta z_j}{z_j}\right|\right\}\right)^{-1} \wedge 1 \qquad (18.7)$$

(with $\theta = 1$ when all ratios vanish), and moves to $(x + \theta\Delta x, w + \theta\Delta w, y + \theta\Delta y, z + \theta\Delta z)$. This is Fig. 18.1 with the shorter step (18.7) the book adopts for its analysis. Along a sequence of iterates, superscripts $^{(k)}$ denote the quantities at the $k$-th iterate, and $\theta^{(k)}$ is the step length computed there.

Formalization targets

Goal: Theorem 18.1 with the explicit constant

If $t > 0$, $M$ is real, and for all $k \le K$ one has $\theta^{(k)} \ge t$, $\|x^{(k)}\|_\infty \le M$, $\|y^{(k)}\|_\infty \le M$, then for all $k \le K$, with $\tilde t = t(1-\delta)$,

$$\|\rho^{(k)}\|_1 \le (1-t)^k\|\rho^{(0)}\|_1, \qquad \|\sigma^{(k)}\|_1 \le (1-t)^k\|\sigma^{(0)}\|_1, \qquad \gamma^{(k)} \le (1-\tilde t)^k \left(\gamma^{(0)} + \frac{M(\|\rho^{(0)}\|_1 + \|\sigma^{(0)}\|_1)}{\delta t}\right).$$

Milestones

The one-step identities for the infeasibilities, $\tilde\rho = (1-\theta)\rho$ (18.8) and $\tilde\sigma = (1-\theta)\sigma$ (18.9); the one-step complementarity estimate

$$\tilde\gamma \le (1 - (1-\delta)\theta)\gamma + M\|\rho\|_1 + M\|\sigma\|_1 \qquad (18.10)$$

under $\|x\|_\infty, \|y\|_\infty \le M$; and the recursion $\gamma^{(k)} \le (1-\tilde t)\gamma^{(k-1)} + M(1-t)^{k-1}(\|\rho^{(0)}\|_1 + \|\sigma^{(0)}\|_1)$ (18.11). Two unnumbered statements complete the picture: every iteration has $0 < \theta \le 1$ and keeps the point strictly positive, and at any strictly positive point the duality gap satisfies $|b^Ty - c^Tx| \le \gamma + \|\sigma\|_1\|x\|_\infty + \|\rho\|_1\|y\|_\infty$ (§18.5.3).

Significance

Theorem 18.1 separates the convergence question for the path-following method into two parts: a rate statement that holds whenever steps stay long and iterates stay bounded, and the remaining question of when those two conditions hold. It also explains an effect seen in practice: the infeasibilities fall by the factor $1 - t$ per iteration while the complementarity, and hence (by the duality-gap estimate) the gap $b^Ty - c^Tx$, falls only by $1 - \tilde t$. The book stresses that the result is partial, because it does not show that the step lengths remain bounded away from zero; that requires modifications of the method and of the starting point that the book does not carry out.

All statements here are proved in the book. The mission's contribution is a machine-checked version, with the constant of the complementarity bound made explicit. Neither Mathlib nor the platform contains a formal proof of this theorem or a formalization of the infeasible-start primal–dual iteration it concerns; the platform's existing path-following result concerns a different, feasible-start short-step method in equality form.

Difficulty

The infeasibility identities are linear and follow from the first two Newton equations. The complementarity is where the Newton system linearizes a bilinear equation, so the new complementarity contains a second-order term $\theta^2(\Delta y^T\rho - \sigma^T\Delta x)$ that has no sign. Bounding it requires relating the size of the step $\theta\Delta$ to the size of the current iterate through the specific form of the step-length rule (18.7); the rule (18.6) of Fig. 18.1, with signed ratios, does not give such a bound. The multi-step estimate then couples two geometric sequences with different rates, and keeping the constant independent of the horizon $K$ is what makes the statement non-trivial.

Formalization scope

Vectors are Fin n → ℝ and Fin m → ℝ, $A$ is a Matrix (Fin m) (Fin n) ℝ, and points and step directions are a structure PDPoint m n with fields x w y z. The sup-norm is ⨆ j, |v j| (the maximum; 0 for an empty vector). The step length is written with the explicit case $\theta = 1$ when all ratios vanish, since Lean's r / 0 = 0 would otherwise give $\theta = 0$. An iteration is a relation between the current point, a step direction and the next point: the current point is strictly positive, the direction is some solution of the Newton system (uniqueness, which the book asserts under a full-rank assumption, is not assumed), and the next point is current $+\ \theta \cdot$ direction. The hypotheses $0 < \delta < 1$, $0 < r < 1$ (pp. 272–273) are stated in every theorem; $M$ is an arbitrary real number and $K$ a natural number. As in the book, the hypotheses of Theorem 18.1 range over $k \le K$, so the iteration from index $K$ is part of the data.

Explicit constants. The book's Theorem 18.1 asserts only "there exists a constant $\bar M < \infty$". Because $K$ is fixed, that existential is satisfied trivially by $\max_{k \le K}\gamma^{(k)}/(1-\tilde t)^k$, and a statement with $\exists \bar M$ would be empty. The goal therefore uses the constant the book's proof establishes (p. 279, last display): $\bar M = \gamma^{(0)} + M(\|\rho^{(0)}\|_1 + \|\sigma^{(0)}\|_1)/(\delta t)$. Eq. (18.11) is stated with the book's $\tilde M = M(\|\rho^{(0)}\|_1 + \|\sigma^{(0)}\|_1)$ written out.

The formalization needs only finite sums, dot products and matrix–vector products from Mathlib; the definitions of the iteration are reusable for other analyses of the same method (Chapters 19–22 of the book). Contributions are welcome for each milestone separately.

Selected references

• R. J. Vanderbei, Linear Programming: Foundations and Extensions, 4th ed., International Series in Operations Research & Management Science 196, Springer, 2014, Chapter 18, pp. 269–283. https://doi.org/10.1007/978-1-4614-7630-6
• S. J. Wright, Primal-Dual Interior-Point Methods, SIAM, 1997. https://doi.org/10.1137/1.9781611971453

Motivation

Interior-point methods solve linear programs by moving through the interior of the feasible region instead of along its edges, as the simplex method does. The methods used in practice are path-following methods: they track a curve, the central path, that runs through the interior of the feasible region and ends at an optimal solution. Before any such method can be analysed, the curve has to exist. This mission formalizes Chapter 17 of R. J. Vanderbei, Linear Programming: Foundations and Extensions (4th ed., Springer 2014), which defines the central path through the logarithmic barrier problem and proves that it exists exactly when the primal and the dual problem both have strictly positive feasible points.

The chapter's results have a short history. Barrier methods for nonlinear programming go back to Fiacco and McCormick (1968). Interest in interior-point methods for linear programming began with Karmarkar (1984), whose projective algorithm does not mention a central path; the connection between Karmarkar's method and the primal–dual central path was found by Megiddo (1989), with central points traced back to Huard (1967) and an extended study of the path by Bayer and Lagarias (1989). The chapter is the textbook entry point to this line of work and the foundation for the path-following algorithm of Chapter 18.

Setting

Let $A$ be a real $m \times n$ matrix, $b \in \mathbb{R}^m$ and $c \in \mathbb{R}^n$. The primal linear program is to maximize $c^T x$ subject to $Ax \le b$, $x \ge 0$; its dual is to minimize $b^T y$ subject to $A^T y \ge c$, $y \ge 0$. With slack variables $w \in \mathbb{R}^m$ and $z \in \mathbb{R}^n$ they read (17.1)

$$Ax + w = b,\ x, w \ge 0 \qquad\text{and}\qquad A^T y - z = c,\ y, z \ge 0.$$

For a vector $\xi$, $\xi > 0$ means that every component is strictly positive. The primal feasible region has nonempty interior when some $(\bar x, \bar w)$ satisfies $A\bar x + \bar w = b$ with $\bar x > 0$, $\bar w > 0$; the dual feasible region has nonempty interior when some $(\bar y, \bar z)$ satisfies $A^T \bar y - \bar z = c$ with $\bar y > 0$, $\bar z > 0$.

For a parameter $\mu > 0$, the barrier function (17.7) is

$$f(x, w) = c^T x + \mu \sum_{j=1}^n \log x_j + \mu \sum_{i=1}^m \log w_i ,$$

and the barrier problem (17.2) is to maximize $f(x, w)$ subject to $Ax + w = b$, over the domain $x > 0$, $w > 0$ where the logarithms are finite. A solution of the barrier problem is a point of that domain at which $f$ attains its maximum over the domain.

Writing $X, Z, Y, W$ for the diagonal matrices carrying $x, z, y, w$ and $e$ for the all-ones vector, the primal–dual central-path system (17.6) is

$$Ax + w = b, \qquad A^T y - z = c, \qquad XZe = \mu e, \qquad YWe = \mu e,$$

with $x, w, y, z > 0$. The last two equations say $x_j z_j = \mu$ and $y_i w_i = \mu$ for all $j$ and $i$. The set of its solutions $(x_\mu, w_\mu, y_\mu, z_\mu)$, $\mu > 0$, is the primal–dual central path.

The chapter also uses one fact from nonlinear programming: for the problem "maximize $f(x)$ subject to $g_i(x) = 0$, $i = 1, \dots, m$", a critical point is a feasible $x^*$ with $\nabla f(x^*) = \sum_i y_i \nabla g_i(x^*)$ for some Lagrange multipliers $y_i$ (17.3), and $H_f(x^*)$ is the Hessian of $f$ at $x^*$.

Formalization targets

Goal: Theorem 17.2 (p. 265)

For each fixed $\mu > 0$,

$$\exists\, (x, w) \text{ solving the barrier problem} \iff \big(\exists\, \bar x, \bar w > 0 : A\bar x + \bar w = b\big) \wedge \big(\exists\, \bar y, \bar z > 0 : A^T\bar y - \bar z = c\big).$$

Both directions are part of the goal. The statement fixes no constants and no rate; it asserts only when the barrier problem is solvable.

Milestones

1. Theorem 17.1 (p. 261), second-order sufficiency under linear constraints: if the constraints are linear, a critical point $x^*$ with $\xi^T H_f(x^*)\xi < 0$ for every $\xi \ne 0$ satisfying $\xi^T \nabla g_i(x^*) = 0$ for all $i$ is a local maximum on the feasible set.
2. Exercise 10.7 (p. 150): if the primal is feasible and its feasible set $\{x : Ax \le b,\ x \ge 0\}$ is bounded, then there are $y > 0$, $z > 0$ with $A^T y - z = c$.
3. Corollary 17.3 (p. 266): if the primal feasible set (or the dual feasible set) has nonempty interior and is bounded, then for each $\mu > 0$ the system (17.6) has exactly one solution with $x, w, y, z > 0$.

The corollary is stronger than the goal in one direction (it adds uniqueness and the dual variables) and weaker in another (it assumes boundedness).

Significance

The result itself. Theorem 17.2 gives an exact criterion for the barrier problem to be solvable for a fixed $\mu$, and Corollary 17.3 turns it into the statement that the central path is a well-defined curve $\mu \mapsto (x_\mu, w_\mu, y_\mu, z_\mu)$ for all $\mu > 0$. Every path-following method, including the one analysed in Chapter 18 of the same book, targets points on this curve; without existence and uniqueness the "target" of an iteration is undefined. The system (17.6) is also the starting point of the primal–dual Newton step.

Formalizing it. These results are classical and proved in the book; none is formalized in the $Ax \le b$, $x \ge 0$ form used here. The platform has the converse fact in the standard form $Ax = b$, $x \ge 0$ (a solution of the central-path conditions minimizes the barrier, Introduction to Linear Optimization), but not existence. A formal proof of Theorem 17.2 and Corollary 17.3 produces a reusable existence theorem for the central path that downstream missions on path-following and self-dual methods can import.

Difficulty

The "if" direction of Theorem 17.2 is an existence claim on a set that is neither closed nor bounded: the domain $x > 0$, $w > 0$ is open, and the feasible region itself may be unbounded, so the obvious appeal to "a continuous function on a compact set attains its maximum" does not apply directly. The example "maximize $0$ subject to $x \ge 0$" (p. 264), whose barrier $\mu \log x$ has no maximum, shows that the dual hypothesis cannot be dropped. The "only if" direction, which the book calls trivial and does not prove, needs first-order conditions at a maximizer over a relatively open set.

Theorem 17.1 needs a second-order Taylor expansion with a remainder that is $o(\|\xi\|^2)$ uniformly along the constraint subspace, not along individual lines. Exercise 10.7 is a theorem of the alternative and is not a consequence of weak duality alone. Uniqueness in Corollary 17.3 requires the positivity of the solution: the equations $x_j z_j = \mu$, $y_i w_i = \mu$ admit sign-flipped solutions.

Formalization scope

Vectors are Fin n → ℝ and Fin m → ℝ; $A$ is a Matrix (Fin m) (Fin n) ℝ. The book's primal–dual pair in $Ax \le b$, $x \ge 0$ form with slacks is used throughout; there are no explicit constants in this chapter.

Conventions committed to:

• "Nonempty interior" means a feasible point with every component strictly positive, as the proof of Theorem 17.2 says. The topological interior of $\{(x, w) : Ax + w = b,\ x, w \ge 0\}$ in $\mathbb{R}^{n+m}$ is empty whenever $m \ge 1$; reading the theorem that way would make its right-hand side always false for $m \ge 1$, and that reading is ruled out.
• The barrier problem is posed over $x > 0$, $w > 0$ explicitly; Real.log returns $0$ at nonpositive arguments and is never evaluated there.
• Solutions of (17.6) are required to be strictly positive, as in Exercise 17.3 (p. 267).
• "Bounded" is Bornology.IsBounded of the feasible set in $\mathbb{R}^n$ (resp. $\mathbb{R}^m$).
• In Theorem 17.1 the constraints are $Gx = \beta$; $f$ is differentiable near $x^*$ with derivative differentiable at $x^*$, and $\xi^T H_f(x^*) \xi$ is the second Fréchet derivative applied to $(\xi, \xi)$. The local maximum is relative to the feasible set.

Infrastructure a complete development needs: attainment of maxima on compact sets (IsCompact.exists_isMaxOn in Mathlib), first-order conditions on relatively open sets, a theorem of the alternative for Exercise 10.7, and concavity facts about the logarithm. A second-order sufficient condition under affine constraints is not in Mathlib and is reusable beyond linear programming. Proofs of any milestone, including the "only if" half of the goal separately, are welcome contributions.

Selected references

• R. J. Vanderbei, Linear Programming: Foundations and Extensions, 4th ed., International Series in Operations Research & Management Science 196, Springer, 2014, Chapter 17 and Exercise 10.7. https://doi.org/10.1007/978-1-4614-7630-6
• A. V. Fiacco and G. P. McCormick, Nonlinear Programming: Sequential Unconstrained Minimization Techniques, Wiley, 1968. https://doi.org/10.1137/1.9781611971316
• N. Karmarkar, A new polynomial-time algorithm for linear programming, Combinatorica 4 (1984), 373–395. https://doi.org/10.1007/BF02579150
• N. Megiddo, Pathways to the optimal set in linear programming, in Progress in Mathematical Programming, Springer, 1989, 131–158. https://doi.org/10.1007/978-1-4613-9617-8_8
• D. A. Bayer and J. C. Lagarias, The nonlinear geometry of linear programming I, II, Transactions of the AMS 314 (1989), 499–526 and 527–581. https://doi.org/10.1090/S0002-9947-1989-1005525-6

Motivation

The r-algorithm of N. Z. Shor is a variable-metric method for minimizing functions that are not differentiable everywhere, such as maxima of smooth functions, penalty functions and Lagrangian dual functions of integer and decomposition problems. It combines a subgradient step with a space dilation along the difference of two successive almost-gradients, and the book reports it competitive with the best variable-metric and conjugate-gradient methods on "gully"-shaped test problems (Shor 1985, §3.6, pp. 68–77). Its convergence theory is limited. Section 3.7 of Shor's book gives a convergence proof only for an idealized version with exact directional minimization, on a class of piecewise smooth functions, and the book itself calls a weakening of the assumptions and rate estimates desirable (p. 85). This mission formalizes that section.

Timeline. Shor and Zhurbenko introduced the r-algorithm in 1971 (Kibernetika, no. 3, 51–59). Shor analysed the convergence of its exact-line-search version in 1975 (Kibernetika, no. 4, 48–53). Section 3.7 of the 1979 Russian monograph, translated in 1985, gives that analysis for piecewise smooth functions. The book states on p. 85 that weaker assumptions and rate estimates would be desirable.

Setting

Let $E_n$ be $n$-dimensional Euclidean space with inner product $(x, y)$, $n \ge 1$. For a unit vector $\xi$ and $\alpha > 0$ the operator of space dilation is $R_\alpha(\xi)x = x + (\alpha - 1)(x, \xi)\xi$.

Widths. For a compact convex set $W$ and a unit vector $\eta$, the width in direction $\eta$ is $d_\eta(W) = \max_{x \in W}(\eta, x) - \min_{x \in W}(\eta, x)$. The width is $d(W) = \min_{\|\eta\| = 1} d_\eta(W)$ and the diameter is $D(W) = \max_{\|\eta\| = 1} d_\eta(W)$. With $e_\rho(W) = \inf_{z \in W}|(\rho, z)|$ put $K_\rho(W) = d_\rho(W)/e_\rho(W)$ (or $+\infty$ when $e_\rho(W) = 0$) and $p(W) = \inf_{\|\rho\| = 1} K_\rho(W)$.

The class $K$. $E_n$ is partitioned into closed sets $\bar D_1, \dots, \bar D_m$, each the closure of its interior, whose interiors are disjoint and homeomorphic to an open ball or open halfspace; $f_i$ is continuously differentiable on an open set containing $\bar D_i$; and $f = f_i$ on $\bar D_i$. At a point $x$, the set of almost-gradients $G_f(x)$ is the set of gradients $\nabla f_i(x)$ of the pieces with $x \in \bar D_i$. For $\delta, \varepsilon > 0$, $\bar P_{\delta,\varepsilon}(x)$ is the closed convex hull of all $G_f(y)$, $\|y - x\| < \delta$, together with the $\varepsilon$-balls around the elements of $G_f(x)$.

The $r(\alpha)$-algorithm. Fix $\alpha > 1$, $\beta = 1/\alpha$. Start from $x_0$, $\tilde g_0 = 0$ and a nonsingular $B_0$. At iteration $k + 1$ choose $g_f(x_k) \in G_f(x_k)$ with $(B_k^* g_f(x_k), \tilde g_k) \le 0$; set $g_k^* = B_k^* g_f(x_k)$, $r_k = g_k^* - \tilde g_k$, $\xi_{k+1} = r_k/\|r_k\|$, $B_{k+1} = B_k R_\beta(\xi_{k+1})$, $\tilde g_{k+1} = R_\beta(\xi_{k+1}) g_k^*$, and $x_{k+1} = x_k - h_{k+1}B_{k+1}\tilde g_{k+1}$ with $h_{k+1} \ge 0$ such that $f$ does not increase along the step and some almost-gradient at $x_{k+1}$ makes a non-acute angle with $\tilde g_{k+1}$ in the transformed metric. The standing assumption is $f(x) \to +\infty$ as $\|x\| \to \infty$ (3.50).

Formalization targets

Goal: convergence to an isolated local minimum (Theorem 3.13)

Assume (3.50) and $\|x_{k+1} - x_k\| \to 0$ (3.52). If $x^*$ is an isolated local minimum, the component of $\{x : f(x^*) \le f(x) \le f(x_0)\}$ containing $x_0$ also contains $x^*$, and no other point $z$ of that component has linearly dependent $G_f(z)$, then

$$\lim_{k \to \infty} x_k = x^* .$$

Milestones

• Lemma 3.2 (p. 80): for $B = SO$ with minimum eigenvalue $\lambda(B)$ of $S$, $\lambda(B)d(W) \le d(BW) \le \lambda(B)D(W)$.
• Lemma 3.3 (p. 80): for $z_1, z_2 \in W$, $0 < \beta \le 1$ and $\gamma = \|z_1 - z_2\|/d(W) \ge 1$,

$$d\Big(R_\beta\big(\tfrac{z_1 - z_2}{\|z_1 - z_2\|}\big)W\Big) \ge \frac{d(W)}{\sqrt{1 + (1-\beta^2)/(\beta^2\gamma^2)}} .$$

• Theorem 3.11 (p. 82): for every $\sqrt[n]{\beta} < v < 1$, $\varepsilon, \delta > 0$ and $r \ge 1$ there is $\bar k > r$ with

$$p\big(\bar P_{\delta,\varepsilon}(x_{\bar k})\big) \ge \sqrt{\frac{v^2\sqrt[n]{\alpha^2} - 1}{\alpha^2 - 1}} .$$

• Theorem 3.12 (p. 84): the level set $\{f = f_\infty\}$, $f_\infty = \lim_k f(x_k)$, contains a point $x^*$ with $G_f(x^*)$ linearly dependent.

Significance

Theorem 3.13 is the only convergence theorem the book gives for the r-algorithm. Theorems 3.11 and 3.12 hold for any function of the class $K$, convex or not, and say that iterates of the exact version cannot stall at a point where the local almost-gradients are linearly independent: linear dependence of $G_f(x)$ is a generalized stationarity condition that includes $0 \in \operatorname{conv} G_f(x)$. Lemmas 3.2 and 3.3 are statements about widths of convex bodies under linear maps and single dilations, and apply to any analysis of space-dilation methods, including the ellipsoid method.

All four results and the goal are proved in the book, the goal as a corollary without a written proof. None of them is formalized. The formalization adds a machine-checked definition of the class $K$ and of the algorithm as a relation on sequences covering every admissible choice of almost-gradients and stepsizes, a verified proof of the book's argument, and a check of the corollary step, which the book leaves to the reader.

Difficulty

The obvious argument for descent methods is that the function value drops by a fixed amount at every step unless the gradient is small. That argument fails here. The r-algorithm can make null steps, with $x_{k+1} = x_k$ while $B_k$ and $\tilde g_k$ change, and its steps are measured in a metric that degenerates as the dilations accumulate ($\det B_k = \beta^k \det B_0$). A small step does not certify near-stationarity, and a large accumulated dilation does not certify progress.

The difficulty is to relate the geometry of the transformed almost-gradient sets $B_k^*\bar P_{\delta,\varepsilon}(x_k)$ to the decay forced on $B_k$. This needs width estimates for convex bodies under linear maps and single dilations, which Mathlib does not have. Theorem 3.12 also needs a uniform lower bound, near a compact level set, on the Gram determinants of the piece gradients, and the book's proof of it is only sketched. Theorem 3.13 is stated as a corollary with no proof at all. Deriving it requires showing that the iterates cannot leave the prescribed component of $\{f(x^*) \le f \le f(x_0)\}$, and that their accumulation points reduce to $x^*$.

Formalization scope

• $E_n$ is EuclideanSpace ℝ (Fin n) with $n \ge 1$; operators are continuous linear maps; $B_k^*$ is the adjoint. $K_\rho$ and $p$ are valued in $[0, +\infty]$, so $K_\rho(W) = +\infty$ is represented exactly. Widths are sInf/sSup of support values, equal to the book's minima and maxima on nonempty compact sets.
• A function of class $K$ is given with its representation (pieces, open domains, smooth $f_i$), and $G_f(x)$ is the set of gradients of the incident pieces, as on p. 79. Linear dependence of $G_f(x)$ is the negation of linear independence of that set.
• A run of the algorithm is a predicate on sequences $x_k, \tilde g_k, B_k, g_f(x_k), h_{k+1}$; every theorem holds for every run. Step (3) divides by $\|r_k\|$, so $r_k \ne 0$ is part of the run: a run that reaches a point where no admissible choice gives $r_k \ne 0$ has no continuation, as in the book. The stepsize $h_{k+1} = \operatorname{argmin}$ is one admissible choice, not the definition.
• The standing assumption (3.50) is a hypothesis of Theorem 3.11 as well as of 3.12 and 3.13, because the book imposes it for the rest of the section on p. 79.
• In Lemma 3.3, $\beta > 0$ is added (the bound divides by $\beta^2$). "Body" means nonempty interior. An isolated local minimum is a local minimum with a neighbourhood containing no other local minimum.
• Runs exist, so the hypotheses are not vacuous. For $f(t) = |t_1| + 2|t_2|$ (four quadrant pieces), started at the minimizer $x_0 = 0$ with null steps, one can choose $g_f(x_{k+1}) = -g_f(x_k)$ at every step, and $0 \notin G_f$ keeps $r_k \neq 0$. A formalization under which no infinite run exists would make every theorem vacuous. The run predicate also forbids the degenerate reading $h_{k+1} < 0$, under which the monotonicity condition (a) would hold vacuously on an empty segment.
• Needed infrastructure: widths of compact convex sets under linear maps (via the singular values of $B$), the effect of a rank-one dilation on widths, determinant bounds for products of dilations, and Gram-determinant continuity. The first three apply to any space-dilation method and are welcome as independent lemmas.

Selected references

• N. Z. Shor, Minimization Methods for Non-Differentiable Functions, Springer Series in Computational Mathematics 3, Springer 1985, §3.7, pp. 77–85. doi:10.1007/978-3-642-82118-9
• N. Z. Shor, N. G. Zhurbenko, A minimization method using the operation of space dilation in the direction of the difference of two successive gradients, Kibernetika (Kiev), no. 3, 51–59, 1971 (listed in the bibliography of Shor 1985; no online copy known).
• N. Z. Shor, The analysis of convergence of a gradient type method with space dilation in the direction of the difference of two successive gradients, Kibernetika (Kiev), no. 4, 48–53, 1975 (listed in the bibliography of Shor 1985; no online copy known).

Motivation

A scenario program replaces an uncertain optimization problem by its worst case over finitely many sampled instances. In its simplest form it reads

$$\min_{\nu\in\mathbb R^{d-1}}\Big[\max_{i=1,\dots,N}\ell(\nu,\delta_i)\Big],$$

where $\ell(\nu,\delta)$ is the cost of a decision $\nu$ when the uncertain parameter takes the value $\delta$, and $\delta_1,\dots,\delta_N$ are independent draws from an unknown probability $\mathbb P$. The classical guarantee of the scenario approach (Campi and Garatti, 2008) bounds the probability that a new instance produces a cost above the optimal value $\ell^*$, and it does so without any knowledge of $\mathbb P$.

That guarantee concerns a single number, $\ell^*$. Two scenario programs with the same $N$ and the same optimal value can look very different at the solution: in one, most sampled costs lie just below $\ell^*$; in the other, they are widely scattered. The costs that do not determine the solution still carry information about how the cost of the chosen decision is distributed on future instances. Carè, Garatti and Campi (2015) showed that this information can be extracted with the same distribution-free character as the classical result, which is the subject of this mission. It is Chapter 8, §8.1 ("Probability box") of Campi and Garatti, Introduction to the Scenario Approach (SIAM/MOS 2018), the third mission of the series formalizing that book.

Timeline:

• 2008: Campi and Garatti prove that the violation of the scenario solution is dominated by a beta distribution $B(d,N-d+1)$, with equality for fully supported problems (doi:10.1137/07069821X).
• 2015: Carè, Garatti and Campi prove that the risks of all empirical costs from index $d$ on have a joint ordered Dirichlet law (doi:10.1137/130928546).
• 2018: the book states the result as Theorem 8.4 and draws the probability box from it.

Setting

Let $\Delta$ be a measurable space with a probability $\mathbb P$, and $\ell:\mathbb R^{d-1}\times\Delta\to\mathbb R$ a cost that is convex in $\nu$ for every $\delta$ (a standing assumption of the book). For a sample $(\delta_1,\dots,\delta_N)$ of independent draws from $\mathbb P$, let $\nu^*$ be the solution of the program above and $\ell^*=\max_i\ell(\nu^*,\delta_i)$ its optimal value.

Empirical costs (Definition 8.1). Sort the costs of the solution on the sampled scenarios in decreasing order, $\ell^*_1\ge\ell^*_2\ge\dots\ge\ell^*_N$; so $\ell^*_1=\ell^*$.

Risk (Definition 8.2). For a decision $\nu$ and a level $\ell$, $R(\nu,\ell)=\mathbb P\{\delta:\ell(\nu,\delta)>\ell\}$. The risk of the $k$-th empirical cost is $R_k=R(\nu^*,\ell^*_k)$, and $R_1\le R_2\le\dots\le R_N$.

Nondegeneracy (Definition 8.3). For every $N\ge d$, with probability $1$, $\ell^*_d\ne\ell^*_{d+1}\ne\dots\ne\ell^*_N$. Costs with index below $d$ are excluded because several scenarios typically attain the maximum at $\nu^*$.

Support constraints and full support (Definitions 5.1 and 5.4). In epigraph form, $\min t$ subject to $t\ge\ell(\nu,\delta_i)$, the constraint of scenario $i$ is a support constraint if removing it lowers the optimal value; the problem is fully supported if for every $m\ge d$ the program with $m$ scenarios has exactly $d$ support constraints with probability $1$.

The ordered Dirichlet distribution with parameters $(d,1,\dots,1)$ is the law on $\{0\le\alpha_d\le\dots\le\alpha_N\le1\}$ with density $\frac{N!}{(d-1)!}\alpha_d^{d-1}$.

Formalization targets

Goal: Theorem 8.4

Under nondegeneracy, for $N\ge d$ and all $\varepsilon_d,\dots,\varepsilon_N$,

$$\mathbb P^N\{R_d\le\varepsilon_d,\dots,R_N\le\varepsilon_N\}=\frac{N!}{(d-1)!}\int_0^{\varepsilon_d}\alpha_d^{d-1}\int_0^{\varepsilon_{d+1}}\!\!\cdots\int_0^{\varepsilon_N}\mathbf 1_{\{0\le\alpha_d\le\dots\le\alpha_N\le1\}}\,\mathrm d\alpha_N\cdots\mathrm d\alpha_d .$$

This is an identity of joint distribution functions, not a bound, and it does not depend on $\ell$ or $\mathbb P$.

Milestones

1. Theorem 3.7 for the min-max program: $\mathbb P^N\{R(\nu^*,\ell^*)>\varepsilon\}\le\sum_{i=0}^{d-1}\binom Ni\varepsilon^i(1-\varepsilon)^{N-i}$ for $\varepsilon\in[0,1]$.
2. Fully supported problems: $\ell^*=\ell^*_d$ with probability $1$.
3. Marginal of $R_d$ (a corollary of the goal): $\mathbb P^N\{R_d\le\varepsilon\}=1-\sum_{i=0}^{d-1}\binom Ni\varepsilon^i(1-\varepsilon)^{N-i}$, the beta law $B(d,N-d+1)$.

Significance

The result. The theorem controls the whole distribution function of the cost $\ell(\nu^*,\delta)$ of the scenario solution on a new instance, not just one quantile of it. Discarding the extreme tails of the laws of $R_d,\dots,R_N$ yields, with a prescribed confidence $1-\beta$, a region (the book's "probability box", Figure 8.2) that contains the entire cumulative distribution function of $\ell(\nu^*,\delta)$, computed from the sample alone. The first marginal recovers the classical Theorem 3.7, since $\ell^*\ge\ell^*_d$ makes the risk of $\ell^*$ at most $R_d$.

Formalizing it. The theorem is proved in Carè, Garatti and Campi (2015); the book states it and gives no proof. No machine-checked version of the scenario approach, of its generalization theorem, or of ordered Dirichlet laws of risks is known to exist. A formalization would produce a checked proof of the distribution-free identity together with the combinatorial and measure-theoretic infrastructure (order statistics of sampled costs, laws of random risks) that the rest of scenario theory reuses. The milestones separate the classical beta bound, which is also the goal of the first mission of this series, from the new exact joint law.

Difficulty

The obvious attempt treats $R_d,\dots,R_N$ as the order statistics of the uniform variables $1-F(\ell(\nu^*,\delta_i))$. That works only for $d=1$, when the decision space is a point and the costs are independent. For $d\ge2$ the decision $\nu^*$ is itself a function of the whole sample, so the sampled costs at $\nu^*$ are neither independent nor identically distributed, and the $d$-th cost is tied to the scenarios that determine the solution. The factor $\alpha_d^{d-1}$ and the constant $N!/(d-1)!$ encode exactly this dependence. Any argument must account for which scenarios are active at $\nu^*$ without assuming full support, since the theorem holds whether or not $\ell^*=\ell^*_d$.

Formalization scope

The decision space is EuclideanSpace ℝ (Fin n) and the book's $d$ is $n+1$; the sample is ω : Fin N → Δ with law Measure.pi (fun _ => P). Empirical costs are read from Tuple.sort with $k$ counted from $1$; risks are real numbers (P {δ | c < ℓ ν δ}).toReal. The right-hand side of the goal is a Lebesgue integral over the box $\prod_k[0,\varepsilon_k]$ intersected with the ordered simplex, stated for all real $\varepsilon_k$. The solution map $\omega\mapsto\nu^*$ is a hypothesis-constrained function, never an arbitrary map.

Implicit hypotheses of the book pinned down in the binders:

• $\ell(\cdot,\delta)$ is convex for every $\delta$ (p. 6).
• Existence and uniqueness of the solution of the program for every sample size $m\ge1$ and every sample; the book's Assumption 3.6 says "every $m$", but the program with no scenario has no solution.
• Nondegeneracy for every sample size $m\ge d$, not only for the $N$ of the theorem, as Definition 8.3 is written.
• Joint measurability of $(\nu,\delta)\mapsto\ell(\nu,\delta)$ and measurability of the solution map (measurability is glossed over in the book, p. 6 footnote 1 and p. 33).
• $N\ge d$, and $\varepsilon\in[0,1]$ in the binomial-form statements.

A statement in which the solution is an arbitrary measurable map, or in which the nondegeneracy or existence hypothesis is unsatisfiable, would make the goal vacuous; the hypotheses here are met, for example, by $\ell(\nu,\delta)=\|\nu-\delta\|^2$ with a continuous law on $\mathbb R^{n}$, and for $d=1$ by any cost independent of $\nu$ with an atomless law.

A complete development needs: the scenario approach generalization theorem (reusable across this series), laws of order statistics of i.i.d. uniform variables, and the combinatorics of support sets of convex min-max programs. Proofs of the milestones, of the $d=1$ case of the goal, and of auxiliary facts about kthLargest are all welcome.

Selected references

• M. C. Campi, S. Garatti, Introduction to the Scenario Approach, MOS-SIAM Series on Optimization 26, SIAM, 2018. doi:10.1137/1.9781611975444
• A. Carè, S. Garatti, M. C. Campi, Scenario min-max optimization and the risk of empirical costs, SIAM J. Optim. 25(4):2061–2080, 2015. doi:10.1137/130928546
• M. C. Campi, S. Garatti, The exact feasibility of randomized solutions of uncertain convex programs, SIAM J. Optim. 19:1211–1230, 2008. doi:10.1137/07069821X

Motivation

Small and medium-sized suppliers that sell to large retailers typically ship first and are paid weeks or months later. Until the invoice is paid, the supplier's cash is trapped in accounts receivable, while production has to be financed up front. Factoring (selling or pledging the receivable to a financial intermediary for immediate cash) and reverse factoring (a buyer-initiated program in which the retailer approves the invoice and the supplier is paid early by the retailer's bank) are the standard remedies, and reverse factoring programs of large buyers are often coupled with an extension of payment terms. Which of these instruments a supplier should use, and how the answer depends on the credit ratings of the two firms, is the question addressed by Kouvelis and Xu, "A Supply Chain Theory of Factoring and Reverse Factoring", Management Science 67(10):6071–6088, 2021 (doi:10.1287/mnsc.2020.3788). The paper embeds each financing scheme in the pull wholesale-price game of Cachon (2004) and compares the resulting equilibria.

Setting

A retailer (the Stackelberg leader) offers a wholesale price $w$ to a capital-constrained supplier (the follower), who then chooses a production quantity $q \ge 0$ and carries the inventory; the retail price $p$ is fixed, and $0 < c < p$ is the unit production cost. Demand $D \ge 0$ has density $f$, complementary distribution function $\bar F$, finite mean, $f > 0$ on $[0, \mathbb Z]$ ($\mathbb Z \le +\infty$), and a strictly increasing failure rate $z = f/\bar F$ (strict IFR). Expected sales are $S(q) = \int_0^q \bar F(\xi)\,d\xi$, and $k(q) = S(q)/\bar F(q)$.

Each firm $j \in \{s, r\}$ has a credit rating $C_j \in (C_{\min}, C_{\max})$. The default probability $\rho_j = \rho(C_j) \in [0,1]$ is strictly decreasing in the rating, and the interest-rate premium $\eta_j = \eta(C_j) > 0$ is decreasing. Production takes a lead time $t_1$, payment follows after a term $t_2$, and the supplier faces a Poisson liquidity shock of rate $\lambda_s$ during production (the retailer, rate $\lambda_r$). Under a scheme $i \in \{\mathcal F, \mathcal N, \mathcal R\}$ (recourse, non-recourse, reverse factoring, the last with a payment extension $\tau \ge 0$), the supplier's expected profit is

$$\pi_i(q; w) = (1-\rho_s)\big(\Lambda_i e^{-\lambda_s t_1} w S(q) - c q e^{\eta_s t_1}\big),$$

with coefficients $\Lambda_{\mathcal F} = (1-\rho_r) + (1-\rho_s) - e^{\eta_s t_2}$, $\Lambda_{\mathcal N} = e^{-\eta_r t_2}(1-\rho_r)$, $\Lambda_{\mathcal R} = e^{-\eta_r(t_2+\tau)}$, and effective unit cost $c_i = c\,e^{(\eta_s+\lambda_s)t_1}/\Lambda_i$. The retailer's profit is $e^{-\lambda_s t_1}(1-\rho_r)(p-w)S(q)$ under factoring and carries the extra factor $2 - e^{-\lambda_r \tau}$ under reverse factoring.

Formalization targets

Goal: Proposition 5

For a given $\tau \ge 0$, let $\mathbb C_3$ be the unique supplier rating with $\Lambda_{\mathcal F} = \Lambda_{\mathcal R}$, $\mathbb C_1$ the one with $\Lambda_{\mathcal F} = \Lambda_{\mathcal N}$, and $\mathbb C_{\mathcal F}, \mathbb C_{\mathcal N}, \mathbb C_{\mathcal R}$ the feasibility thresholds ($c_i = p$). With all three schemes available,

$$e^{-\eta_r\tau} < 1-\rho_r:\quad \mathcal F \text{ adopted} \iff C_s > \mathbb C_{\mathcal F}\vee\mathbb C_1,\qquad \mathcal N \text{ adopted} \iff \mathbb C_{\mathcal N} < C_s \le \mathbb C_1;$$ $$1-\rho_r \le e^{-\eta_r\tau}:\quad \mathcal F \text{ adopted} \iff C_s > \mathbb C_{\mathcal F}\vee\mathbb C_3,\qquad \mathcal R \text{ adopted} \iff \mathbb C_{\mathcal R} < C_s \le \mathbb C_3.$$

Milestones

1. Proposition 2 (p. 6078): recourse factoring is feasible iff $C_s > \mathbb C_{\mathcal F}$, and then the unique equilibrium solves $p\bar F(q) = c_{\mathcal F}[1+z(q)k(q)]$, $w = c_{\mathcal F}/\bar F(q)$.
2. Proposition 3 (p. 6079): the same for non-recourse factoring with $c_{\mathcal N}$.
3. §5.2 (p. 6082, Proposition EC.2): the same for reverse factoring with $c_{\mathcal R}(\tau)$.
4. §4.4 (p. 6080, Lemma EC.2): between feasible schemes, the supplier's equilibrium profits are ordered as the $\Lambda_i$.
5. Proposition 4 (p. 6080): the choice between $\mathcal F$ and $\mathcal N$ alone.

A follow-on item states Corollary 1(i) (p. 6081): when the retailer's rating is at least the supplier's, non-recourse factoring strictly dominates recourse factoring whenever the latter is feasible.

Significance

Proposition 5 is the paper's main prediction: non-recourse factoring is the choice of medium-rated suppliers, recourse factoring of highly rated ones, and reverse factoring with a fixed payment extension is chosen only by low-to-medium-rated suppliers and only if the retailer's rating is not too high relative to the extension. The paper uses it to explain the documented reluctance of suppliers to join reverse-factoring programs that come with extended payment terms (Wuttke, Rosenzweig and Heese 2019; Corsten 2010, a teaching case), and it is the starting point of the paper's Proposition 6 on the retailer's optimal payment extension, which is the subject of the second mission of this series.

The results are proved in the paper's Online Appendix B. To the best of available knowledge none of them, and no model of supply chain finance, has a machine-checked proof. A formal development would also supply a reusable treatment of the pull wholesale-price game under a strictly IFR demand with an arbitrary effective unit cost, of which Propositions 2, 3 and EC.2 are three instances.

Difficulty

The case analysis of Proposition 5 is short once the milestones are available; the substance lies in them. Proposition 2 requires solving a bilevel problem: the follower's best response has to be identified for every wholesale price, including prices at which the supplier produces nothing, and the leader's objective has to be shown to have a single maximizer. The natural route through the first-order condition needs the IFR property to exclude multiple stationary points, and it has to cope with a support end $\mathbb Z$ that may be finite (where $\bar F$ vanishes and $k$, $z$ lose meaning) or infinite. Lemma EC.2 requires the equilibrium supplier profit as an explicit function of the effective cost and a monotonicity argument for it; comparing the $\Lambda_i$ directly, without the equilibrium, does not prove it. The proofs themselves are not in the main text.

Formalization scope

Demand is a probability measure μ on ℝ with a density f, the support end Z : EReal, and the complementary distribution function Fbar μ x = μ((x, ∞)). The model data, the coefficients coef, the effective costs cF, cN, cR, the two profit functions, best responses, equilibria, feasibility and adoption are definitions; all theorems are for arbitrary data satisfying Params.Valid and DemandModel. Readings fixed by the formalization:

• "Feasible": some wholesale price $w \ge 0$ with a supplier best response gives the retailer strictly positive profit. It is not defined as $c_i < p$, which is what Propositions 2 and 3 prove.
• "Adopted" / "should be adopted": the scheme is feasible and its equilibrium supplier profit is best among the feasible available schemes, with ties resolved $\mathcal R \succ \mathcal N \succ \mathcal F$, the order forced by the paper's boundaries ($C_s = \mathbb C_1$, $C_s = \mathbb C_3$, $C_r = \mathbb C_2$). Profits are compared, not the $\Lambda_i$.
• "The unique value of $C_s$ that satisfies …": a point of $(C_{\min}, C_{\max})$ satisfying the equation and the only such point. The equation $c_{\mathcal F} = p$ is cross-multiplied, $c\,e^{(\eta_s+\lambda_s)t_1} = p\Lambda_{\mathcal F}$, since $\Lambda_{\mathcal F}$ may be $\le 0$.
• "$C_r > \mathbb C_2$" is read through the paper's equivalence $\tau > -\eta_r^{-1}\ln(1-\rho_r)$, i.e. $e^{-\eta_r\tau} < 1-\rho_r$: the stated assumptions give no single crossing in $C_r$.
• "The unique equilibrium can be derived from …": the equilibrium exists and is unique, and it is exactly the solution of the system with $q \in (0, \mathbb Z)$; this part is stated under feasibility.
• "Increasing/decreasing" is weak unless stated (p. 6075); $\rho$ is strictly decreasing, $\eta$ weakly.
• Continuity of $f$ is read on $[0,\mathbb Z]$, and $\mathbb Z$ as the upper end of the support.
• Additions: wholesale prices are nonnegative; the retailer's reverse-factoring profit is the last line of the p. 6082 display (the middle line has a stray factor $w$).
• Corollary 1(i): "higher credit rating" is read as $C_r \ge C_s$, "dominates" as a strictly larger equilibrium supplier profit whenever recourse is feasible.

A formalization in which feasibility or adoption is defined by the inequalities being proved ($c_i < p$, largest $\Lambda_i$) makes the propositions trivial and is ruled out by the definitions above. The profit functions (7), (10) and the p. 6082 display are taken as the model; the derivation from bank and factor pricing (Eqs. (1), (5), (6), Lemma 1) is not formalized, and the unstated assumptions of Table EC.2 are not included. Contributions of reusable lemmas on $S$, $k$ and $z$ under strict IFR are welcome.

Selected references

• P. Kouvelis and F. Xu, A Supply Chain Theory of Factoring and Reverse Factoring, Management Science 67(10):6071–6088, 2021. https://doi.org/10.1287/mnsc.2020.3788
• G. P. Cachon, The Allocation of Inventory Risk in a Supply Chain: Push, Pull, and Advance-Purchase Discount Contracts, Management Science 50(2):222–238, 2004. https://doi.org/10.1287/mnsc.1030.0189
• D. A. Wuttke, E. Rosenzweig, H. S. Heese, An Empirical Analysis of Supply Chain Finance Adoption, Journal of Operations Management 65(3):242–261, 2019.
• P. Kouvelis and W. Zhao, Who Should Finance the Supply Chain? Impact of Credit Ratings on Supply Chain Decisions, Manufacturing & Service Operations Management 20(1):19–35, 2018.

Motivation

Vehicle routing problems ask for minimum-cost vehicle routes that deliver goods from a depot to a set of customers. In the split-delivery vehicle routing problem with time windows (SDVRPTW) a customer's demand may be served by several vehicles, and each customer must be visited inside a prescribed time window. Allowing split deliveries matters in practice: for the variant without time windows, Dror and Trudeau (1989) showed empirically that splitting can save substantially, and Archetti, Savelsbergh and Speranza (2006) proved the savings can reach 50%.

Exact methods for the SDVRPTW are branch-and-price algorithms, and they rely on structural properties of optimal solutions to prune the search. Desaulniers (Operations Research 58(1), 2010) collects these properties in §2 of his paper and uses the strongest one, Corollary 2, as a family of valid inequalities (constraint (7)) in his branch-and-price-and-cut method.

Timeline (as reported by Desaulniers 2010, §§1–2).

• 1989–1990. Dror and Trudeau prove, for the SDVRP without time windows and with the triangle inequality, that some optimal solution has no two routes sharing more than one split customer.
• 2006. Gendreau, Dejax, Feillet and Gueguen observe that the property holds with time windows, derive the arc-based Corollary 1, and remark that elementary routes suffice (Remark 1).
• 2010. Desaulniers strengthens Corollary 1 to pairs of reverse arcs (Corollary 2) and exploits it as cutting planes.

Setting

An instance has $n$ customers $\mathcal N$, a start depot $0$ and an end depot $n+1$ (the same location at the beginning and end of the planning horizon), and the following data: a vehicle capacity $Q > 0$; a demand $d_i > 0$ for each customer, which may exceed $Q$; a time window $[e_v, l_v]$ for each node, shared by the two depot copies; nonnegative travel times $t_{vw}$, which include the service time at $v$; and nonnegative costs $c_{vw}$. The arc set $\mathcal A$ contains the idle arc $(0, n+1)$ and every arc $(v, w)$, $v \ne w$, with $e_v + t_{vw} \le l_w$. The triangle inequality $t_{vx} \le t_{vw} + t_{wx}$, $c_{vx} \le c_{vw} + c_{wx}$ is assumed throughout.

A route is a walk $0 \to v_1 \to \dots \to v_m \to n+1$ along arcs of $\mathcal A$, customers possibly repeated, with service start times inside the time windows that respect travel times (waiting is allowed), and nonnegative quantities delivered at its visits whose total is at most $Q$. Its cost is the sum of its arc costs. A solution is a finite family of routes, one per vehicle, with no bound on their number; it is feasible if every customer $i$ receives in total at least $d_i$, and optimal if it is feasible and no feasible solution costs less. For customers $i, j$, let $x_{ij}$ be the number of times arc $(i, j)$ is traversed, summed over all routes of a solution, and let $\mathcal A(\mathcal N) = \mathcal A \cap (\mathcal N \times \mathcal N)$.

Formalization targets

All four statements assume the triangle inequality and that the instance has a feasible solution, and assert the existence of an optimal solution with a structural property.

Goal: Corollary 2 (p. 181)

$$\exists \text{ optimal solution with } x_{ij} + x_{ji} \le 1 \quad \text{for all } (i,j) \in \mathcal A(\mathcal N).$$

The two arcs of a pair of reverse customer arcs are used at most once in total. This is the form constraint (7) of the paper gives to the corollary.

Milestones

1. Remark 1. Some optimal solution has only elementary routes: no route visits a customer twice.
2. Theorem 1. Some optimal solution has no two distinct routes with two customers in common.
3. Corollary 1. Some optimal solution has $x_{ij} \le 1$ for all $(i, j) \in \mathcal A(\mathcal N)$.

Significance

The result. Corollary 2 turns a property of optimal solutions into linear inequalities on arc-flow variables. Adding them to the arc-flow formulation cuts off fractional points of its linear relaxation while keeping an optimal integer solution, which is how Desaulniers uses them. Remark 1 justifies pricing only elementary routes in column generation. Theorem 1 is the combinatorial fact underneath both corollaries and is reused throughout the split-delivery literature.

Formalizing it. The four results are proved in the literature (Dror–Trudeau; Gendreau et al. 2006); Desaulniers states them without proof. No machine-checked version exists. This mission produces a reusable Lean model of the SDVRPTW (instances, arc sets, feasible routes with schedules and delivery patterns, solutions, optimality) and checked proofs of these properties, including the existence of an optimal solution, which the paper takes for granted.

Difficulty

The obvious argument is local: take an optimal solution that violates the property, shift quantities between two routes, remove a visit, and shortcut. Three points make this less routine than it sounds. First, the statements are existential: each exchange must keep the solution optimal and not reintroduce a violation already removed, so one needs a termination measure that decreases under every exchange (Corollary 2 needs elementarity and the Theorem 1 property simultaneously, not two separate optimal solutions). Second, removing a visit is feasible only because the arc set is defined by time windows: the shortcut arc $(v, w)$ must be shown to exist from the schedule and the triangle inequality on travel times, and the new schedule must be built explicitly. Third, the feasible set is infinite (real quantities, unbounded walks, unbounded number of vehicles), so the existence of an optimal solution is itself a statement to prove, not a hypothesis.

Formalization scope

Namespace SplitDeliveryVRPTW.Known. Nodes are the inductive type Node n (start, cust i for i : Fin n, finish). All quantities, times and costs are real. The arc set is exactly the set the paper defines (read as "if and only if"), with arcs into the start depot and out of the end depot excluded. The triangle inequality for $t$ is imposed on pairwise distinct nodes and for $c$ on triples of arcs, where the paper's data is defined. A route is a customer list (repetitions allowed) with a time function over path positions and a quantity function over visits. A solution is an indexed family Fin m → Route I, so identical routes may appear twice. Demand satisfaction uses $\ge$, as constraint (2) does. The per-vehicle bound $\min\{d_i, Q\}$ of constraint (14) is omitted because it changes neither the feasible route patterns nor the costs.

Explicit readings of imprecise phrases:

• "split customer", "in common" (Theorem 1) are undefined in the paper. A customer two distinct routes visit is split, and the routes have it in common. This visit-based reading is at least as strong as a delivery-based one.
• Corollary 2's wording "at most one arc in set $\mathcal A^*_{ij}$ appears at most once" is read through constraint (7): the total number of traversals of the arcs of $\mathcal A^*_{ij}$ is at most one. It is stated in the equivalent form free of the choice of $\mathcal A^*(\mathcal N)$.
• "there exists an optimal solution" is conditional on feasibility, which is the hypothesis added.

Ruled-out trivializations: routes are not restricted to elementary walks (that would make Remark 1 definitional), solutions are not sets (that would forbid duplicate routes), optimality compares against solutions with any number of routes and any visit pattern, "split" is never counted over visits by the same route, capacity and time windows are part of route feasibility, and deliveries occur only at visits.

Useful infrastructure, reusable beyond this mission: shortcut lemmas for feasible routes (removing a visit), exchange lemmas between two routes, and existence of an optimum for split-delivery routing. Contributions of any of these as separate lemmas are welcome.

Selected references

• G. Desaulniers, Branch-and-Price-and-Cut for the Split-Delivery Vehicle Routing Problem with Time Windows, Operations Research 58(1):179–192, 2010. https://doi.org/10.1287/opre.1090.0713
• M. Dror, P. Trudeau, Savings by split delivery routing, Transportation Science 23(2):141–145, 1989. https://doi.org/10.1287/trsc.23.2.141
• M. Dror, P. Trudeau, Split delivery routing, Naval Research Logistics 37(3):383–402, 1990. https://doi.org/10.1002/nav.3800370304
• M. Gendreau, P. Dejax, D. Feillet, C. Gueguen, Vehicle routing with time windows and split deliveries, Technical Report 2006-851, Laboratoire Informatique d'Avignon, 2006.
• C. Archetti, M. W. P. Savelsbergh, M. G. Speranza, Worst-case analysis for split delivery vehicle routing problems, Transportation Science 40(2):226–234, 2006. https://doi.org/10.1287/trsc.1050.0117

Motivation

The linear-quadratic regulator (LQR) is the basic problem of optimal control: steer a linear system $\dot x = Ax + Bu$ so as to minimize an integrated quadratic cost of state and input. When the full state is measured, the optimal law is a static linear feedback $u = -Kx$ obtained from an algebraic Riccati equation. When only an output $y = Cx$ is measured, the optimal static output feedback $u = -Ky$ has no closed form, and the problem is nonconvex. In both cases the cost can be written as a function $f(K)$ of the gain alone, which makes it natural to minimize $f$ by gradient methods directly in the space of gains. This view goes back to Levine and Athans (1970) for output feedback and was revived in reinforcement learning by Fazel, Ge, Kakade and Mesbahi (2018), who proved global convergence of policy gradient for discrete-time state-feedback LQR.

Fatkhullin and Polyak (arXiv:2004.09875v2, SIAM J. Control Optim. 2021) carry out this program for the continuous-time problem. They show that $f$ is coercive, that it is $L$-smooth on every sublevel set, and that under state feedback it satisfies a gradient-domination (Łojasiewicz–Polyak) inequality. From these properties they derive convergence of the continuous gradient flow and of the discrete gradient method. This mission formalizes the continuous part, Theorem 4.1.

Setting

Fix real matrices $A\in\mathbb R^{n\times n}$, $B\in\mathbb R^{n\times m}$, $C\in\mathbb R^{r\times n}$ and symmetric positive definite weights $Q\in\mathbb R^{n\times n}$, $R\in\mathbb R^{m\times m}$ and covariance $\Sigma\in\mathbb R^{n\times n}$. A gain is $K\in\mathbb R^{m\times r}$, and the closed-loop matrix is $A_K = A - BKC$. A square matrix is Hurwitz if all its eigenvalues have negative real part. The stabilizing set is

$$\mathcal S = \{K : A_K \text{ is Hurwitz}\}.$$

For $K\in\mathcal S$, let $X(K)$ and $Y(K)$ be the unique solutions of the Lyapunov equations

$$A_K^\top X + XA_K + C^\top K^\top RKC + Q = 0,\qquad A_KY + YA_K^\top + \Sigma = 0 .$$

The LQR cost is $f(K) = \mathrm{Tr}(X(K)\Sigma)$. It equals the expected infinite-horizon quadratic cost when the initial state has covariance $\Sigma$. Its gradient with respect to the Frobenius inner product $\langle M,N\rangle = \mathrm{Tr}(M^\top N)$ is

$$\nabla f(K) = 2\,(RKC - B^\top X(K))\,Y(K)\,C^\top .$$

Given a known stabilizing gain $K_0\in\mathcal S$, the sublevel set is $\mathcal S_0 = \{K\in\mathcal S : f(K)\le f(K_0)\}$. The gradient flow is the ODE

$$\dot K(t) = -\nabla f(K(t)),\qquad K(0) = K_0 .$$

State feedback (SLQR) is the case $C = I$; there $f$ is written $f_S$, and $K_*$ denotes a minimizer of $f_S$ on $\mathcal S$. Norms: $\|\cdot\|_F$ is the Frobenius norm and $\|\cdot\|$ the spectral norm. $\lambda_1(M)$ is the smallest eigenvalue of a symmetric $M$.

Formalization targets

Goal: Theorem 4.1 for state feedback

For $C=I$, let $L>0$ be a Lipschitz constant of $\nabla f$ on $\mathcal S_0$ and let $\mu$ be the constant (3.11),

$$\mu = \frac{\lambda_1(R)\lambda_1^2(\Sigma)\lambda_1(Q)}{8f_S(K_*)\big(\|A\| + \|B\|^2 f_S(K_0)/(\lambda_1(\Sigma)\lambda_1(R))\big)^2}.$$

Then the flow has a solution on $[0,\infty)$, and every solution stays in $\mathcal S_0$, has $f(K_t)$ nonincreasing, has $\nabla f(K_t)\to0$ and $\min_{0\le t\le T}\|\nabla f(K_t)\|_F^2\le f(K_0)/T$, and

$$\|K_t - K_*\|_F \le \frac{\sqrt{2L(f(K_0)-f(K_*))}}{\mu}\,e^{-\mu t}\qquad(t\ge0).$$

Milestones

In attack order:

1. Existence of a minimizer (Corollary 3.10).
2. The gradient formula (Lemma 3.11).
3. Existence of a Lipschitz constant of $\nabla f$ on $\mathcal S_0$ (Theorem 3.15, qualitative).
4. The gradient-domination inequality $\tfrac12\|\nabla f_S(K)\|_F^2\ge\mu(f_S(K)-f_S(K_*))$ on $\mathcal S_0$ (Theorem 3.17).
5. The energy identity $\tfrac{d}{dt}f(K_t) = -\|\nabla f(K_t)\|_F^2$.
6. The integral bound of Appendix D.1.
7. The output-feedback part of Theorem 4.1 (general $C$ with $\operatorname{rank} C = r$): existence, invariance of $\mathcal S_0$, monotonicity and (4.2).

Significance

The theorem certifies a simple, model-based procedure: start from any stabilizing gain and follow the negative gradient of the cost. The iterate never loses stability, and it converges exponentially to the optimal gain. The same happens for state feedback even though neither $f$ nor $\mathcal S$ is convex. For output feedback it guarantees convergence to stationarity with an explicit $O(1/T)$ rate. The continuous-time result is the template for the discrete gradient method of the same paper (Theorem 4.2) and for policy-gradient analyses of LQR in learning-based control.

The result is proved in the paper. The paper writes out the proof of (4.2) in Appendix D.1. It describes the whole proof as a replica of Theorems 8 and 9 of Polyak 1963, which treat functions satisfying the Łojasiewicz–Polyak inequality on the whole space. No machine-checked version of any of these statements is known. A formal proof would give the first verified treatment of LQR as an optimization problem over gains. It needs Lyapunov equations, the stabilizing set and the cost, together with the analytic facts on which the policy-gradient literature rests: coercivity, smoothness on sublevel sets, and gradient domination.

Difficulty

The obvious argument is the textbook convergence proof for gradient flows of smooth functions satisfying the Łojasiewicz–Polyak inequality. That proof assumes the function is defined and smooth on the whole space. Here $f$ lives only on the open set $\mathcal S$, blows up at its boundary, is unbounded on $\mathcal S$ and is not globally $L$-smooth. The Łojasiewicz–Polyak inequality holds only on $\mathcal S_0$, with a constant that depends on $K_0$. The work therefore lies in keeping the trajectory inside $\mathcal S_0$ and showing it exists for all time. This requires compactness of $\mathcal S_0$ (coercivity of $f$) and a positive distance from $\mathcal S_0$ to the boundary of $\mathcal S$. Local ODE existence alone gives neither. The explicit constant $\mu$ additionally requires lower bounds on solutions of Lyapunov equations and an upper bound on $\|K\|_F$ over $\mathcal S_0$.

Formalization scope

Matrices are Matrix (Fin p) (Fin q) ℝ. Hurwitz means every complex eigenvalue of the complexified matrix has negative real part. $X(K)$ and $Y(K)$ are the unique solutions of their Lyapunov equations. Where the solution is not unique (possible only off $\mathcal S$) they take the value $0$, and no statement uses these off-$\mathcal S$ values. $\nabla f$ is defined by formula (3.3), and Lemma 3.11 states that it is the Fréchet derivative with respect to the Frobenius inner product. A solution of the flow is a curve with $K(0)=K_0$ that stays in $\mathcal S$ for $t\ge0$. Its entries have the prescribed derivative within $[0,\infty)$. Global existence and invariance of $\mathcal S_0$ are conclusions, never hypotheses. Assuming a global solution, or a solution confined to $\mathcal S_0$, would trivialize the theorem, and that encoding is excluded. The case $C=I$ is a matrix argument $C$ with $C=1$. $\lambda_1$ is the minimum of the (unsorted) eigenvalues, and the spectral norm is $\sqrt{\lambda_{\max}(M^\top M)}$. "Monotone decreasing" is rendered as nonincreasing. $\min_{0\le t\le T}$ is rendered as the existence of a point of $[0,T]$ where the bound holds.

The standing assumptions of p. 3 are hypotheses wherever they apply: $K_0\in\mathcal S$, $Q,R,\Sigma\succ0$, $\operatorname{rank}C=r$, and $B\ne0$. Deviations from the printed paper:

• The paper's explicit Lipschitz constant (3.8) is false as printed; the planning notes record counterexamples for both output and state feedback. The goal therefore takes $L$ as any Lipschitz constant of $\nabla f$ on $\mathcal S_0$, which is the paper's own definition of $L$-smoothness (§3.6). Theorem 3.15 enters only in its qualitative form, as the existence of such an $L$.
• The set-builder for $\mathcal S$ on p. 3 writes $\mathbb R^{m\times n}$; gains are $m\times r$.
• The optimal gain $K_*$ is a hypothesis (its existence is Corollary 3.10), and $f(K_*)>0$ is not assumed.

A complete development needs the following: Lyapunov equations and their solution theory (uniqueness, positivity, integral representation), continuity and differentiability of $K\mapsto X(K)$ on $\mathcal S$, coercivity of $f$, global existence for ODEs confined to a compact invariant set, and the Łojasiewicz–Polyak argument for flows on a subset. The Lyapunov and stabilizing-set layer is reusable across control missions. Proofs of any milestone, and of auxiliary lemmas such as Lemmas 3.8, A.5, C.1–C.3, are welcome.

Selected references

• I. Fatkhullin, B. Polyak, Optimizing Static Linear Feedback: Gradient Method, SIAM J. Control Optim. 59 (2021); preprint arXiv:2004.09875v2. https://arxiv.org/abs/2004.09875
• M. Fazel, R. Ge, S. Kakade, M. Mesbahi, Global Convergence of Policy Gradient Methods for the Linear Quadratic Regulator, ICML 2018. https://arxiv.org/abs/1801.05039
• W. S. Levine, M. Athans, On the determination of the optimal constant output feedback gains for linear multivariable systems, IEEE Trans. Automat. Control 15 (1970). https://doi.org/10.1109/TAC.1970.1099363
• H. Karimi, J. Nutini, M. Schmidt, Linear Convergence of Gradient and Proximal-Gradient Methods Under the Polyak–Łojasiewicz Condition, ECML PKDD 2016. https://arxiv.org/abs/1608.04636
• B. T. Polyak, Gradient methods for the minimisation of functionals, USSR Comput. Math. Math. Phys. 3 (1963), 864–878. https://doi.org/10.1016/0041-5553(63)90382-3