Operations Research

Motivation

Multistage stochastic programs are solved on a finite scenario tree: a discrete probability distribution whose support points are paths of a random process. Realistic trees have far too many scenarios for the resulting optimization problem, so practitioners reduce the tree, keeping $n$ of its $N$ scenarios and redistributing the probability of the deleted ones. The reduction should keep the reduced distribution as close as possible to the original one in a probability metric that controls the optimal value of the stochastic program (Dupačová, Gröwe-Kuska, Römisch, Math. Program. 95 (2003)).

Choosing the best $n$ scenarios is a set-covering problem and NP-hard, and the algorithms used in practice (backward reduction, fast forward selection) are heuristics without error guarantees. Heitsch and Römisch (2003) therefore derived test instances with an exactly known optimum: regular binary and ternary scenario trees, for which the minimal reduction distance has a closed form once $n$ is not too small. This mission formalizes the ternary case, Proposition 3.2 of that paper. The binary case (Proposition 3.1) is a separate mission of the same series.

Setting

Fix a depth $K \in \mathbb{N}$ and branch widths $\delta^1, \dots, \delta^K \ge 0$, with $\delta^0 = 0$. A regular ternary scenario tree has $N = 3^K$ scenarios, one for each index tuple $(i_1, \dots, i_K) \in \{1, 2, 3\}^K$, where $i_k$ is the successor chosen at level $k$. Choosing successor $i_k$ adds the increment $\delta^k_{i_k} = (i_k - 2)\,\delta^k \in \{-\delta^k, 0, \delta^k\}$, and scenario $i$ is the vector $\omega_i = (\omega_i^0, \dots, \omega_i^K) \in \mathbb{R}^{K+1}$ with

$$\omega_i^k = \sum_{j=0}^{k} \delta^j_{i_j}, \qquad k = 0, \dots, K \quad \text{(eq. (19))}.$$

All scenarios have probability $p_i = 1/N$. The distance between scenarios is the maximum norm $c(\omega_i, \omega_j) = \|\omega_i - \omega_j\|_\infty = \max_{0 \le k \le K} |\omega_i^k - \omega_j^k|$.

Deleting the scenarios of an index set $J$ and moving each deleted scenario's probability to a nearest kept scenario costs the reduction distance

$$D_J = \sum_{i \in J} p_i \min_{j \notin J} \|\omega_i - \omega_j\|_\infty \quad \text{(eq. (8))},$$

which by Theorem 2.1 of the paper is the optimal transport-type distance between the original distribution and the best distribution supported on the kept scenarios. The minimal reduction distance to $n$ scenarios is $D^{min}_n = \min\{D_J : \#J = N - n\}$.

Formalization targets

Goal: Proposition 3.2 (7/9-solution)

Let $K \ge 3$ and let $k_0 \in \arg\min_{1 \le k \le K} \delta^k$ with $k_0 \le K - 2$ and $\max\{\delta^{k_0+1}, \delta^{k_0+2}\} \le 2\delta^{k_0}$. Then any two distinct scenarios are at distance at least $\delta^{k_0}$; there is a set of $\tfrac79 N$ scenarios each paired with a scenario outside it at distance exactly $\delta^{k_0}$; and for each $n \in \mathbb{N}$ with $\tfrac29 N \le n < N$,

$$D^{min}_n = \min\{D_J : \#J = N - n\} = \frac{N - n}{N}\,\delta^{k_0} \quad \text{(eq. (21))},$$

with the minimum attained.

Milestones

1. Distinct scenarios satisfy $\|\omega_i - \omega_j\|_\infty \ge \delta^{k_0}$.
2. Every $J$ with $\#J = N - n$ has $D_J \ge \frac{N-n}{N}\delta^{k_0}$.
3. The index set $I_{**}$ of the proof has $\#I_{**} = \tfrac29 N$, and its complement $J_{**}$ has $\tfrac79 N$ elements.
4. Every $j \in J_{**}$ has a partner $i \in I_{**}$ with $\|\omega_i - \omega_j\|_\infty = \delta^{k_0}$.
5. Example 4.2: for $K = 6$ and $(\delta^1, \dots, \delta^6) = (0.7, 0.9, 1.2, 1.5, 2.6, 3.3)$, $D^{min}_n = 0.7\,\frac{N-n}{N}$ for $162 \le n < 729$.

Significance

The result gives an exact optimal value for an NP-hard reduction problem on an infinite family of instances. Heitsch and Römisch use it in their numerical section to measure how far the heuristics' reduced trees are from optimal (Examples 4.1 and 4.2 are the binary and ternary test trees of that study). A closed form of this kind is also the only way to certify that a heuristic is exactly optimal on some instances rather than only competitive with other heuristics.

The proposition is proved in the paper, but the published proof of its central step is one sentence: "Similarly as in Proposition 3.1 it can be shown that there exists an index $i \in I_{**}$ for each $j \in J_{**}$ …". A formal proof supplies that case analysis, which is absent from the literature. The formalization also settles two points the printed statement leaves loose (see Formalization scope): the count of pairs at distance $\delta^{k_0}$, and the definition of $I_{**}$ when some widths vanish. No machine-checked version of this result or of the reduction distance $D_J$ is known to exist.

Difficulty

The lower bound is routine: two distinct scenarios first differ at some level $l$, where their coordinates differ by $\delta^l$ or $2\delta^l$. The substance is attainment: one must exhibit, for every $n \ge \tfrac29 N$, a kept set of size $n$ whose every deleted scenario lies at distance exactly $\delta^{k_0}$ from some kept one. The obvious candidate, keeping the scenarios that take the middle branch at level $k_0$, has every other scenario at distance exactly $\delta^{k_0}$ from a kept one, but it keeps $\tfrac13 N$ scenarios and so covers only $n \ge \tfrac13 N$. Going down to $\tfrac29 N$ kept scenarios forces a deleted scenario and its partner to differ at more than one level, and since coordinates are running sums the differences at levels $k_0+1$ and $k_0+2$ accumulate on top of the one at level $k_0$. The paper's proof of this step is not written out, and it depends on the widths of the two levels below $k_0$: the hypothesis $\max\{\delta^{k_0+1}, \delta^{k_0+2}\} \le 2\delta^{k_0}$ is essential, and the result is false without it (for $K = 3$ and $(\delta^1, \delta^2, \delta^3) = (1, 3, 3)$ one has $D^{min}_6 = 31/27$, not $7/9$).

Formalization scope

• A scenario is an index tuple $\sigma : \mathrm{Fin}\,K \to \mathrm{Fin}\,3$; $\sigma(r)$ is the successor at paper level $r + 1$, with $\mathrm{Fin}\,3$ values $0, 1, 2$ standing for the paper's $i = 1, 2, 3$. The widths are $\delta : \mathbb{N} \to \mathbb{R}$, of which only $\delta(1), \dots, \delta(K)$ are used; the standing assumption $\delta^k \in \mathbb{R}_+$ (p. 196) is the hypothesis $\delta(k) \ge 0$ for $1 \le k \le K$. Scenarios live in $\mathrm{Fin}(K+1) \to \mathbb{R}$, whose Mathlib norm is the maximum norm. Probabilities are uniform, $1/3^K$.
• $D_J$ is defined for a general finite index set, probabilities and cost, with the inner minimum a Finset.inf' over the complement of $J$; the complement must be nonempty, so no default value arises. $D^{min}_n$ is stated as IsLeast of the set of all values $D_J$ with $\#J = N - n$: the goal asserts both the lower bound for every $J$ and attainment by some $J$. A formalization that exhibits a single $J$ with $D_J = \frac{N-n}{N}\delta^{k_0}$, or states an infimum without attainment, or drops any hypothesis on $k_0$, is a different (and in the last case false) statement.
• $\tfrac29 N \le n$ is written $2 \cdot 3^K \le 9n$, and $\tfrac79 N$ as $7 \cdot 3^{K-2}$.
• Pairs. As printed, "there are $\tfrac79 N$ distinct pairs of scenarios such that the distance between the members of each pair is exactly $\delta^{k_0}$" is false as an exact count (for $K = 3$ and $\delta = (1,1,1)$ there are 130 such pairs, not 21). The goal states what the proof constructs: a set $J_{**}$ of $\tfrac79 N$ scenarios, each paired with a scenario outside $J_{**}$ at distance exactly $\delta^{k_0}$.
• $I_{**}$. The paper defines $I_{**}$ by testing whether the increments $\delta^k_{i_k}$ vanish. When one of $\delta^{k_0}, \delta^{k_0+1}, \delta^{k_0+2}$ is $0$ this no longer identifies the middle branch and the count $\tfrac29 N$ fails, although the proposition remains true. The formalization defines $I_{**}$ by branch indices (middle branch versus outer branches), which agrees with the paper whenever these three widths are positive. No positivity hypothesis is added to the goal.
• Welcome contributions: a reusable library for regular scenario trees (first differing level, distance of paths), and the case analysis of milestone 4. The binary mission of this series needs the same lower-bound argument with the constant $2\delta^{k_0}$.

Selected references

• H. Heitsch, W. Römisch, Scenario Reduction Algorithms in Stochastic Programming, Computational Optimization and Applications 24 (2003), 187–206. https://doi.org/10.1023/A:1021805924152
• J. Dupačová, N. Gröwe-Kuska, W. Römisch, Scenario reduction in stochastic programming: An approach using probability metrics, Mathematical Programming 95 (2003), 493–511. https://doi.org/10.1007/s10107-002-0331-0

Why exact reduction distances matter

Multistage stochastic programs are solved on a finite scenario tree, a discrete probability measure whose atoms are paths of a stochastic process. The size of the deterministic equivalent grows with the number of scenarios, so practitioners replace the original measure $P=\sum_{i=1}^N p_i\delta_{\omega_i}$ by a measure supported on $n<N$ of its scenarios. Stability theory for stochastic programs (Dupačová, Gröwe-Kuska and Römisch, Math. Program. 95 (2003), doi:10.1007/s10107-002-0331-0) bounds the change of the optimal value by a probability metric between the two measures, which leads to the optimal scenario reduction problem: choose which $N-n$ scenarios to delete so that this distance is smallest.

That problem is a set-covering problem and is NP-hard, and the algorithms that Heitsch and Römisch study in the same paper (backward reduction, fast forward selection) are heuristics without error guarantees. To test them one needs original measures whose optimal reduction distance is known exactly. Section 3 of Heitsch and Römisch, Scenario Reduction Algorithms in Stochastic Programming, Comput. Optim. Appl. 24 (2003) (doi:10.1023/A:1021805924152) supplies such instances: regular binary and ternary scenario trees, for which the minimal distance to any reduced tree with at least a fixed fraction of the scenarios is an explicit formula. This mission formalizes the binary case, Proposition 3.1.

Setting

Fix a horizon $K\in\mathbb N$ and level parameters $\delta^1,\dots,\delta^K\ge0$, with $\delta^0=0$. A regular binary scenario tree has $N=2^K$ scenarios. A scenario is determined by a branch $i_k\in\{1,2\}$ at every level $k=1,\dots,K$, and is the vector $\omega_i=(\omega_i^0,\dots,\omega_i^K)\in\mathbb R^{K+1}$ with

$$\omega_i^k=\sum_{j=0}^k\delta^j_{i_j},\qquad \delta^j_{i_j}=(2i_j-3)\,\delta^j\in\{-\delta^j,+\delta^j\}.\tag{19}$$

All scenarios start at the root $\omega_i^0=0$ and carry probability $p_i=1/N$. Scenarios are compared in the maximum norm $\|\omega-\tilde\omega\|_\infty=\max_{k=0,\dots,K}|\omega^k-\tilde\omega^k|$.

Deleting the scenarios with indices in $J\subset\{1,\dots,N\}$ and moving each deleted scenario's probability to a nearest kept scenario costs the reduction cost

$$D_J=\sum_{i\in J}p_i\min_{j\notin J}\|\omega_i-\omega_j\|_\infty,\tag{8}$$

which by Theorem 2.1 of the paper is the minimal Kantorovich-type distance between $P$ and a measure supported on the kept scenarios. The minimal reduction distance for $n$ kept scenarios is $D^{min}_n=\min\{D_J:\#J=N-n\}$.

In the Lean development, scenarios are indexed by σ : Fin K → Fin 2 (Fin-index $r$ is tree level $r+1$, value $0$ is the branch $-\delta$, value $1$ is $+\delta$), lev σ k is the branch at level $k$, scenario δ σ : Fin (K+1) → ℝ is $\omega_\sigma$, and redCost δ J hJ is $D_J$.

Formalization targets

Goal: Proposition 3.1 (3/4-solution)

Let $K\ge3$, $k_0\in\arg\min_{1\le k\le K}\delta^k$, $k_0\le K-2$ and $\max\{\delta^{k_0+1},\delta^{k_0+2}\}\le2\delta^{k_0}$. Then any two distinct scenarios are at distance at least $2\delta^{k_0}$; there is a set $J_*$ of $\frac34N$ scenarios each of which has a partner outside $J_*$ at distance exactly $2\delta^{k_0}$; and for every $n\in\mathbb N$ with $\frac N4\le n<N$

$$D^{min}_n=\min\{D_J:\#J=N-n\}=\frac{N-n}{N}\,2\delta^{k_0}.\tag{20}$$

Milestones

1. Two scenarios that first differ at level $l$ are at distance $\ge2\delta^l\ge2\delta^{k_0}$.
2. $D_J\ge\frac{N-n}{N}2\delta^{k_0}$ for every $J$ with $\#J=N-n$.
3. The index set $I_*$ (branch at $k_0$ opposite to the common branch at $k_0+1,k_0+2$) has $\#I_*=N/4$, so $\#J_*=\frac34N$.
4. Every $j\in J_*$ has a partner $i\in I_*$ with $\|\omega_i-\omega_j\|_\infty=2\delta^{k_0}$.
5. Example 4.1: for $K=10$ and the paper's parameters, Proposition 3.1 applies with $k_0=1$ and $D^{min}_n=\frac{N-n}{N}$ for $256\le n<1024$.

Significance

The result. Proposition 3.1 gives the exact optimum of an NP-hard combinatorial problem on an explicit, parametrized family of instances of every size $N=2^K$. Section 4 of the paper uses it (Example 4.1, $N=1024$) as ground truth for the relative accuracy of backward reduction and fast forward selection. Without it, the quality of a heuristic reduction on a large tree could only be compared with other heuristics or with lower bounds.

The formalization. The proposition is proved in the paper; to our knowledge no machine-checked version exists. A formal proof certifies the benchmark values, and the statement also corrects the printed text in two places. First, "there are $\frac34N$ distinct pairs … at distance exactly $2\delta^{k_0}$" is false as an exact count (for $K=3$, $\delta=(1,1,1)$ there are twenty such pairs, not six), so the goal states "at least", in the form the proof exhibits. Second, the paper's sign-based definition of $I_*$ degenerates when some $\delta^k=0$, although the proposition still holds; the mission defines $I_*$ by branch indices.

Difficulty

The lower bound is a direct computation. The content is the matching upper bound: a set $J$ of the prescribed size for which every deleted scenario has a kept scenario at the minimal possible distance. A natural first attempt pairs scenarios that differ only at level $k_0$. That handles only half of the scenarios with a single partner each, and it cannot reach $\frac34N$ deleted scenarios. Once scenarios differ at more than one level, their maximum-norm distance is a maximum of several partial sums, and keeping all of them at most $2\delta^{k_0}$ is exactly where the hypothesis $\max\{\delta^{k_0+1},\delta^{k_0+2}\}\le2\delta^{k_0}$ enters. Without it, eq. (20) fails: for $K=3$, $\delta=(1,3,3)$ and $n=2$ the true minimum is $\frac52$, not $\frac32$. The passage from $n=N/4$ to general $n\ge N/4$ also needs care: the deleted set must shrink while each remaining deleted scenario keeps its partner among the kept ones.

Formalization scope

• The index type is Fin K → Fin 2, which has exactly $2^K$ elements. The paper's $(K+1)$-tuple has a level-0 entry with no choice, so it is dropped, and the vector $\omega$ keeps its $K+1$ coordinates with $\omega^0=0$.
• The parameters are δ : ℕ → ℝ, with $\delta^k\ge0$ and the arg min stated for $k=1,\dots,K$ only. $\delta^k>0$ is not assumed, because the paper allows $\delta^k\in\mathbb R_+$ and the proposition holds with zeros.
• The cost is Mathlib's norm on Fin (K+1) → ℝ, which is the maximum norm, and $p_i=1/2^K$ is written out.
• $D_J$ requires a nonempty set of kept scenarios, and its inner minimum is a finite Finset.inf'.
• $D^{min}_n$ is stated with IsLeast over the set of attained values $D_J$, $\#J=N-n$, so both the lower bound and attainment are part of the goal. A proof of $D_{J}\le\frac{N-n}{N}2\delta^{k_0}$ for a single exhibited $J$, or a real infimum without attainment, does not prove the goal.
• "$\frac N4\le n$" is written $2^K\le4n$, and $\frac34N$ is $3\cdot2^{K-2}$.
• The pairs claim is "at least $\frac34N$ pairs", expressed as a set $J_*$ of that size with a partner outside $J_*$ for every member.
• All hypotheses on $k_0$ appear in the goal; dropping any of them makes (20) false.

Useful infrastructure: sup-norm lemmas for Fin n → ℝ (pi_norm_le_iff_of_nonneg, norm_le_pi_norm), Finset.inf' lemmas, and counting functions Fin K → Fin 2 with prescribed values (Fintype.card_fun, Fintype.card_pi). The tree and reduction-cost definitions are shared in spirit with the ternary-tree mission of this series (Proposition 3.2), and a proof whose structure transfers to $d=3$ is welcome. Contributions of proofs of the milestones individually, in any order, are welcome.

Selected references

• H. Heitsch, W. Römisch, Scenario Reduction Algorithms in Stochastic Programming, Computational Optimization and Applications 24 (2003), 187–206. doi:10.1023/A:1021805924152
• J. Dupačová, N. Gröwe-Kuska, W. Römisch, Scenario reduction in stochastic programming: An approach using probability metrics, Mathematical Programming 95 (2003), 493–511. doi:10.1007/s10107-002-0331-0

Why reduce scenarios

Multistage and two-stage stochastic programs are solved numerically on a discrete probability distribution: a finite set of scenarios $\omega_1,\dots,\omega_N$ with probabilities $p_1,\dots,p_N$. The size of the resulting optimization problem grows with $N$, and scenario sets produced by sampling or by historical data are often far too large to be solved directly. Scenario reduction replaces the original distribution by one supported on a small subset of the scenarios, chosen so that the optimal value and solutions of the stochastic program change as little as possible.

Stability theory for stochastic programs (Rachev and Römisch, 2002) shows that this change is controlled by a probability distance of Fortet–Mourier type, which for discrete measures is bounded by the value of a transportation problem. Dupačová, Gröwe-Kuska and Römisch (2003) turned this into a combinatorial problem and proposed greedy backward and forward heuristics. Heitsch and Römisch (2003) gave faster versions of both heuristics; the forward version, fast forward selection, is the subject of this mission. Implementations of these reduction heuristics are distributed with the GAMS modelling system (SCENRED) and are used in energy and finance applications of stochastic programming.

Setting

Let $E$ be a finite-dimensional real vector space with a norm $\|\cdot\|$, let $\omega_0\in E$, and let $h:[0,\infty)\to[0,\infty)$ be continuous and nondecreasing with $h(0)=0$. The cost between two points of $E$ is

$$c(\omega,\tilde\omega)=\max\bigl\{1,\,h(\|\omega-\omega_0\|),\,h(\|\tilde\omega-\omega_0\|)\bigr\}\,\|\omega-\tilde\omega\| .$$

It is nonnegative, symmetric, and zero on the diagonal.

The original distribution is $P=\sum_{i=1}^N p_i\delta_{\omega_i}$ with $p_i>0$ and $\sum_i p_i=1$. Deleting the scenarios in a set $J\subset\{1,\dots,N\}$ and assigning new weights $q_j\ge 0$, $\sum_{j\notin J}q_j=1$, to the kept ones gives $Q=\sum_{j\notin J}q_j\delta_{\omega_j}$. The distance $D(J;q)$ between $P$ and $Q$ is the optimal value of the transportation problem

$$D(J;q)=\min\Bigl\{\sum_{i=1}^N\sum_{j\notin J}c(\omega_i,\omega_j)\eta_{ij}:\ \eta_{ij}\ge 0,\ \sum_{i}\eta_{ij}=q_j,\ \sum_{j\notin J}\eta_{ij}=p_i\Bigr\}.$$

The reduction cost of deleting $J$ is

$$D_J=\sum_{i\in J}p_i\min_{j\notin J}c(\omega_i,\omega_j),$$

and the optimal reduction problem (8) minimizes $D_J$ over all $J$ with $\#J=N-n$, where $n$ is the number of scenarios to keep.

Forward selection builds the kept set greedily. With $J^{[0]}=\{1,\dots,N\}$ and $J^{[i]}=\{1,\dots,N\}\setminus\{u_1,\dots,u_i\}$, it chooses

$$u_i\in\arg\min_{u\in J^{[i-1]}}D_{J^{[i-1]}\setminus\{u\}},\qquad i=1,\dots,n. \tag{16}$$

Fast forward selection (Algorithm 2.4) computes the same choices through an updated cost matrix: $c^{[1]}_{ku}=c(\omega_k,\omega_u)$, $c^{[i]}_{ku}=\min\{c^{[i-1]}_{ku},c^{[i-1]}_{ku_{i-1}}\}$, $z^{[i]}_u=\sum_{k\in J^{[i-1]}\setminus\{u\}}p_kc^{[i]}_{ku}$, and $u_i\in\arg\min_{u\in J^{[i-1]}}z^{[i]}_u$.

Formalization targets

Goal: Theorem 2.5

For $1\le n\le N$ and every run $u_1,\dots,u_n$ of Algorithm 2.4, with any tie-breaking in the arg min,

$$u_i\ \text{satisfies (16)}\quad\text{and}\quad z^{[i]}_{u_i}=D_{J^{[i]}}\qquad(i=1,\dots,n).$$

Milestones

1. Theorem 2.1 (redistribution). For $J$ with at least one kept scenario, $D_J=\min_q D(J;q)$, and the minimum is attained at $\bar q_j=p_j+\sum_{i\in J,\,j(i)=j}p_i$ for every choice of nearest kept scenarios $j(i)$.
2. Eq. (10). $D_{\{1,\dots,N\}\setminus\{u\}}=\sum_{i=1}^Np_ic(\omega_i,\omega_u)$, so (8) with $\#J=N-1$ is problem (10).
3. Eq. (12). The sum $lb$ of the $N-n$ smallest single-deletion costs $p_l\min_{j\neq l}c(\omega_l,\omega_j)$, taken in the greedy order (11), is at most $D_J$ for every $J$ with $\#J=N-n$.
4. Optimality condition (p. 191). If each $l_i$ has a nearest other scenario outside $\{l_1,\dots,l_{N-n}\}\setminus\{l_i\}$, then $\{l_1,\dots,l_{N-n}\}$ solves (8).
5. Eq. (17), unrolled recursion. For any index sequence, $c^{[i]}_{ku}=\min_{j\notin J^{[i-1]}\setminus\{u\}}c(\omega_k,\omega_j)$ for $u\in J^{[i-1]}$.
6. Eq. (17), conclusion. For any index sequence, $z^{[i]}_u=D_{J^{[i-1]}\setminus\{u\}}$ for $u\in J^{[i-1]}$.

Significance

Theorem 2.5 certifies that the cheap update of Algorithm 2.4 (one pairwise minimum per matrix entry and step) produces exactly the greedy forward selection defined through the reduction costs, and that the running objective $z^{[i]}_{u_i}$ is the reduction cost of the scenarios deleted so far. Combined with Theorem 2.1, $z^{[i]}_{u_i}$ is the optimal transportation distance between $P$ and the best measure on the kept scenarios, which is the quantity practitioners monitor to decide how many scenarios to keep. The lower bound (12) and the optimality condition give a posteriori quality certificates for any reduced set.

All results of this mission are proved in the paper or in the works it cites (Dupačová et al., 2003); none is open. To the best of a platform search, none has been machine-checked. The mission provides a verified specification of a widely deployed algorithm, a formal link between a combinatorial set-covering objective and a finite transportation problem, and definitions (reduction cost, transportation plans with a partially free target marginal, greedy runs with arbitrary tie-breaking) reusable by the regular-tree missions of this series and by later scenario-tree construction papers.

Difficulty

The mathematics is elementary; the difficulty is bookkeeping. The recursion for $c^{[i]}$ refers to the previous step's column $u_{i-1}$, which itself was updated, so unrolling it to a minimum over $\{u,u_1,\dots,u_{i-1}\}$ is an induction on $i$ in which the index sets $J^{[i]}$, the 1-based step counter and the complement structure all move together. The natural first attempt, identifying $c^{[i]}_{ku}$ with the minimum over the complement of $J^{[i]}$, is off by one step: the correct set is the complement of $J^{[i-1]}\setminus\{u\}$, which contains $u$ itself. For Theorem 2.1 the lower bound requires using that $c(\omega_i,\omega_i)=0$ for kept scenarios and that every plan ships all of $p_i$ somewhere outside $J$; the attainment part requires constructing the plan explicitly from the choice $j(\cdot)$, including scenarios for which several kept scenarios are equally near.

Formalization scope

• Scenarios are ω : Fin N → E with E a finite-dimensional real normed space; the paper's closed set $\Omega\subset\mathbb R^s$ plays no role beyond containing the scenarios and is omitted. Scenarios need not be distinct.
• $h$ is a function ℝ → ℝ with the paper's assumptions imposed on $[0,\infty)$ (IsGrowthFunction); every theorem carries them, together with $p_i>0$ and $\sum_ip_i=1$.
• The functions $f_0$, $g$ and the stochastic program (1)–(2) that motivate $c$ appear in no statement.
• $D(J;q)$ is the paper's finite transportation problem (p. 188), not the Kantorovich functional on measures. Weights $q$ and plans $\eta$ are indexed by all of $\{1,\dots,N\}$ with entries at deleted indices fixed to $0$.
• $D_J$ requires a proof that the complement of $J$ is nonempty; minima are Finset.inf', never a real infimum with a default value.
• Algorithm 2.4 is a relation on sequences u : ℕ → Fin N with 1-based steps. $c^{[i]}$ is the printed recursion, extended to all indices; runs are any sequences satisfying the arg-min conditions, so every tie-breaking rule is covered.
• The paper's standing restriction $n<N$ is relaxed to $n\le N$ in Theorem 2.5; the statement remains true at $n=N$.
• A trivializing formalization is ruled out: defining $c^{[i]}$ or $z^{[i]}$ directly as the minimum over the selected set or as $D_{J^{[i-1]}\setminus\{u\}}$ would make Theorem 2.5 hold by definition, and proving it for one fixed tie-breaking rule would prove less than the paper; neither is done.
• Proofs of the milestones, alternative proofs of Theorem 2.1 via LP duality, and a verified executable implementation of Algorithm 2.4 are all welcome.

Selected references

• H. Heitsch, W. Römisch, Scenario Reduction Algorithms in Stochastic Programming, Computational Optimization and Applications 24 (2003), 187–206. https://doi.org/10.1023/A:1021805924152
• J. Dupačová, N. Gröwe-Kuska, W. Römisch, Scenario reduction in stochastic programming: an approach using probability metrics, Mathematical Programming 95 (2003), 493–511. https://doi.org/10.1007/s10107-002-0331-0
• S. T. Rachev, W. Römisch, Quantitative stability in stochastic programming: the method of probability metrics, Mathematics of Operations Research 27 (2002), 792–818. https://doi.org/10.1287/moor.27.4.792.304

Motivation

A monotone variational inequality asks for a point at which a monotone operator and a convex function are in equilibrium. It unifies convex minimization (where $F$ is a gradient), convex–concave saddle-point problems (where $F$ is the skew gradient of a Lagrangian), Nash equilibria of monotone games, and complementarity problems in economics and traffic assignment. In operations research, first-order methods for such problems are the workhorse behind large-scale saddle-point formulations of linear and conic programs, where only one operator evaluation and one projection or proximal step per iteration are affordable.

The classical method for Lipschitz monotone operators is Korpelevich's extragradient method (1976) and its proximal variant, Tseng's forward–backward–forward method (2000); both need two evaluations of $F$ per iteration. The reflected projected gradient method of Malitsky (SIAM J. Optim., 2015) uses one evaluation of $F$ but evaluates it at $2z^k-z^{k-1}$, a point that may lie outside the domain of $g$. Malitsky's Golden Ratio Algorithm (GRAAL), introduced in Golden Ratio Algorithms for Variational Inequalities (preprint 2018; published in Mathematical Programming, doi:10.1007/s10107-019-01416-w), uses one evaluation of $F$, always at a feasible point, and one proximal step per iteration. Its fixed-step version, Theorem 1 of that paper, is the subject of this mission; the explicit, adaptive-step version (Theorem 2) is a separate mission of this series.

Setting

Let $\mathcal E$ be a finite-dimensional real inner product space with inner product $\langle\cdot,\cdot\rangle$ and norm $\|\cdot\| = \sqrt{\langle\cdot,\cdot\rangle}$. Let $g:\mathcal E\to(-\infty,+\infty]$ and write $\operatorname{dom} g = \{x : g(x)<+\infty\}$. Let $F:\operatorname{dom} g\to\mathcal E$. The variational inequality is

$$\text{find } z^*\in\mathcal E \quad\text{such that}\quad \langle F(z^*), z-z^*\rangle + g(z)-g(z^*)\ \ge\ 0\qquad \forall z\in\mathcal E. \tag{1}$$

The standing assumptions are:

• (C1) the solution set $S$ of (1) is nonempty;
• (C2) $g$ is proper (never $-\infty$, finite somewhere), convex, and lower semicontinuous;
• (C3) $F$ is monotone: $\langle F(u)-F(v),u-v\rangle\ge0$ for all $u,v\in\operatorname{dom} g$.

The proximal operator of $g$ is $\operatorname{prox}_g(z) = \operatorname{argmin}_x\{g(x)+\tfrac12\|x-z\|^2\}$. Let $\varphi = \frac{\sqrt5+1}{2}$ be the golden ratio, so that $\varphi^2 = 1+\varphi$. For a step $\lambda>0$ and arbitrary starting points $z^1,\bar z^0\in\mathcal E$, the Golden Ratio Algorithm generates, for $k\ge1$,

$$\bar z^k = \frac{(\varphi-1)z^k + \bar z^{k-1}}{\varphi},\qquad z^{k+1} = \operatorname{prox}_{\lambda g}\big(\bar z^k - \lambda F(z^k)\big). \tag{6}$$

The first line is a convex combination of the newest iterate and the previous average; the second is a forward–backward step taken from the average rather than from $z^k$.

Formalization targets

Goal: Theorem 1

If $F$ is $L$-Lipschitz on $\operatorname{dom} g$ ($L>0$), (C1)–(C3) hold, and $\lambda\in\big(0,\frac{\varphi}{2L}\big]$, then there is $z^*\in S$ with

$$z^k\to z^*\qquad\text{and}\qquad \bar z^k\to z^*\qquad(k\to\infty).$$

Both sequences converge, to one and the same solution. The goal is stated with the paper's exact step range; no rate is claimed, as the paper claims none.

Milestones

1. Eq. (4), the prox-inequality: for proper convex lsc $g$,

$$\bar x = \operatorname{prox}_g z \iff \langle\bar x - z, x-\bar x\rangle\ge g(\bar x)-g(x)\quad\forall x\in\mathcal E.$$

2. Eq. (12), an identity using only the averaging step of (6): for every point $z^*$,

$$\|z^{k+1}-z^*\|^2 = (1+\varphi)\|\bar z^{k+1}-z^*\|^2-\varphi\|\bar z^k-z^*\|^2+\tfrac1\varphi\|z^{k+1}-\bar z^k\|^2 .$$

3. Eq. (14), the energy inequality: for $z^*\in S$ and $k\ge2$,

$$(1+\varphi)\|\bar z^{k+1}-z^*\|^2+\tfrac\varphi2\|z^{k+1}-z^k\|^2\le(1+\varphi)\|\bar z^k-z^*\|^2+\tfrac\varphi2\|z^k-z^{k-1}\|^2-\varphi\|z^k-\bar z^k\|^2 .$$

4. Lemma 1 (Bauschke–Combettes, Theorem 5.5): a sequence that is Fejér monotone with respect to a nonempty set $C$ and whose cluster points all lie in $C$ converges to a point of $C$.

Significance

The result. Theorem 1 shows that monotone variational inequalities with a Lipschitz operator can be solved with one operator evaluation and one proximal step per iteration, with $F$ evaluated only at points of $\operatorname{dom} g$, where it is defined. This matters when $F$ is expensive (a large matrix–vector product, a simulation) or undefined outside the feasible set (for instance an operator involving $\log x$ on the positive orthant). The analysis also explains the constant: the averaging weight $\varphi$ is the largest $c$ with $1/c\ge c-1$, and the step bound $\varphi/(2L)$ follows from it. The fixed-step analysis is the template for the explicit, adaptive-step EGRAAL of the same paper (Theorem 2), which needs only local Lipschitz continuity of $F$.

The formalization. The theorem has a published proof, and no machine-checked version of it or of GRAAL is known. Mathlib contains the golden ratio, Lipschitz conditions, lower semicontinuity and cluster points, but no proximal operator of an extended-valued function, no prox-inequality and no Fejér-monotonicity convergence lemma. This mission produces those pieces and a complete convergence proof for a first-order VI method, which are reusable for projected gradient, forward–backward, extragradient and reflected-gradient analyses.

Difficulty

The naive approach, to show that $\|z^k-z^*\|$ decreases, fails: GRAAL is not Fejér monotone in $z^k$, because the forward step is taken from the average $\bar z^k$ and uses $F(z^k)$ rather than $F$ at the new point. The quantity that decreases is an energy mixing $\|\bar z^k-z^*\|^2$ with the successive difference $\|z^k-z^{k-1}\|^2$, and both the averaging identity and the Lipschitz estimate must produce matching coefficients for the cross terms to cancel. The energy inequality alone gives only boundedness and vanishing successive differences; convergence of the whole sequence, and the fact that the limit solves (1) when $g$ is merely lower semicontinuous and extended-valued, is a separate step. On the formal side, $g$ takes the value $+\infty$, so the prox-inequality and the variational inequality must be handled in extended arithmetic without letting $\infty-\infty$ decide anything.

Formalization scope

• $\mathcal E$ is a type E with [NormedAddCommGroup E] [InnerProductSpace ℝ E] [FiniteDimensional ℝ E].
• $g$ is E → EReal. (C2) is IsProperConvexLSC g: never $\bot$, somewhere finite, convex epigraph $\{(x,t)\in E\times\mathbb R: g(x)\le t\}$, and LowerSemicontinuous g on all of E. $\operatorname{dom} g$ is effDom g = {x | g x ≠ ⊤}.
• $F$ is a total function E → E; monotonicity and the Lipschitz bound $\|F(u)-F(v)\|\le L\|u-v\|$ are required on effDom g only. The step range is 0 < λ, λ ≤ φ / (2 * L) with 0 < L and φ = Real.goldenRatio.
• $S$ is solutionSet g F: points of effDom g satisfying (1) for every $z\in E$, evaluated in EReal.
• The proximal step is the argmin predicate IsProxPoint (fun x => λ * g x) w z⁺, not a choice function, so no junk value is involved. A run of (6) is IsGRAALRun g F λ z zbar on sequences ℕ → E indexed as in the paper: $z^1$ and $\bar z^0$ are free and the entry $z^0$ is unused.
• The conclusion is ∃ zs ∈ solutionSet g F, Tendsto z atTop (𝓝 zs) ∧ Tendsto zbar atTop (𝓝 zs).

The hypotheses of the goal are jointly satisfiable, so the theorem is not vacuous: for $g\equiv0$ and $F\equiv0$ every point is a solution and constant sequences form a run of (6); a formalization under which IsGRAALRun has no instances, or in which $S$ may be empty, is ruled out. Two hypotheses are added to printed statements and flagged in their notes: $C\neq\emptyset$ in Lemma 1, which is false without it, and $z^1\in\operatorname{dom} g$ in Eq. (14), needed at $k=2$ because the paper's $F$ is only defined on $\operatorname{dom} g$.

Welcome contributions: existence and uniqueness of the proximal point of a proper convex lsc function in finite dimensions; the prox-inequality; Fejér-monotonicity lemmas; the energy inequality; and the final convergence argument. The prox and Fejér infrastructure is independent of the golden ratio and is shared with the second mission of this series.

Selected references

• Y. Malitsky, Golden Ratio Algorithms for Variational Inequalities, preprint, Optimization Online 6598, 2018. https://optimization-online.org/wp-content/uploads/2018/05/6598.pdf ; published in Mathematical Programming. https://doi.org/10.1007/s10107-019-01416-w
• H. H. Bauschke, P. L. Combettes, Convex Analysis and Monotone Operator Theory in Hilbert Spaces, Springer, 2011 (2nd ed. 2017). https://doi.org/10.1007/978-3-319-48311-5
• G. M. Korpelevich, The extragradient method for finding saddle points and other problems, Ekonomika i Matematicheskie Metody 12 (1976) 747–756.
• P. Tseng, A modified forward–backward splitting method for maximal monotone mappings, SIAM J. Control Optim. 38 (2000) 431–446. https://doi.org/10.1137/S0363012998338806
• Y. Malitsky, Projected reflected gradient methods for monotone variational inequalities, SIAM J. Optim. 25 (2015) 502–520. https://doi.org/10.1137/14097238X

Motivation

Pattern search methods minimize a function $f:\mathbb R^n\to\mathbb R$ by comparing values of $f$ at points of a structured set of trial points, without evaluating or approximating derivatives. Coordinate search and the method of Hooke and Jeeves (Hooke–Jeeves 1961) are the classical members of the family. Such methods remain in use when derivatives are unavailable, unreliable or expensive, for instance when $f$ is the output of a simulation, and practical problems of this kind usually carry simple bounds on the variables.

Torczon (SIAM J. Optim. 1997) gave a global convergence theory for pattern search on unconstrained problems: under compactness of the level set and continuous differentiability of $f$, $\liminf_k\|\nabla f(x_k)\|=0$, and under stronger hypotheses $\lim_k\|\nabla f(x_k)\|=0$. Lewis and Torczon extended this theory to bound constrained problems (ICASE Report 96-20, 1996; SIAM J. Optim. 1999). The extension is not automatic: the paper exhibits a pattern search method for unconstrained problems (Box's evolutionary operation with factorial designs) that fails on bound constrained ones, and identifies the structural condition on the pattern that restores convergence.

Timeline.

• 1961: Hooke and Jeeves introduce "direct search" pattern methods.
• 1987–1988: Calamai and Moré (Math. Program. 1987) and Conn, Gould and Toint (SIAM J. Numer. Anal. 1988) develop the projected-gradient stationarity theory for bound and linear constraints, for methods that use derivatives.
• 1997: Torczon proves global convergence of generalized pattern search for unconstrained problems.
• 1996/1999: Lewis and Torczon prove the bound constrained theory formalized here.

Setting

The problem is

$$\min f(x)\quad\text{subject to}\quad \ell\le x\le u,$$

with $\ell,u$ vectors of extended reals and $\ell_j<u_j$ for every $j$; $\ell_j=-\infty$ or $u_j=+\infty$ is allowed. The feasible region is $\Omega=\{x:\ell\le x\le u\}$, $P$ is the coordinatewise projection onto $\Omega$, $g=\nabla f$, and $L_\Omega(y)=\{x\in\Omega:f(x)\le f(y)\}$ is the feasible level set. A stationary point is an $x\in\Omega$ with $\langle g(x),z-x\rangle\ge0$ for all $z\in\Omega$. The stationarity measure is

$$q(x)=P\bigl(x-g(x)\bigr)-x,$$

which vanishes exactly at stationary points.

A generalized pattern search method is fixed by a nonsingular basis matrix $B\in\mathbb R^{n\times n}$, a finite set $\mathcal M$ of nonsingular integer matrices, a rational $\tau>1$, an integer $w_0<0$ and nonnegative integers $w_1,\dots,w_L$. At iteration $k$ the generating matrix is $C_k=[M_k\ \ {-M_k}\ \ L_k]=[\Gamma_k\ \ L_k]$ with $M_k\in\mathcal M$, $L_k$ an integer matrix containing a zero column, and $BM_k$ diagonal. A trial step is $\Delta_kBc$ for a column $c$ of $C_k$. The step $s_k$ is a trial step with $x_k+s_k\in\Omega$, and it must decrease $f$ whenever some feasible trial step from the core $\Delta_kB\Gamma_k$ does. The iterate moves, $x_{k+1}=x_k+s_k$, exactly when $f(x_k+s_k)<f(x_k)$. The step length $\Delta_k$ is multiplied by $\theta=\tau^{w_0}<1$ after an unsuccessful iteration and by some $\tau^{w_i}\ge1$ after a successful one. The Strong Hypotheses additionally require $f(x_k+s_k)$ to be no larger than the best feasible core trial value whenever that value is below $f(x_k)$.

Formalization targets

Goal: Theorem 3.3

If $L_\Omega(x_0)$ is compact, $f$ is continuously differentiable, the columns of the $C_k$ are uniformly bounded, $\Delta_k\to0$, and the Strong Hypotheses hold, then

$$\lim_{k\to\infty}\|q(x_k)\|=0 .$$

Milestones

• Lemma 2.1, Theorem 2.2, Lemma 2.3: the unconstrained results the paper recalls from Torczon (1997): nonzero steps have length at least $\zeta_*\Delta_k$; the iterates lie on the translated lattice $x_0+\beta^{r_{LB}}\alpha^{-r_{UB}}\Delta_0B\,\mathbb Z^n$ (with $\tau=\beta/\alpha$); bounded columns give $\Delta_k\ge\psi_*\|s_k^i\|$.
• Proposition 3.1 (6), (8): $\|q(x)\|\le\|g(x)\|$, and $x$ is stationary iff $q(x)=0$.
• All iterates lie in $L_\Omega(x_0)$ (§4, p. 10).
• Propositions 4.1–4.3: a descent estimate along short steep directions; a feasible core step with $g_k^Ts\le-n^{-1/2}\|q_k\|\|s\|$ whenever $q_k\ne0$ and the step length is small; a uniform $\delta$ (and, under the Strong Hypotheses, a $\sigma$) with $f(x_{k+1})\le f(x_k)-\sigma\|q(x_k)\|\|s_k\|$ when $\Delta_k<\delta$ and $\|q(x_k)\|>\eta$.
• Corollary 4.4 and Theorem 4.5: $\liminf\|q(x_k)\|\ne0$ keeps $\Delta_k$ bounded away from zero, whereas compactness alone forces $\liminf\Delta_k=0$.
• Theorem 3.2: $\liminf_k\|q(x_k)\|=0$.

Significance

Theorem 3.2 shows that a method which never computes a gradient still has a subsequence approaching first-order stationarity for the bound constrained problem, even though it cannot enforce a sufficient decrease condition measured by the projected gradient. Theorem 3.3 upgrades this to the whole sequence, so every limit point of the iterates is a KKT point. These results justify the bound constrained variants of coordinate search and Hooke–Jeeves discussed in §5 of the paper, and they are the template for the later theory of pattern search under linear constraints and generating set search.

The results are proved in the paper, and three of the milestones are proved in Torczon (1997). None of them has a machine-checked proof. Formalizing them produces a Lean model of generalized pattern search (patterns, exploratory moves, step-length updates) that later missions on direct search, mesh adaptive direct search or linearly constrained pattern search can reuse, and checks the details the paper handles briefly: the lattice argument, the feasibility of the chosen coordinate step, and uniform constants.

Difficulty

The obvious argument copies the unconstrained proof with $\nabla f$ replaced by $q$. The step that fails is the existence of a good trial step: in the unconstrained case some pattern direction makes an acute angle with $-\nabla f(x_k)$, but near the boundary of $\Omega$ that direction may leave the feasible region, and a feasible direction may not be a descent direction. For a general pattern no uniform choice exists, and the paper's counterexample in §5.2 shows convergence can fail. The diagonality of $BM_k$ is what makes the pattern contain coordinate directions, one of which is both feasible and a descent direction of quality $n^{-1/2}\|q_k\|$ (Proposition 4.2). The second difficulty is Theorem 4.5, which uses no derivatives: it rests on the rationality of $\tau$ and the integrality of the $C_k$, which confine the iterates to a lattice that meets the compact set $L_\Omega(x_0)$ in finitely many points.

Formalization scope

Points are EuclideanSpace ℝ (Fin n) with the Euclidean norm; the bounds are Fin n → EReal with the hypothesis $\ell_j<u_j$ for all $j$, so infinite bounds are allowed as in the paper. Paper coordinates $1,\dots,n$ are Lean's Fin n. The gradient is Mathlib's gradient f. A run of the method is a structure of sequences $(x_k,\Delta_k,s_k,M_k,L_k)$ together with a predicate IsGPSRun that encodes §2.1–§2.4 clause by clause; the parameter $m\ge1$ is the number of columns of $L_k$ (the paper's $p-2n$). $\tau$ is rational and the $C_k$ are integer matrices, as the lattice argument requires. "$\min\{f(x_k+y):\dots\}<f(x_k)$" over the finite set of feasible core trial points is encoded as "some feasible core trial step strictly decreases $f$". $\liminf$ statements are encoded with ∃ᶠ, not Filter.liminf. The modulus of continuity $\omega$ is not formed as a real supremum; Proposition 4.1 takes an explicit radius $\delta>\|d\|$.

Standing assumptions and every departure from the page:

1. Smoothness. The page assumes $f$ continuously differentiable on $L_\Omega(x_0)$. The mission assumes $f$ is $C^1$ on an open set $U\supseteq\Omega$. The proofs evaluate $\nabla f$ along segments to trial points that lie in $\Omega$ but generally outside $L_\Omega(x_0)$, and $L_\Omega(x_0)$ may have empty interior, so the page's hypothesis does not define what the proofs use.
2. Strong Hypothesis 3. The page prints $f(x_k+s_k)<\min\{\cdots\}$. No core step can satisfy the strict form, which would exclude coordinate search, which the paper says satisfies it. The mission uses $\le$, the form of Torczon (1997) and the one the proof of Proposition 4.3 uses. The theorem with $\le$ implies the one with $<$.
3. The nonemptiness of $\{w_1,\dots,w_L\}$ is made explicit.
4. Proposition 3.1 (7) is omitted: its $P(g(x))$ is a projected gradient the paper does not define.
5. Proposition 4.2 quantifies over every step length below $\nu_k$, because $\nu_k$ does not depend on $\Delta_k$.

The run predicate is not vacuous: an explicit run of coordinate search on $f(x)=x$ over $[0,\infty)$ satisfies IsGPSRun, the Strong Hypotheses, bounded columns, $\Delta_k\to0$ and compactness of $L_\Omega(x_0)$. That check is proved in Lean without sorry, so the goal cannot be closed by exhibiting an unsatisfiable hypothesis.

A complete development needs the mean value theorem along segments, uniform continuity of $\nabla f$ near the compact set $L_\Omega(x_0)$, finiteness of a discrete lattice inside a compact set, and elementary facts about the coordinatewise projection. The projection and lattice lemmas are reusable beyond this mission. Proofs of any milestone, alternative proofs, and a general statement of Proposition 3.1 for closed convex $\Omega$ are welcome.

Selected references

• R. M. Lewis and V. Torczon, Pattern Search Algorithms for Bound Constrained Minimization, ICASE Report No. 96-20 (NASA CR-198306), 1996; SIAM J. Optim. 9(4):1082–1099, 1999. https://doi.org/10.1137/S1052623496300507
• V. Torczon, On the Convergence of Pattern Search Algorithms, SIAM J. Optim. 7(1):1–25, 1997. https://doi.org/10.1137/S1052623493250780
• P. H. Calamai and J. J. Moré, Projected Gradient Methods for Linearly Constrained Problems, Math. Program. 39:93–116, 1987. https://doi.org/10.1007/BF02592073
• A. R. Conn, N. I. M. Gould and P. L. Toint, Global Convergence of a Class of Trust Region Algorithms for Optimization with Simple Bounds, SIAM J. Numer. Anal. 25(2):433–460, 1988. https://doi.org/10.1137/0725029
• R. Hooke and T. A. Jeeves, "Direct Search" Solution of Numerical and Statistical Problems, J. ACM 8(2):212–229, 1961. https://doi.org/10.1145/321062.321069

Motivation

The maximum-flow problem asks how much of a commodity can be sent from a source to a sink through a network whose edges carry capacities. It is a basic model of operations research. Transportation, scheduling, bipartite matching and image segmentation reduce to it, and it is the inner step of many combinatorial algorithms.

Goldberg and Tarjan introduced the push–relabel (preflow) method in A New Approach to the Maximum-Flow Problem (J. ACM 35(4), 1988). Ford–Fulkerson-type algorithms augment along whole source–sink paths. The push–relabel method instead moves excess flow across single edges, guided by integer distance labels on the vertices. Whatever order its local operations are applied in, it is correct and performs $O(n^2 m)$ of them (§3 of the paper). Section 4 shows that one particular order, processing the active vertices first-in, first-out, cuts the dominant term, the number of nonsaturating pushes, to $O(n^3)$. The method and its FIFO and highest-label variants are the standard practical maximum-flow codes.

Timeline:

• 1956: Ford and Fulkerson, augmenting paths and max-flow min-cut.
• 1970–72: Dinic, and Edmonds and Karp, give polynomial augmenting-path bounds.
• 1974: Karzanov introduces preflows and obtains $O(n^3)$.
• 1982: Shiloach and Vishkin give a parallel $O(n^2 \log n)$ preflow algorithm with a first-in, first-out flavour.
• 1988: Goldberg and Tarjan, the generic push–relabel method, the FIFO bound of this mission, and $O(nm \log(n^2/m))$ with dynamic trees.

Setting

A flow network has a finite vertex set $V$ with $n = |V|$, a capacity $c(v,w) \ge 0$ on every ordered pair, a source $s$ and a sink $t \neq s$. The edges are the pairs with $c(v,w) > 0$, and there are no loops. A preflow is a function $f$ on vertex pairs with $f(v,w) \le c(v,w)$ and $f(v,w) = -f(w,v)$. Its excess $e(v) = \sum_u f(u,v)$ must be nonnegative at every $v \neq s$. The residual capacity is $r_f(v,w) = c(v,w) - f(v,w)$. A labeling $d$ assigns each vertex a value in $\mathbb{N} \cup \{\infty\}$. A vertex $v \notin \{s,t\}$ is active if $d(v) < \infty$ and $e(v) > 0$.

The two basic operations (Fig. 1 of the paper) are:

• push$(v,w)$, applicable when $v$ is active, $r_f(v,w) > 0$ and $d(v) = d(w)+1$. It sends $\delta = \min(e(v), r_f(v,w))$ from $v$ to $w$. The push is saturating if $r_f(v,w) = 0$ afterwards and nonsaturating otherwise.
• relabel$(v)$, applicable when $v$ is active and $d(v) \le d(w)$ for every residual edge $(v,w)$. It sets $d(v) \leftarrow \min\{d(w)+1 : r_f(v,w) > 0\}$.

The algorithm starts by saturating every edge leaving $s$, with $d(s) = n$ and $d(v) = 0$ for $v \neq s$.

In the first-in, first-out algorithm (§4), each vertex $v$ scans a fixed list $L(v)$ of its neighbours through a current edge. The push/relabel$(v)$ operation pushes through the current edge if possible. Otherwise it advances the current edge, or, at the end of the list, returns to the first edge and relabels $v$. Active vertices wait in a queue $Q$, initially $\{v \in V - \{s,t\} : c(s,v) > 0\}$. The discharge operation removes the front vertex $v$ and repeats push/relabel$(v)$ until $e(v) = 0$ or $d(v)$ increases. Every vertex that becomes active meanwhile is appended to $Q$, and $v$ is appended too if it is still active. Passes over the queue are defined inductively. Pass 1 consists of the discharges of the initially queued vertices. Pass $i+1$ consists of the discharges of vertices added during pass $i$.

Formalization targets

Goal: Corollary 4.4 (p. 931)

For every network, every edge-list order, every initial queue order, and every run of the FIFO algorithm,

$$\#\{\text{nonsaturating pushes}\} \le 4n^3 .$$

The constant is the printed one.

Milestones

• Lemma 4.1 (p. 929): the push/relabel operation relabels only when relabeling is applicable.
• Lemma 3.5 (p. 926): from any vertex with positive excess, the source is reachable in the residual graph.
• Lemma 3.7 (p. 927): at any time, $d(v) \le 2n-1$ for every vertex.
• Lemma 3.8 (p. 927): at most $2n-1$ relabelings per vertex and at most $(2n-1)(n-2) < 2n^2$ in total.
• Lemma 4.3 (p. 930): at most $4n^2$ passes over the queue.

Significance

Corollary 4.4 is the combinatorial core of Theorem 4.5, which states that the FIFO algorithm runs in $O(n^3)$ time. Theorem 4.2 shows that the remaining work of the implementation is $O(nm)$ plus constant time per nonsaturating push. The bound of Corollary 4.4 is therefore what separates the $O(n^3)$ FIFO method from the $O(n^2 m)$ bound of the generic method, which matters on dense networks. The same pass-counting argument is reused for the parallel algorithm of §6 and underlies later analyses of highest-label and wave variants.

The results are proved in the paper. Formalizing them adds an analysis of a push–relabel algorithm, which the platform does not yet have. Its existing network-flow material states max-flow min-cut and Ford–Fulkerson termination in an arc-based model with nonnegative flows (the Introduction to Linear Optimization missions). The mission builds a precise operational model of the FIFO implementation, with edge lists, current edges and a queue carrying pass numbers, and states an explicit operation count for it. A companion mission in this series treats the generic algorithm's correctness and its $(2n-1)(n-2) + 2nm + 4n^2m$ operation bound.

Difficulty

The obvious argument is the potential-function count of §3, over the sum of the labels of active vertices. It yields only $4n^2 m$ and does not use the queue discipline at all. The $4n^3$ bound has to charge nonsaturating pushes to passes over the queue, and then bound the number of passes by the total growth of the labels. Neither step is visible in the generic algorithm, because both depend on the order in which vertices are processed.

Making this rigorous requires invariants of the implementation that the paper uses silently:

• a vertex is in the queue exactly when it is active, and at most once;
• pass numbers are nondecreasing along the queue;
• current edges only move forward between relabelings.

Lemma 4.1 in particular depends on the current-edge scan: an edge passed over earlier stays inadmissible until $v$ is relabeled.

Formalization scope

The Lean development works in namespace GoldbergTarjan.FIFO. Vertices form a type V with [Fintype V] [DecidableEq V], and $n$ is Fintype.card V. Capacities are c : V → V → ℝ with c ≥ 0 and c v v = 0. Flows are antisymmetric real functions on all ordered pairs, not nonnegative arc flows. Excess is computed from the flow. Labels are in ℕ∞, and the empty minimum in relabel is ⊤.

The state of the algorithm consists of the flow, the labels, the current-edge index cur v into the edge list L v, and the queue Q : List (V × ℕ), each entry tagged with its pass number. Push/relabel (Fig. 3) is a total function, and a discharge (Fig. 4) is a relation carrying the number of push/relabel operations it performs. A run consists of the states S 0, …, S K with S 0 the initial state and consecutive states related by one discharge. The printed variant of Fig. 4, which stops as soon as $v$ is relabeled, is the one formalized. Counts are natural numbers over all push/relabel operations of all discharges. The number of passes is the largest pass tag of a discharged entry.

All constants are explicit, exactly as printed:

• $2n-1$ (Lemmas 3.7, 3.8);
• $(2n-1)(n-2)$ and $2n^2$ (Lemma 3.8);
• $4n^2$ (Lemma 4.3);
• $4n^3$ (Corollary 4.4).

No asymptotic notation is used, and no $m \ge n-1$ assumption is made.

A model without current edges, where relabeling happens whenever no push applies, would make Lemma 4.1 vacuous and change the algorithm. Pass numbers that are not propagated by the "added during pass $i$" rule would make the pass count arbitrary. Both are ruled out by the definitions. A sorry-free check, outside the proposal, exhibits a three-vertex network with two legal discharges, two passes and no nonsaturating push, so the run hypotheses are satisfiable.

Reusable beyond this mission are the network, preflow, push and relabel definitions and Lemma 3.5, which is about an arbitrary preflow. Contributions welcome: invariants of FIFO runs (preflow, valid labeling, queue = active set, cur within bounds), proofs of the milestones, and the reduction of Corollary 4.4 to Lemma 4.3.

Selected references

• A. V. Goldberg, R. E. Tarjan, A New Approach to the Maximum-Flow Problem, Journal of the ACM 35(4):921–940, 1988. https://doi.org/10.1145/48014.61051
• A. V. Karzanov, Determining the maximal flow in a network by the method of preflows, Soviet Math. Doklady 15:434–437, 1974.
• Y. Shiloach, U. Vishkin, An O(n² log n) parallel max-flow algorithm, Journal of Algorithms 3(2):128–146, 1982. https://doi.org/10.1016/0196-6774(82)90013-X
• L. R. Ford, D. R. Fulkerson, Maximal flow through a network, Canadian Journal of Mathematics 8:399–404, 1956. https://doi.org/10.4153/CJM-1956-045-5
• J. Edmonds, R. M. Karp, Theoretical improvements in algorithmic efficiency for network flow problems, Journal of the ACM 19(2):248–264, 1972. https://doi.org/10.1145/321694.321699