All missions

Motivation

Quasi-Newton methods minimize a smooth function $f:\mathbb R^n\to\mathbb R$ using gradients only, while building up an approximation to the inverse of the Hessian matrix from the changes of gradient they observe. They are the default algorithms of unconstrained nonlinear optimization: BFGS and its limited-memory variant L-BFGS sit inside most numerical optimization libraries and most large-scale fitting codes in statistics and machine learning.

The family starts with William C. Davidon's Argonne report ANL-5990 of 1959, Variable Metric Method for Minimization, rejected by a journal in 1957 and published only in 1991 as the first article of the SIAM Journal on Optimization, with a preface by the author (DOI 10.1137/0801001). The report's body describes a flowchart algorithm; its Appendix describes "a simplified method embodying some of the ideas" of that algorithm, in one and a half pages. That simplified method is what is now called the symmetric rank-one (SR1) update, and the Appendix contains its quadratic-termination property: on a quadratic, after at most $n$ steps the trial matrix equals the inverse Hessian, and the next step lands at the minimum.

Timeline.

• 1952: Hestenes and Stiefel's conjugate gradient method minimizes a strictly convex quadratic in at most $n$ steps (J. Res. NBS 49); Davidon cites it as [2].
• 1959: Davidon's ANL-5990, including the Appendix method treated here.
• 1963: Fletcher and Powell simplify the body's method into the DFP update and prove its quadratic termination with exact line searches (Comput. J. 6).
• 1967: Broyden's survey of quasi-Newton methods discusses the symmetric rank-one update (Math. Comp. 21); textbooks later analyse it under the name SR1.
• 1980: Nocedal's limited-storage BFGS, whose termination on quadratics with exact line searches is a separate Prove2Me mission (Math. Comp. 35).

Setting

Vectors are elements of $\mathbb R^n$ and matrices are real $n\times n$. Fix a symmetric positive definite matrix $G$ and a point $\xi$, and let

$$f(x)=\tfrac12\,(x-\xi)^{\mathsf T}G\,(x-\xi)+c ,$$

a quadratic with constant Hessian $G$ and minimum point $\xi$. Its gradient is $\nabla(x)=G(x-\xi)$, display (6) of the Appendix.

The method keeps a point $x_k$ and a symmetric trial matrix $H_k$, which plays the role of $G^{-1}$. Write $\nabla_k=\nabla(x_k)$. One iteration is

$$x_{k+1}=x_k-H_k\nabla_k \quad\text{(2)},\qquad H_{k+1}=H_k+a_k\,(H_k\nabla_{k+1})(H_k\nabla_{k+1})^{\mathsf T}\quad\text{(3)},$$

a full step with no line search followed by a symmetric rank-one update. With the two scalars

$$N_k=\nabla_{k+1}^{\mathsf T}H_k\nabla_{k+1},\qquad M_k=\nabla_{k+1}^{\mathsf T}H_k\nabla_k ,$$

footnote 2 of the Appendix fixes the coefficient on a quadratic to $a_k=(M_k-N_k)^{-1}$, which makes the paper's quality measure $\Delta$ of (4) vanish. Write $S_k=x_{k+1}-x_k$ for the step and $D_k=\nabla_{k+1}-\nabla_k$ for the change of gradient. The paper writes $x$, $H$, $\nabla$ for the current quantities and $x^{+}$, $H^{+}$, $\nabla^{+}$ for the updated ones, and uses the letter $N$ both for the number of variables and for the scalar $N_k$; here the number of variables is $n$.

Formalization targets

Goal: termination in $n$ steps

If $H_0$ is symmetric and $M_k\ne N_k$ for every $k<n$, then

$$H_n=G^{-1}\qquad\text{and}\qquad x_{n+1}=\xi .$$

This is the sentence after (8) on p. 17: "After no more than $N$ steps (for which $\Delta=0$), $H$ will equal $G^{-1}$ and the following step will be to the exact minimum."

Milestones

1. Display (6): the gradient of $f$ is $G(x-\xi)$.
2. Display (8): with $a=(M-N)^{-1}$ and $M\ne N$, $H^{+}GS=H^{+}D=S$.
3. Display (7): if $HGu=u$ then $H^{+}Gu=u$, for any coefficient $a$.
4. The sentence "so that $S$ becomes another such eigenvector", read along the run: if $M_i\ne N_i$ for $i<k$, then $H_kGS_j=S_j$ for every $j<k$.

Further statements of the Appendix

• Condition 1 (p. 16): if $H$ is positive definite and $R^{-1}\le\det H^{+}/\det H\le R$ with $R>1$, then $H^{+}$ is positive definite.
• The constraint remark (p. 17): for any function and any coefficients, if $H_0$ is symmetric and $H_0b=0$, then every step is perpendicular to $b$ and $b\cdot x_k$ is conserved.
• Table (5) (p. 16): for general (non-quadratic) $f$, the coefficient $a$ that minimizes $\Delta$ subject to condition 1, and the minimum value, on each of five ranges of $M$.

Significance

The termination theorem says that on a quadratic the rank-one update recovers the exact inverse Hessian from $n$ gradient differences, without line searches and without positive definiteness of the trial matrix, and then takes the Newton step. This is the property that justifies the name "variable metric": the metric $H_k$ converges to the one that makes the problem trivial. The secant equation (8) and the persistence of earlier secant equations (7) are the template for the analysis of every later quasi-Newton update, and the SR1 update remains in use in trust-region methods because it does not force positive definiteness.

The result is classical and its proof is short; what this mission adds is a machine-checked version stated in the paper's own form. In particular the hypothesis is the paper's: only that each of the first $n$ coefficients is defined. Textbook statements of SR1 termination often assume in addition that the steps are linearly independent; the paper does not, and the formal goal does not either. No Lean formalization of SR1, DFP or BFGS termination was found in Mathlib or on Prove2Me. On Prove2Me, the related open targets are L-BFGS termination (Nocedal 1980) and conjugate-gradient termination (Hestenes–Stiefel 1952); the definitions here are independent of theirs.

Difficulty

Each update changes $H$ by a rank-one term, and nothing in a single step forces $H_n$ to be $G^{-1}$: the coefficient $a_k$ depends on the current iterate and the update could, a priori, destroy what earlier steps achieved, as general rank-one updates do. The obvious count — $n$ steps, $n$ conditions — does not by itself give $n$ independent conditions on $H_n$; the steps produced by the method are not chosen to be independent, and the only hypothesis available is that the denominators $M_k-N_k$ are nonzero. The difficulty is to get the full conclusion from that hypothesis alone, without assuming independence, positive definiteness of $H_0$, or a line search.

Formalization scope

Vectors are Fin n → ℝ with 0-based indices, matrices Matrix (Fin n) (Fin n) ℝ; $u^{\mathsf T}Hv$ is u ⬝ᵥ (H *ᵥ v) and $(Hg)(Hg)^{\mathsf T}$ is vecMulVec (H *ᵥ g) (H *ᵥ g). The definition file DavidonVM.Termination.Method provides quad, grad, Nval, Mval, the run predicate IsRunWith for an arbitrary gradient field and coefficients, and IsQuadRun for the quadratic with $a_k=(M_k-N_k)^{-1}$; DavidonVM.Termination.Delta provides the update and the quantity $\Delta$ of (4).

Committed conventions, each a reading of something the paper leaves implicit:

• The function is an explicit quadratic with constant symmetric positive definite $G$ ("in the neighborhood of a minimum"); the run's gradient field is $G(x-\xi)$, and milestone 1 ties it to the quadratic. The goal uses positive definiteness; milestones use only symmetry where that suffices; the secant milestone needs nothing on $G$.
• Only $H_0$ is assumed symmetric (the paper, p. 4: "It is to be symmetric"); positive definiteness of $H_0$ is not assumed.
• "Steps for which $\Delta=0$" is read through footnote 2 as $M_k\ne N_k$. Lean's real inverse gives $0^{-1}=0$, so without this hypothesis the update would silently leave $H$ unchanged.
• The stopping rule "provided that $N$ is greater than some preassigned $\varepsilon$" is not modelled; the run continues for every $k$, and the goal looks at steps $0,\dots,n$ only.
• $G^{-1}$ is Mathlib's matrix inverse, the true inverse since $G$ is positive definite.

A statement of the goal with the extra hypothesis "$S_0,\dots,S_{n-1}$ are linearly independent", or with $H_k G S_j=S_j$ assumed, would be a weaker theorem than the paper's and is ruled out. The case $n=0$ is trivial; the hypotheses are jointly satisfiable for $n=1$ (checked locally).

Welcome contributions: proofs of the four milestones and the goal; general facts about rank-one updates (determinant and positive definiteness, Sherman–Morrison form of the inverse) that the further statements need and that are reusable for DFP, BFGS and trust-region analyses.

Selected references

• W. C. Davidon, Variable Metric Method for Minimization, SIAM J. Optim. 1(1) (1991) 1–17 (Argonne report ANL-5990, 1959). https://doi.org/10.1137/0801001
• M. R. Hestenes and E. Stiefel, Methods of conjugate gradients for solving linear systems, J. Res. Nat. Bur. Standards 49 (1952) 409–436. https://doi.org/10.6028/jres.049.044
• R. Fletcher and M. J. D. Powell, A rapidly convergent descent method for minimization, Comput. J. 6 (1963) 163–168. https://doi.org/10.1093/comjnl/6.2.163
• C. G. Broyden, Quasi-Newton methods and their application to function minimisation, Math. Comp. 21 (1967) 368–381. https://doi.org/10.1090/S0025-5718-1967-0224273-2
• J. Nocedal, Updating quasi-Newton matrices with limited storage, Math. Comp. 35 (1980) 773–782. https://doi.org/10.1090/S0025-5718-1980-0572855-7

Motivation

Bounded-confidence opinion dynamics model a population in which each individual revises an opinion, represented by a real number, by looking only at opinions close to its own. The best-known example is the discrete-time model of Krause (1997), also called the Hegselmann–Krause model, in which every agent moves to the average of the opinions within distance 1 of its own. Such models are studied in the social sciences as models of consensus formation and polarization, and in control theory as prototypes of multiagent systems with state-dependent interaction topology: the communication graph is not given in advance but is determined by the state itself, which is what rendezvous algorithms for mobile robots and flocking models share (Lorenz 2007; Olfati-Saber, Fax, Murray 2007).

Blondel, Hendrickx and Tsitsiklis (SIAM J. Control Optim. 48 (2010)) study the continuous-time, symmetric counterpart of Krause's model. Each opinion is continuously attracted by every other opinion differing from it by less than 1, with an intensity proportional to the difference. Simulations show convergence to clusters, groups of agents sharing a common value, with distinct clusters at distance at least 1 and typically close to 2. This mission formalizes the first of the paper's results: convergence to clusters at least 1 apart for finitely many agents.

Timeline.

• 1997: Krause introduces the discrete-time bounded-confidence model.
• 2003: Jadbabaie, Lin and Morse prove convergence results for nearest-neighbour coordination with switching topologies (IEEE TAC 48).
• 2009: Blondel, Hendrickx and Tsitsiklis analyse Krause's discrete-time model and the observed intercluster distance of about 2 (IEEE TAC 54).
• 2010: the present paper. For the continuous-time model with n agents it proves that every proper solution converges to an equilibrium with clusters at least 1 apart (Theorem 2). It also proves that almost every initial condition is proper (Theorem 1, proof sketched) and gives a continuum-agent analysis with a sharper intercluster bound.

Setting

There are $n$ agents, labelled $1,\dots,n$. Agent $i$ holds a real opinion $x_i(t)$ at each time $t\ge 0$, and $x_i$ is a continuous function of time. Agent $j$ is a neighbour of agent $i$ at time $t$ when $|x_i(t)-x_j(t)|<1$; the inequality is strict, and every agent is its own neighbour. The intended dynamics are

$$\dot x_i(t)=\sum_{j:\,|x_i(t)-x_j(t)|<1}\bigl(x_j(t)-x_i(t)\bigr).\tag{1.1}$$

The right-hand side jumps whenever a pair of opinions crosses distance 1, so (1.1) usually has no differentiable solution. The model is therefore the integral equation

$$x_i(t)=x_i(0)+\int_0^t\sum_{j:\,|x_i(\tau)-x_j(\tau)|<1}\bigl(x_j(\tau)-x_i(\tau)\bigr)\,d\tau\qquad(t\ge0,\ i=1,\dots,n).\tag{2.1}$$

Solutions of (2.1) need not be unique. The two-agent initial condition $(-\tfrac12,\tfrac12)$ is a fixed point, but the agents can also start attracting each other immediately. An initial condition $\tilde x\in\mathbb R^n$ is proper if (a) (2.1) has exactly one solution $x$ with $x(0)=\tilde x$; (b) the set of times at which $x$ is not differentiable is at most countable and has no accumulation point; and (c) agents that meet stay together: $x_i(t)=x_j(t)$ implies $x_i(t')=x_j(t')$ for all $t'\ge t$. That solution is then a proper solution.

The set of equilibria is

$$F=\bigl\{\tilde s\in\mathbb R^n:\ \text{for all } i,j,\ \tilde s_i=\tilde s_j\ \text{or}\ |\tilde s_i-\tilde s_j|\ge 1\bigr\}.$$

The average opinion is $\bar x(t)=\frac1n\sum_i x_i(t)$, and $V(x(t))=\sum_i\bigl(x_i(t)-\bar x(t)\bigr)^2$ is the sum of squared differences from it.

Formalization targets

Goal: Theorem 2

Every proper solution $x$ of (2.1) converges to a limit in $F$:

$$\exists\,x^*\in F:\qquad \lim_{t\to\infty}x(t)=x^*,$$

that is, all limits exist and two agents with different limits end at distance at least 1.

Milestones, in the order the argument uses them

1. Order preservation (§2.1, p. 5218): if $x_i(t)\ge x_j(t)$ for some $t$, then $x_i(t')\ge x_j(t')$ for every $t'\ge t$.
2. Proposition 1 (p. 5218): $\bar x$ is constant, $V(x(t))$ is nonincreasing, and outside a countable set of times $dV/dt$ is negative when $x(t)\notin F$ and zero when $x(t)\in F$.
3. Monotone partial sums (2.3) (p. 5219): for a sorted initial condition and every $k$, outside a countable set of times the derivative of $\sum_{i\le k}x_i(t)$ equals $\sum_{i\le k}\sum_{j>k:\,|x_i-x_j|<1}(x_j-x_i)\ge 0$, and $\sum_{i\le k}x_i(t)$ is nondecreasing in $t$.
4. Convergence of each opinion (p. 5219): every $x_i(t)$ has a limit as $t\to\infty$.

Significance

The result. Theorem 2 establishes the clustering behaviour seen in simulations for every proper solution. The paper's Theorem 1 states that almost every initial condition is proper, so the conclusion covers almost every initial condition. The lower bound 1 on the distance between clusters is the baseline that the paper's continuum analysis improves to $1+\min\{W_A,W_B\}/\max\{W_A,W_B\}$. Average preservation and the decrease of $V$ (Proposition 1) are the structural facts that the continuum model also has.

Formalization. Theorem 2 is proved in the paper; no machine-checked version exists on Prove2Me or, as far as is known, elsewhere. The mission produces a Lean definition of solutions of a discontinuous integral equation with a state-dependent neighbour graph, and a convergence proof for it. The remaining work is to formalize the known argument, from the integral equation to the existence and location of the limit.

Difficulty

The obvious approach treats (1.1) as an ODE and differentiates along trajectories. That fails because the right-hand side is discontinuous in the state. A solution is only an integral solution, it may fail to be differentiable, and derivative identities such as (2.3) hold only outside an exceptional set of times. Every step from a derivative sign to a monotonicity or limit statement therefore has to go through the integral equation, not the differential one. A second obstacle is that pairs of agents at distance exactly 1 sit on the discontinuity of the interaction rule. This is where uniqueness fails (footnote 3 of the paper) and why conditions (a)–(c) are hypotheses. A decreasing Lyapunov function $V$ alone does not give convergence of the trajectory to a single point, so Proposition 1 does not by itself settle Theorem 2.

Formalization scope

Agents are Fin n. Time is ℝ, and only values at $t\ge 0$ enter. A trajectory is x : ℝ → Fin n → ℝ. A solution of (2.1) is continuous on $[0,\infty)$ and satisfies (2.1) for every $t\ge 0$ and every $i$, with the neighbour set Finset.univ.filter (fun j => |x τ i - x τ j| < 1). The integrand is also required to be interval-integrable on $[0,t]$; for a continuous trajectory this always holds. Explicit readings of the paper's phrases:

• "unique solution" in (a) is uniqueness among all solutions of (2.1), with no sortedness or differentiability assumed;
• condition (b) is "for every $T$, the non-differentiability times in $[0,T]$ form a finite set", which is equivalent for subsets of $[0,\infty)$;
• "converges to a limit" is Tendsto x atTop (𝓝 xstar) in Fin n → ℝ, with real time;
• "with the exception of a countable set of times" is the existence of a countable set $S$ outside which the derivative of $V(x(t))$ exists and has the stated sign;
• "constant" and "nonincreasing" are on $[0,\infty)$;
• Proposition 1 assumes $n\ge1$ so that $\bar x$ is meaningful.

The paper's sortedness convention ("we assume … that the components of proper initial conditions are sorted") is a hypothesis only of the partial-sum milestone. Theorem 2 is stated for every proper solution. A formalization that dropped condition (a), restricted solutions to sorted, differentiable or eventually constant trajectories, or replaced (2.1) by the differential equation would prove a different, and in parts trivial, statement. These are excluded.

The work needs interval integrals and the fundamental theorem of calculus for integrands with jumps, monotonicity from almost-everywhere derivative signs, and bounded monotone convergence. The definition of solutions of (2.1) and the equilibrium set is reusable for the stability results of the same paper (Theorem 3), which are not part of this mission. Contributions are welcome at every milestone, as are alternative convergence proofs, such as the paper's reference [9], that avoid symmetry.

Selected references

• V. D. Blondel, J. M. Hendrickx, J. N. Tsitsiklis, Continuous-time average-preserving opinion dynamics with opinion-dependent communications, SIAM J. Control Optim. 48(8), 2010, 5214–5240. https://doi.org/10.1137/090766188
• V. D. Blondel, J. M. Hendrickx, J. N. Tsitsiklis, On Krause's multi-agent consensus model with state-dependent connectivity, IEEE Trans. Automat. Control 54(11), 2009, 2586–2597. https://doi.org/10.1109/TAC.2009.2031211
• A. Jadbabaie, J. Lin, A. S. Morse, Coordination of groups of mobile autonomous agents using nearest neighbor rules, IEEE Trans. Automat. Control 48(6), 2003, 988–1001. https://doi.org/10.1109/TAC.2003.812781
• J. Lorenz, Continuous opinion dynamics under bounded confidence: a survey, Internat. J. Modern Phys. C 18(12), 2007, 1819–1838. https://doi.org/10.1142/S0129183107011789
• R. Olfati-Saber, J. A. Fax, R. M. Murray, Consensus and cooperation in networked multi-agent systems, Proc. IEEE 95(1), 2007, 215–233. https://doi.org/10.1109/JPROC.2006.887293
• U. Krause, Soziale Dynamiken mit vielen Interakteuren. Eine Problemskizze, in Modellierung und Simulation von Dynamiken mit vielen interagierenden Akteuren, Universität Bremen, 1997, 37–51.

Motivation

Many optimization problems in machine learning, signal processing and distributed computation have an objective that is a sum of a large number $m$ of component functions: one per data point, per sensor, or per processor. Evaluating the whole sum, or a subgradient of it, at every iteration is expensive when $m$ is large. Incremental methods process one component at a time, cycling through the components or sampling them, and update the iterate after each one. Incremental gradient and subgradient methods go back to the training of neural networks and to least squares; incremental proximal methods replace the subgradient step of a component by a proximal step, which is more stable and handles nonsmooth terms exactly.

Bertsekas (LIDS-P-2847, 2010, rev. 2011; Math. Program. 129 (2011)) proposed a unified family that mixes the two: each component is split as $F_i=f_i+h_i$, a proximal step is taken on $f_i$ and a subgradient step on $h_i$, followed by a projection onto the constraint set. This mission formalizes the convergence analysis of that family under the cyclic order of component selection.

Timeline. Incremental subgradient methods were analysed by Nedić and Bertsekas (SIAM J. Optim. 12 (2001)) for cyclic and randomized orders. Incremental proximal methods were studied by Bertsekas (Optimization for Machine Learning, 2011, ch. 4). The 2011 paper combines the two analyses; its cyclic-order results are Propositions 3–6.

Setting

Let $\mathbb R^n$ carry the Euclidean norm $\|\cdot\|$ and inner product $x'y$. The data are a nonempty closed convex set $X\subseteq\mathbb R^n$ and real-valued convex functions $f_i,h_i:\mathbb R^n\to\mathbb R$, $i=1,\dots,m$. The problem is

$$\text{minimize } F(x)=\sum_{i=1}^m F_i(x),\qquad F_i=f_i+h_i,\qquad\text{subject to } x\in X,$$

with optimal value $F^*=\inf_{x\in X}F(x)$, possibly $-\infty$, and optimal set $X^*=\{x^*\in X:F(x^*)=F^*\}$, possibly empty.

A subgradient $\tilde\nabla f(x)$ of a convex $f$ at $x$ is a vector $g$ with $f(y)\ge f(x)+g'(y-x)$ for all $y$; $P_X(u)$ is the point of $X$ nearest to $u$. At iteration $k$ one component $i_k$ is processed with stepsize $\alpha_k>0$. The three iterations are

$$\begin{aligned} &(19)\quad z_k=P_X\big(x_k-\alpha_k\tilde\nabla f_{i_k}(z_k)\big),&& x_{k+1}=P_X\big(z_k-\alpha_k\tilde\nabla h_{i_k}(z_k)\big),\\ &(20)\quad z_k=x_k-\alpha_k\tilde\nabla f_{i_k}(z_k),&& x_{k+1}=P_X\big(z_k-\alpha_k\tilde\nabla h_{i_k}(z_k)\big),\\ &(21)\quad z_k=x_k-\alpha_k\tilde\nabla h_{i_k}(x_k),&& x_{k+1}=P_X\big(z_k-\alpha_k\tilde\nabla f_{i_k}(x_{k+1})\big). \end{aligned}$$

The $f$-subgradient is taken at the new point, which makes the $f$-part a proximal step (Proposition 1). In the cyclic order $i_k=(k\bmod m)+1$; a block of $m$ consecutive iterations is a cycle, and $\alpha_k$ is constant within a cycle.

The analysis assumes that the subgradients the iteration uses are bounded by a constant $c$ and that, at each cycle start $k>0$, the function values at the intermediate points of the cycle differ from those at $x_k$ by at most $c$ times the distance (Assumption 1 for (19), (20); Assumption 2 for (21)). Both hold, for example, when all $f_i$, $h_i$ are polyhedral or Lipschitz.

Formalization targets

Goal: Proposition 6

If $\alpha_k\to0$ and $\sum_k\alpha_k=\infty$, then

$$\liminf_{k\to\infty}F(x_k)=F^*,$$

and if moreover $X^*\neq\emptyset$ and $\sum_k\alpha_k^2<\infty$, then $x_k\to x^*$ for some $x^*\in X^*$. Both parts are one item, as in the paper.

Milestones

1. Proposition 1: the proximal iteration over $X$ equals $P_X(x_k-\alpha_k\tilde\nabla f(x_{k+1}))$, and $\|x_{k+1}-y\|^2\le\|x_k-y\|^2-2\alpha_k(f(x_{k+1})-f(y))-\|x_k-x_{k+1}\|^2$ for $y\in X$.
2. Eq. (30) and Eq. (37): one-iteration estimates for (19)/(20) and for (21).
3. Proposition 3: over one cycle, $\|x_{k+m}-y\|^2\le\|x_k-y\|^2-2\alpha_k(F(x_k)-F(y))+\alpha_k^2\beta m^2c^2$, with $\beta=\frac1m+4$ for (19), (20) and $\beta=\frac5m+4$ for (21).
4. Footnote 5: nonnegative sequences with $Y_{k+1}\le Y_k-Z_k+W_k$ and $\sum W_k<\infty$ have $Y_k$ convergent.
5. Proof of Proposition 6: $\|x_{k+1}-x_k\|\le2\alpha_kc$.
6. Proposition 4 (constant stepsize $\alpha$): $\liminf F(x_k)=F^*$ if $F^*=-\infty$, and $\liminf F(x_k)\le F^*+\alpha\beta m^2c^2/2$ otherwise.

Significance

Proposition 6 is the exact-convergence guarantee for the cyclic incremental subgradient-proximal methods: with a diminishing, non-summable stepsize they reach the optimal value whatever the starting point, and with square-summable stepsizes the iterates converge to a single optimal solution. It covers the pure incremental subgradient method ($f_i=0$) and the pure incremental proximal method ($h_i=0$) as special cases, and it is the template for the randomized-order analysis of the companion mission. Proposition 4 quantifies what a constant stepsize achieves.

These results are proved in the paper. None of them, nor the underlying estimates, is formalized on the platform or, to our knowledge, in Mathlib. The mission produces machine-checked versions of the proximal-step lemma for extended-real-valued functions on a constrained set, of the cycle estimate (27) with its explicit constants, and of the deterministic supermartingale convergence lemma, all reusable for other incremental and proximal methods.

Difficulty

The obvious argument treats one cycle as one step of a subgradient method on $F$. It fails because within a cycle the subgradients are taken at the intermediate points $z_{k+j-1}$ or $x_{k+j-1}$, not at $x_k$, and the proximal steps use subgradients at the points they produce. The error between $F_j(x_k)$ and $F_j$ at those points must be controlled through the drift of the iterates during the cycle, and the bookkeeping differs between (19)/(20) and (21), which is why $\beta$ differs. For the second part of Proposition 6, convergence of the values in the $\liminf$ sense does not by itself give convergence of the iterates, and $X^*$ may contain more than one point.

Formalization scope

• $\mathbb R^n$ is EuclideanSpace ℝ (Fin n). Components are Fin m (0-based), so $i_k$ is k % m, and the paper's $j=1,\dots,m$ with $z_{k+j-1}$ becomes $j=0,\dots,m-1$ with z (k + j). $m\ge1$ is assumed (NeZero m).
• $P_X$ is a predicate (IsProj: a nearest point of $X$), not a choice function.
• Each iteration is a relation Step between $x_k$, $z_k$, $x_{k+1}$ and the two subgradients it uses, with the $f$-subgradient at the point the paper names. A run records these subgradients, so Assumptions 1–2 bound exactly them. The assumptions are not replaced by Lipschitz continuity, which is only a sufficient condition (p. 9).
• $F^*$ and $\liminf F(x_k)$ are taken in EReal; $F^*=-\infty$ is allowed and no boundedness is assumed.
• The standing assumptions of p. 4 (X nonempty closed convex, $f_i,h_i$ convex), positive stepsizes, constancy within a cycle and Assumptions 1–2 (stated "throughout this section", p. 8) are explicit hypotheses of every statement, including the goal, whose sentence does not repeat them.
• The starting point $x_0$ is arbitrary, as in the paper. Assumptions 1–2 make the cycle-start conditions only for $k>0$, and the bound $\|x_{k+1}-x_k\|\le2\alpha_kc$ needs $x_k\in X$, which holds from $k=1$ on. Proposition 3 and the step bound are therefore stated for $k>0$ and $k\ge1$; this does not affect Propositions 4 and 6.
• Proposition 1 is stated for extended-real-valued closed proper convex $f$, as in the paper, using the published class $\Gamma_0$ (MoreauProx.Characterization.GammaZero) and its subgradients, with relative interiors as Mathlib's intrinsicInterior.

The goal cannot be satisfied vacuously: a sorry-free check confirms that the hypotheses hold together (for $n=1$, $m=2$, $X=\mathbb R$, $\alpha_k=1/(\lfloor k/2\rfloor+1)$, $X^*\neq\emptyset$). Restricting to Lipschitz components, to bounded-below $F$, or to $x_0\in X$ would state a weaker theorem and is ruled out.

Proofs of any item are welcome, as are general lemmas (the projection theorem as a variational inequality, the subdifferential sum rule under the relative-interior condition) that the proofs need.

Selected references

• D. P. Bertsekas, Incremental Proximal Methods for Large Scale Convex Optimization, Report LIDS-P-2847, MIT, 2010 (rev. March 2011); Math. Program. 129(2) (2011) 163–195. https://doi.org/10.1007/s10107-011-0472-0
• A. Nedić, D. P. Bertsekas, Incremental Subgradient Methods for Nondifferentiable Optimization, SIAM J. Optim. 12(1) (2001) 109–138. https://doi.org/10.1137/S1052623499362111
• D. P. Bertsekas, Incremental Gradient, Subgradient, and Proximal Methods for Convex Optimization: A Survey, in Optimization for Machine Learning, MIT Press, 2011. https://arxiv.org/abs/1507.01030
• D. P. Bertsekas, A. Nedić, A. E. Ozdaglar, Convex Analysis and Optimization, Athena Scientific, 2003.

Motivation

Balanced truncation is the standard method for reducing the state dimension of a linear control system while keeping its input–output behaviour: states that are both hard to reach and hard to observe are discarded. Its justification rests on two Gramians, the solutions $P$ and $Q$ of a pair of Lyapunov equations, whose quadratic forms measure how much input energy a state costs to reach and how much output energy it produces. For linear systems the facts are classical: the image of the controllability Gramian is the reachable subspace, and the kernel of the observability Gramian is the unobservable subspace.

Bilinear systems, in which the input multiplies the state, are the simplest nonlinear extension of linear control systems. Their generalized Lyapunov equations have the same form as the Lyapunov equations of linear stochastic systems with multiplicative noise, which is why Benner and Damm treat the two classes together. Balanced truncation was extended to them using the solutions of generalized Lyapunov equations (Al-Baiyat and Bettayeb, 1993; Gray and Mesko, 1998, both as cited by Benner and Damm). Benner and Damm, in Lyapunov Equations, Energy Functionals, and Model Order Reduction of Bilinear and Stochastic Systems (SIAM J. Control Optim. 49(2), 2011, doi:10.1137/09075041X), examine which of the linear-system energy interpretations survive for bilinear systems. Their Theorem 3.1 is the structural result that does survive: the kernels of the Gramians still identify states that are unreachable from the origin and states that produce no output, so these states can be eliminated without changing the transfer behaviour.

Setting

Fix $n, m, p \ge 0$ and real matrices $A, N_1,\dots,N_m\in\mathbb R^{n\times n}$, $B\in\mathbb R^{n\times m}$, $C\in\mathbb R^{p\times n}$. The bilinear control system (3.1)–(3.2) is

$$\dot x = Ax + \sum_{j=1}^m N_j x\,u_j + Bu,\qquad y = Cx,$$

with state $x(t)\in\mathbb R^n$, input $u(t) = [u_1(t),\dots,u_m(t)]^T\in\mathbb R^m$ and output $y(t)\in\mathbb R^p$. The homogeneous system (3.5) drops the term $Bu$: $\dot x = Ax + \sum_j N_j x\,u_j$, $y = Cx$. The same $m$ indexes the bilinear terms and the columns of $B$. Write $x(t, x_0, u)$ for a solution with input $u$ and initial state $x(0) = x_0$, and $y(t,x_0,u) = Cx(t,x_0,u)$.

An admissible input is a function $u:[0,\infty)\to\mathbb R^m$ that is integrable on every interval $[0,T]$. A solution with $x(0) = x_0$ is a continuous $x:[0,\infty)\to\mathbb R^n$ satisfying

$$x(t) = x_0 + \int_0^t \Big(Ax(s) + \sum_{j=1}^m u_j(s)N_jx(s) + Bu(s)\Big)\,ds\qquad (t\ge 0).$$

The generalized Lyapunov equations (3.6) are

$$AP + PA^T + \sum_{j=1}^m N_jPN_j^T = -BB^T,\qquad A^TQ + QA + \sum_{j=1}^m N_j^TQN_j = -C^TC .$$

A solution $P$ of the first is a controllability Gramian and a solution $Q$ of the second an observability Gramian. "$P\ge 0$" means $P$ is symmetric nonnegative definite; $\mathrm{Im}\,P$ and $\mathrm{Ker}\,P$ are its image and kernel.

Formalization targets

Goal: Theorem 3.1 (p. 695)

(a) If $P\ge 0$ solves the first equation of (3.6), then for every admissible input $u$

$$x(t,0,u)\in\mathrm{Im}\,P\qquad\text{for all } t\ge 0 .$$

(b) If $Q\ge 0$ solves the second equation of (3.6) and $x_0\in\mathrm{Ker}\,Q$, then for the homogeneous system (3.5) and every admissible input $u$

$$y(t,x_0,u) = Cx(t,x_0,u) = 0\qquad\text{for all } t\ge 0 .$$

Milestones (the steps of the paper's proof, pp. 695–696)

1. $N_j^T\,\mathrm{Ker}\,P\subset\mathrm{Ker}\,P\subset\mathrm{Ker}\,B^T$.
2. $A^T\,\mathrm{Ker}\,P\subset\mathrm{Ker}\,P$.
3. $\mathrm{Im}\,P$ is invariant under the vector field of (3.1): $z\in\mathrm{Im}\,P$ implies $Az + \sum_j w_jN_jz + Bw\in\mathrm{Im}\,P$ for every $w\in\mathbb R^m$.
4. $N_j\,\mathrm{Ker}\,Q\subset\mathrm{Ker}\,Q\subset\mathrm{Ker}\,C$ and $A\,\mathrm{Ker}\,Q\subset\mathrm{Ker}\,Q$.
5. $\mathrm{Ker}\,Q$ is invariant under the vector field of (3.5).

Significance

The result. Part (a) says that a state outside $\mathrm{Im}\,P$ cannot be reached from the origin by any input, so its controllability energy $E_c$ is infinite; part (b) says that a state in $\mathrm{Ker}\,Q$ is indistinguishable from the origin at the output, so its observability energy $E_o$ is zero. Together they justify the first step of balanced truncation for bilinear systems (Remark 3.2 of the paper): the system can be restricted to $\mathrm{Im}\,P$ and projected along $\mathrm{Ker}\,Q$ without changing its input–output map, after which the Gramians are nonsingular and can be balanced. The theorem needs neither stability of $A$ nor uniqueness of the Gramians, so it applies to every nonnegative definite solution of (3.6). The rest of §3 of the paper shows that the finer, quantitative energy estimates of the linear case fail for bilinear systems; this qualitative statement is what remains exact.

Formalizing it. The theorem is proved in the paper; no machine-checked version is known. A formal proof yields a reusable invariance principle: a subspace that is invariant under every vector field $z\mapsto Az + \sum_j w_jN_jz + Bw$ is invariant under the flow of the input-driven system, for merely integrable inputs. It also yields the linear-algebra lemmas relating kernels of nonnegative definite solutions of generalized Lyapunov equations to the system matrices, which are used throughout model reduction and stochastic stability theory.

Difficulty

The algebraic part (milestones 1, 2, 4) rests on the fact that for a nonnegative definite matrix $M$, $v^TMv = 0$ forces $Mv = 0$, applied to the terms of the Lyapunov equation tested against a kernel vector. The passage from the vector field to trajectories is where care is needed. The paper's sentence "$\dot x(t)\in\mathrm{Im}\,P$ whenever $x(t)\in\mathrm{Im}\,P$, so $\mathrm{Im}\,P$ is invariant" is not by itself an argument: tangency of a vector field to a subspace implies invariance only together with a uniqueness or Gronwall-type estimate, and here the field depends on time through an input that is only integrable, so $x$ is merely absolutely continuous and $\dot x$ exists only almost everywhere. A pointwise "the derivative points into the subspace" argument therefore does not apply as stated, and the step needs an estimate on the component of $x(t)$ orthogonal to the subspace.

Formalization scope

The Lean development uses states Fin n → ℝ, matrices Matrix (Fin n) (Fin n) ℝ, the matrix–vector product *ᵥ, and indices $j = 1,\dots,m$ as Fin m. $\mathrm{Im}\,P$ and $\mathrm{Ker}\,P$ are LinearMap.range and LinearMap.ker of Matrix.toLin' P; "$P\ge 0$" is Matrix.PosSemidef, which includes symmetry. Admissible inputs are integrable on every $[0,T]$ (this contains $L^2[0,\infty[$, the input space of the paper's energy functionals). A solution is a function continuous on $[0,\infty)$ whose integrand is integrable on every $[0,t]$ and which satisfies the integral equation for all $t\ge 0$; no uniqueness of solutions is presupposed, and the theorem is asserted for every such solution.

Two reading decisions are fixed. First, part (b) is stated for every admissible input $u$, not only for $u = 0$ as printed in the theorem: the paper's proof concludes "$x(t,x_0,u)\in\mathrm{Ker}\,Q$ for all $t\ge0$, implying $y(t,x_0,u) = 0$", and the theorem's clause "i.e., $E_o(x_0) = 0$" (with $E_o$ a maximum over inputs) needs this form. The $\forall u$ statement contains the printed one. Second, the clauses "i.e., $E_c(x_0) = \infty$" and "i.e., $E_o(x_0) = 0$", which restate the conclusions through the energy functionals (3.3)–(3.4), are not formalized.

The statements must not be trivialized: assuming $P$ positive definite makes (a) empty ($\mathrm{Im}\,P = \mathbb R^n$), assuming $A$ stable or the Gramians unique adds hypotheses the paper does not have, and a solution notion requiring $C^1$ trajectories for integrable inputs would make the solution set empty. None of these is used. A local check confirms the hypotheses have content: for $A = -I_2$, $N_1 = 0$, $B = e_1$, the matrix $P = \mathrm{diag}(1/2, 0)$ is a nonnegative definite solution of (3.6) with $e_2\notin\mathrm{Im}\,P$.

A complete development needs Mathlib's facts on nonnegative definite matrices, interval integrals of vector-valued functions, and a Gronwall inequality with an integrable coefficient. The flow-invariance lemma for input-driven linear-in-state systems is reusable beyond this mission. Contributions welcome: proofs of the five milestones, and the passage from vector-field invariance to trajectory invariance.

Selected references

• P. Benner, T. Damm, Lyapunov Equations, Energy Functionals, and Model Order Reduction of Bilinear and Stochastic Systems, SIAM J. Control Optim. 49(2), 686–711, 2011. https://doi.org/10.1137/09075041X
• S. A. Al-Baiyat, M. Bettayeb, A new model reduction scheme for k-power bilinear systems, Proceedings of the 32nd IEEE Conference on Decision and Control, San Antonio, 1993, pp. 22–27 (reference [1] of Benner–Damm).
• W. S. Gray, J. Mesko, Energy functions and algebraic Gramians for bilinear systems, Preprints of the 4th IFAC Nonlinear Control Systems Design Symposium, Enschede, 1998, pp. 103–108 (reference [25] of Benner–Damm).
• B. C. Moore, Principal component analysis in linear systems: Controllability, observability, and model reduction, IEEE Trans. Automat. Control 26(1), 17–32, 1981. https://doi.org/10.1109/TAC.1981.1102568