Optimization

Motivation

Interior point methods for a linear program $\min\{c^\top x : Ax\ge b\}$ with $A\in\mathbb R^{m\times n}$ follow the central path of the logarithmic barrier $-\sum_i\log s_i$, where $s=Ax-b$ is the slack vector. Renegar's path-following analysis (1988) gives $O(\sqrt m\,L)$ iterations, and for decades this was the best bound for methods whose iterations cost a linear system solve. Vaidya's volumetric barrier $-\log\det(A^\top S^{-2}A)$ and the hybrid volumetric barriers of Vaidya and of Anstreicher (references [45] and [2] of the paper) reached $O((m\,\mathrm{rank}(A))^{1/4}L)$ iterations at the price of more expensive linear algebra. Nesterov and Nemirovski showed that a universal barrier gives $O(\sqrt n\,L)$ iterations, but that barrier cannot be evaluated efficiently.

Lee and Sidford (FOCS 2014; full version arXiv:1312.6677) obtained $\tilde O(\sqrt{\mathrm{rank}(A)}\,L)$ iterations, each costing $\tilde O(1)$ linear system solves, by following a weighted central path whose weights are recomputed from the slacks. The weights come from a weight function $g$, defined as the minimizer of a regularized D-optimal-design problem. This mission is about that weight function and the theorem (Theorem 1 of the paper) certifying its properties. The companion mission, Path-Finding Methods for Linear Programming I, formalizes the path-following framework (Theorem 5 of §IV.C) that consumes these properties.

Setting

Fix $A\in\mathbb R^{m\times n}$ with full column rank, $\mathrm{rank}(A)=n$, and $1\le n<m$. For vectors $s,w\in\mathbb R^m_{>0}$ write $S=\mathrm{diag}(s)$, $W=\mathrm{diag}(w)$, $W^\alpha=\mathrm{diag}(w_i^\alpha)$, and $A_s=S^{-1}A$. For a matrix $M$ let $\|v\|_M=\sqrt{v^\top Mv}$.

Projection matrix and slack sensitivity (Definition 2, p. 428). The projection matrix is $P_{S^{-1}A}(w)=W^{1/2}S^{-1}A\,(A^\top S^{-1}WS^{-1}A)^{-1}A^\top S^{-1}W^{1/2}$, and the slack sensitivity is

$$\gamma(s,w)=\max_{i\in[m]}\big\|W^{-1/2}\mathbb 1_i\big\|_{P_{S^{-1}A}(w)} .$$

Weight function (Definition 4, p. 428). A map $g:\mathbb R^m_{>0}\to\mathbb R^m_{>0}$ is a weight function with constants $c_1,c_\gamma,c_r$ if it is differentiable and, for every $s>0$, with $G(s)=\mathrm{diag}(g(s))$, $G'(s)$ the Jacobian of $g$ at $s$, and $\|y\|_{G(s)}=\sqrt{\sum_ig_i(s)y_i^2}$:

1. Size: $\|g(s)\|_1\le c_1$;
2. Slack sensitivity: $c_\gamma\ge1$ and $\gamma(s,g(s))\le c_\gamma$;
3. Step consistency: $c_r\ge1$ and for all $r\ge c_r$, $y\in\mathbb R^m$: $\|(I+r^{-1}G^{-1}G'S)y\|_{G(s)}\le\|y\|_{G(s)}$ and $\|y+r^{-1}G^{-1}G'Sy\|_\infty\le\|y\|_\infty+c_r\|y\|_{G(s)}$;
4. Uniformity: $\|g(s)\|_\infty\le2$.

The regularized objective (6), p. 429. For $\alpha,\beta\in\mathbb R$,

$$\hat f(s,w)=\mathbb 1^\top w-\frac1\alpha\log\det\big(A_s^\top W^\alpha A_s\big)-\beta\sum_{i\in[m]}\log w_i ,\qquad g(s)=\arg\min_{w\in\mathbb R^m_{>0}}\hat f(s,w).$$

At $\alpha=1,\beta=0$ this is the D-optimal design problem, dual to computing the John ellipsoid of the polytope $\{y:|[A(y-x)]_i|\le s_i\}$ (§V.B).

Formalization targets

Goal: Theorem 1 (Properties of Weight Function), §V.A, p. 429

With

$$\alpha=1-\Big(\log_2\frac{2m}{\mathrm{rank}(A)}\Big)^{-1},\qquad \beta=\frac{\mathrm{rank}(A)}{2m},$$

the objective $\hat f(s,\cdot)$ has a unique minimizer over $\mathbb R^m_{>0}$ for every $s>0$, and the resulting $g$ is a weight function with

$$c_1(g)=2\,\mathrm{rank}(A),\qquad c_\gamma(g)=2,\qquad c_r(g)=2\log_2\frac{2m}{\mathrm{rank}(A)} .$$

Milestones: the three bullets of Theorem 1

• Size: every minimizer $w$ of $\hat f(s,\cdot)$ satisfies $\|w\|_1\le2\,\mathrm{rank}(A)$.
• Slack sensitivity: every minimizer $w$ satisfies $\gamma(s,w)\le2$.
• Step consistency: any map $g$ selecting a minimizer at every $s>0$ is differentiable on $\mathbb R^m_{>0}$ and satisfies the two step-consistency inequalities for every $r\ge2\log_2\frac{2m}{\mathrm{rank}(A)}$.

A supporting (non-milestone) item states the existence and uniqueness of the minimizer on its own.

Significance

The result. Theorem 1 is the input that turns the weighted path-following framework into an $\tilde O(\sqrt{\mathrm{rank}(A)}\,L)$-iteration method: the framework needs $O(c_\gamma^{-1}c_r^{-3}c_1^{-1/2})$-sized steps in $t$ (p. 428), and Theorem 1 makes that $\tilde\Omega(1/\sqrt{\mathrm{rank}(A)})$. The step consistency bound is what allows the weights to be recomputed after each Newton step without losing centrality. The same construction underlies later work on Lewis-weight barriers and on fast approximate John ellipsoids and maximum flow (§VIII of the paper).

Formalizing it. The theorem is proved in the full version of the paper (arXiv:1312.6677); the FOCS extended abstract contains no proofs. No part of it has a machine-checked proof. A complete formalization would give a verified account of leverage-score calculus (sums of leverage scores equal the rank; derivatives of projection matrices), of the convexity of $w\mapsto-\log\det(A^\top W^\alpha A)$ for $\alpha\in(0,1)$, and of differentiability of an argmin via the implicit function theorem, none of which is currently packaged in Mathlib in this form.

Difficulty

Size and slack sensitivity are statements about the minimizer, which is only characterized implicitly; they require precise matrix calculus for $\log\det(A_s^\top W^\alpha A_s)$ and a comparison between the matrices $A^\top WA$ (which defines $\gamma$) and $A^\top W^\alpha A$ (which defines $g$). The specific values of $\alpha$ and $\beta$ matter here: the unregularized choice $\alpha=1$, $\beta=0$ makes the problem degenerate (p. 429).

The hard part is step consistency. The Jacobian $G'$ of an argmin is available only implicitly, as the solution of a linear system obtained by differentiating the optimality condition. A bound on $\|G'\|$ that depends on $m$ is easy to get and useless: the theorem needs the operator norm of $I+r^{-1}G^{-1}G'S$ in the $G(s)$-norm to be at most $1$ as soon as $r$ exceeds $2\log_2(2m/\mathrm{rank}(A))$, and an $\ell_\infty$ bound with only an additive $c_r\|y\|_{G(s)}$ loss.

Existence and differentiability of the minimizer are conclusions, not hypotheses. The minimization is over an open orthant on which the objective is not obviously coercive or strictly convex for $\alpha<1$, and differentiability of $g$ requires the Hessian of $\hat f$ at the minimizer to be invertible.

Formalization scope

Vectors are Fin m → ℝ, matrices Matrix (Fin m) (Fin n) ℝ; inverses are Matrix.inv, $\log\det$ is Real.log (Matrix.det …), $w_i^\alpha$ is Real.rpow, $\log_2$ is Real.logb 2, the Jacobian is fderiv ℝ g s, and $\|\cdot\|_\infty$ is Mathlib's sup norm on Fin m → ℝ.

Conventions and pinned hypotheses:

• Full column rank A.rank = n is assumed in every theorem. The paper never states it, but without it $A_s^\top W^\alpha A_s$ is singular and every formula is undefined (in Lean, Matrix.inv and Real.log would return junk $0$).
• $1\le n<m$. $\beta=\mathrm{rank}(A)/(2m)$ and $\log_2(2m/\mathrm{rank}(A))$ need $\mathrm{rank}(A)\ge1$; at $m=\mathrm{rank}(A)$ the page's $\alpha$ is $0$ and $1/\alpha$ in (6) is undefined.
• Reading of $\alpha$: the exponent $-1$ is the reciprocal of $\log_2\frac{2m}{\mathrm{rank}(A)}$, giving $\alpha\in(0,1)$.
• Size is an upper bound $\|g(s)\|_1\le c_1$ (the paper's weight function has $\|g(s)\|_1=\tfrac32\mathrm{rank}(A)$, while Theorem 1 reports $c_1=2\,\mathrm{rank}(A)$).
• The first step-consistency bullet (an operator-norm bound) is stated for every vector $y$.
• $g$ is any map $\mathbb R^m\to\mathbb R^m$ whose value at each positive $s$ minimizes $\hat f(s,\cdot)$ over $\mathbb R^m_{>0}$. Only its values on the open orthant matter. The goal also asserts that such minimizers exist and are unique, so it is not vacuous.

Ruling out trivializations: the goal does not assume $g$ to be a weight function or to be differentiable, and it does not replace $g$ by an arbitrary weight function; differentiability is a conclusion (a predicate using fderiv without it would make step consistency hold vacuously wherever $g$ fails to be differentiable).

Useful infrastructure, reusable beyond this mission: leverage scores and their sum; derivatives of $w\mapsto\log\det(A^\top WA)$ and of projection matrices; convexity of $-\log\det(A^\top W^\alpha A)$ in $w$ (related to the published ConvexOptimization.log_det_concaveOn); differentiability of the argmin of a strictly convex smooth function. Contributions of these as separate theorems are welcome, as is a proof of any single bullet of Theorem 1.

Selected references

• Y. T. Lee, A. Sidford, Path Finding Methods for Linear Programming: Solving Linear Programs in Õ(√rank) Iterations and Faster Algorithms for Maximum Flow, FOCS 2014, pp. 424–433. https://doi.org/10.1109/FOCS.2014.52
• Y. T. Lee, A. Sidford, Path Finding I: Solving Linear Programs with Õ(√rank) Linear System Solves, arXiv, 2013. https://arxiv.org/abs/1312.6677
• J. Renegar, A polynomial-time algorithm, based on Newton's method, for linear programming, Mathematical Programming 40 (1988). https://doi.org/10.1007/BF01580724

Motivation

Interior point methods solve a linear program by following a central path: a curve of minimizers of a penalized objective that trades off cost against distance from the boundary of the feasible region. The classical analysis of path following with the logarithmic barrier needs $O(\sqrt{m}\,L)$ iterations for a program with $m$ constraints, where $L$ is the bit complexity of the input (Renegar 1988). For programs with many more constraints than variables, $m$ can be far larger than the dimension $n$ or the rank of the constraint matrix, and the $\sqrt m$ factor is then the bottleneck.

Lee and Sidford (FOCS 2014) reduce the iteration count to $\tilde O(\sqrt{\mathrm{rank}(A)}\,L)$ by following a weighted central path in which each constraint carries its own positive weight, and the weights are re-computed as the algorithm moves. Their improved maximum-flow algorithm is an application of the same method.

Timeline. Karmarkar (1984) gave the first polynomial-time interior point method for linear programming. Renegar (1988) showed that path following with the logarithmic barrier needs $O(\sqrt m L)$ iterations. Nesterov and Nemirovskii (1994) showed that a universal self-concordant barrier yields $O(\sqrt n L)$ iterations, but that barrier is not known to be efficiently computable. Lee and Sidford (2014) achieved $\tilde O(\sqrt{\mathrm{rank}(A)}L)$ iterations, each reducible to $\tilde O(1)$ linear-system solves.

This mission covers the first half of that framework (§IV of the paper): the weighted central path, the weighted Newton step, and the centering theorem that shows a single step followed by re-weighting makes constant-factor progress.

Setting

Let $A\in\mathbb R^{m\times n}$, $b\in\mathbb R^m$, $c\in\mathbb R^n$, and consider the linear program

$$\min_{x\in\mathbb R^n:\ Ax\ge b} c^Tx .$$

The slack of a point $x$ is $s(x)=Ax-b$, and the interior is $S^0=\{x : Ax>b\}$, the points with all slacks strictly positive. For a path parameter $t$ and a vector of positive weights $w\in\mathbb R^m_{>0}$, the weighted penalized objective is

$$f_t(x,w)=t\,c^Tx-\sum_{i=1}^m w_i\log s(x)_i .$$

A pair $(x,w)$ is feasible if $x\in S^0$ and $w>0$.

Write $S_x=\mathrm{diag}(s(x))$, $W=\mathrm{diag}(w)$ and $\|v\|_M=\sqrt{v^TMv}$. The Newton step and the centrality are

$$\vec h_t(x,w)=\big(A^TS_x^{-1}WS_x^{-1}A\big)^{-1}\big(tc-A^TS_x^{-1}w\big),\qquad \delta_t(x,w)=\big\|\vec h_t(x,w)\big\|_{A^TS_x^{-1}WS_x^{-1}A}.$$

The matrix $A^TS_x^{-1}WS_x^{-1}A$ is the Hessian of $f_t$ in $x$, and $tc-A^TS_x^{-1}w$ is its gradient; $\delta_t(x,w)=0$ exactly when $x$ minimizes $f_t(\cdot,w)$.

For slacks $s$ and weights $w$ the projection matrix is $P_{S^{-1}A}(w)=W^{1/2}S^{-1}A(A^TS^{-1}WS^{-1}A)^{-1}A^TS^{-1}W^{1/2}$ and the slack sensitivity is

$$\gamma(s,w)=\max_{i\in[m]}\big\|W^{-1/2}\vec 1_i\big\|_{P_{S^{-1}A}(w)} .$$

A weight function (Definition 4) is a differentiable map $\vec g:\mathbb R^m_{>0}\to\mathbb R^m_{>0}$ from slacks to weights with constants $c_1$ (size, a bound on $\|\vec g(s)\|_1$), $c_\gamma\ge1$ (slack sensitivity, $\gamma(s,\vec g(s))\le c_\gamma$), $c_r\ge1$ (step consistency, two inequalities on the Jacobian $G'(s)$ of $\vec g$ that hold for every $r\ge c_r$), and uniformity $\|\vec g(s)\|_\infty\le2$.

Formalization targets

Goal: Theorem 5 (Centering with Weights), §IV.C

Let $\vec g$ be a weight function for $A$ with constants $c_1,c_\gamma,c_r$, let $x^{(old)}\in S^0$, $s^{(old)}=s(x^{(old)})$, and

$$x^{(new)}=x^{(old)}-\frac{1}{1+c_r}\,\vec h_t\big(x^{(old)},\vec g(s^{(old)})\big).$$

If $\delta_t(x^{(old)},\vec g(s^{(old)}))\le\frac{1}{100c_\gamma c_r^2}$, then $x^{(new)}\in S^0$ and

$$\delta_t\big(x^{(new)},\vec g(s^{(new)})\big)\le\Big(1-\frac{1}{4c_r}\Big)\delta_t\big(x^{(old)},\vec g(s^{(old)})\big).$$

The theorem is stated for every weight function, not for the specific one constructed in §V of the paper; that construction is the subject of a separate mission.

Milestone: Lemma 3 (Split Newton Step), §IV.B

For feasible $(x^{(old)},w^{(old)})$ and $r\ge0$, the split step $x^{(new)}=x^{(old)}-\frac1{1+r}\vec h_t$, $w^{(new)}=w^{(old)}+\frac r{1+r}W_{(old)}S_{(old)}^{-1}A\vec h_t$ satisfies, whenever $\delta_t\le\frac1{8\gamma}$,

$$\delta_t\big(x^{(new)},w^{(new)}\big)\le\frac{2}{1+r}\,\gamma\,\delta_t^2,$$

with $\gamma=\gamma(s(x^{(old)}),w^{(old)})$, and the new pair is feasible.

Milestone: Lemma 1, §IV.B

For feasible $(x,w)$ and $\alpha,t\ge0$:

$$\delta_{(1+\alpha)t}(x,w)\le(1+\alpha)\delta_t(x,w)+\alpha\sqrt{\|w\|_1}.$$

Significance

Theorem 5 is the centering half of the weighted path-following method. Combined with Lemma 1, it shows that the path parameter can be doubled, while staying close to the weighted central path, in a number of steps of the form (5) controlled by $c_\gamma$, $c_r$ and $\sqrt{c_1}$. The paper then constructs (§V, Theorem 1) a weight function with $c_1=2\,\mathrm{rank}(A)$, $c_\gamma=2$ and $c_r$ logarithmic in $m/\mathrm{rank}(A)$, which yields the $\tilde O(\sqrt{\mathrm{rank}(A)})$ iteration bound. The theorem isolates exactly which properties of a weighting scheme are needed, so it applies to any weight function satisfying Definition 4.

The FOCS extended abstract states these results without proofs; the proofs are in the arXiv full version (arXiv:1312.6677). The results are proved on paper. No machine-checked formalization of weighted path following, or of the Lee–Sidford framework, is known. A formal proof would check the constants $\frac1{100}$, $\frac1{4}$, $\frac18$ and $\frac2{1+r}$ as stated in the extended abstract, and would produce reusable Lean infrastructure for Newton steps of barrier functions with explicit matrix formulas.

Difficulty

The standard analysis of Newton's method on a self-concordant barrier gives quadratic convergence of centrality for a fixed barrier. Here the barrier changes during the step: the weights are reset to $\vec g(s(x^{(new)}))$, so the new centrality is measured with respect to a different Hessian and a different gradient. The obvious argument, analysing the step at fixed weights and then treating the re-weighting as a small perturbation, does not give a contraction factor independent of $m$: without control of how $\vec g$ reacts to changes in the slacks, the re-weighting can undo the progress of the step. The step-consistency conditions of Definition 4 are the only hypotheses that control this reaction, and they are pointwise bounds on the Jacobian of $\vec g$, while the step moves the slacks by a finite amount.

Formalization scope

Vectors are Fin n → ℝ and Fin m → ℝ, matrices Matrix (Fin m) (Fin n) ℝ, and products are Matrix.mulVec and dotProduct. $S^{-1}$ is the diagonal matrix of reciprocals, $W^{\pm1/2}$ the diagonal matrices of $\sqrt{w_i}^{\pm1}$, and $\|v\|_M=\sqrt{v^TMv}$. The Newton step and centrality are defined by the explicit formulas (3) and (4), not by derivatives of $f_t$; the centrality uses the Hessian-norm form of (4). The Jacobian $G'(s)$ is the Fréchet derivative fderiv ℝ g s, and $\|\cdot\|_\infty$ is Mathlib's sup norm.

Conventions fixed where the paper is silent:

1. Full column rank. Every theorem assumes A.rank = n. The paper uses $(A^TS_x^{-1}WS_x^{-1}A)^{-1}$ without comment; the inverse exists for positive slacks and weights exactly when $A$ has full column rank. Lean's matrix inverse is $0$ on singular matrices, which would make $\vec h_t$, $\delta_t$ and $\gamma$ vanish and every statement trivially true; the rank hypothesis rules this trivializing reading out.
2. Size as an upper bound. Definition 4's "$c_1(\vec g)=\|\vec g(s)\|_1$" is read as $\|\vec g(s)\|_1\le c_1$ for all $s>0$ (the paper's own weight function reports a $c_1$ above its $\ell_1$ norm). $c_1$ does not enter Theorem 5.
3. Operator norm. Step consistency's first bullet is written as $\|(I+r^{-1}G^{-1}G'S)y\|_{G(s)}\le\|y\|_{G(s)}$ for all $y$.
4. Lemma 3's $r$ ranges over $r\ge0$, and $\gamma(x,w)$ means $\gamma(s(x),w)$.
5. Feasibility of the new point is part of the conclusion of Lemma 3 and Theorem 5, since the page's conclusion evaluates quantities defined only on the interior.
6. Maximum over $[m]$ is a supremum over Fin m (attained for $m\ge1$, equal to $0$ for $m=0$).
7. The path parameter $t$ is unrestricted in Theorem 5 and Lemma 3, as on the page; Lemma 1 assumes $t\ge0$ as the page does.

A complete development needs basic facts about weighted norms and the projection matrix $P_{S^{-1}A}(w)$, spectral comparison of the matrices $A^TS^{-1}WS^{-1}A$ for nearby slacks and weights, and calculus for vector-valued maps on the positive orthant. The weighted-norm and projection-matrix material is reusable for any interior point analysis. Proofs of the milestones, alternative arguments, and sharper constants are welcome.

Selected references

• Y. T. Lee, A. Sidford, Path Finding Methods for Linear Programming: Solving Linear Programs in Õ(√rank) Iterations and Faster Algorithms for Maximum Flow, FOCS 2014, pp. 424–433. https://doi.org/10.1109/FOCS.2014.52
• Y. T. Lee, A. Sidford, Path Finding I: Solving Linear Programs with Õ(√rank) Linear System Solves, arXiv:1312.6677, 2013. https://arxiv.org/abs/1312.6677
• J. Renegar, A polynomial-time algorithm, based on Newton's method, for linear programming, Mathematical Programming 40, 1988, pp. 59–93. https://doi.org/10.1007/BF01580724
• N. Karmarkar, A new polynomial-time algorithm for linear programming, Combinatorica 4, 1984, pp. 373–395. https://doi.org/10.1007/BF02579150
• Y. Nesterov, A. Nemirovskii, Interior-Point Polynomial Algorithms in Convex Programming, SIAM, 1994. https://doi.org/10.1137/1.9781611970791

Motivation

A large engineering or construction project is a set of jobs with precedence constraints, and most jobs can be finished faster at a higher cost (overtime, more crews, faster equipment). Planners want to know, for every possible project duration, the cheapest way to meet it. The resulting trade-off between duration and direct cost is what management compares with overhead, penalties and market losses when it picks a schedule.

J. E. Kelley, Jr. and M. R. Walker introduced the critical-path method (CPM) in 1959, from work at du Pont and Remington Rand (Kelley and Walker 1959). Alongside the critical-path computation, they modelled each job's cost as a linear function of its duration and posed the choice of durations as a parametric linear program. They stated that its optimal value, as a function of the project duration $\lambda$, is a non-increasing, piecewise linear, convex function, which they called the project cost curve. The 1959 paper gives no proof and defers the detailed development to a separate paper (Kelley 1961). Fulkerson (1961) gave a network-flow algorithm that computes the curve. Time–cost trade-off analysis ("crashing") has been a standard part of project management since then.

Setting

A project network has events labelled $0, 1, \dots, n$ with $n \ge 1$. Event $0$ is the origin and event $n$ the terminus. A finite set $P$ of jobs is given, each an ordered pair $(i,j)$: an arrow from event $i$ to event $j$. As in the paper, labels increase along arrows ($i < j$ for every $(i,j) \in P$), the origin precedes every event, and the terminus follows every event.

For job durations $y = (y_{ij})$, the earliest event times are given by recursion (1):

$$t_0^{(0)} = 0,\qquad t_j^{(0)} = \max\,[\,y_{ij} + t_i^{(0)} \mid i<j,\ (i,j)\in P\,],\quad 1\le j\le n,$$

and $t_n^{(0)}(y)$ is the earliest project completion time.

Each job has a crash duration $d_{ij}$ and a normal duration $D_{ij}$ with $0 \le d_{ij} \le D_{ij}$, and a linear job cost $a_{ij}y_{ij} + b_{ij}$ with $a_{ij} \le 0$, $b_{ij} \ge 0$. The project (direct) cost is

$$\text{(7)}\qquad \sum_{(i,j)\in P} (a_{ij} y_{ij} + b_{ij}).$$

A schedule for $\lambda$ is a pair $(y,t)$ with

$$\text{(5)}\ d_{ij}\le y_{ij}\le D_{ij},\qquad \text{(8)}\ y_{ij}\le t_j-t_i\quad ((i,j)\in P),\qquad \text{(9)}\ t_0=0,\ t_n=\lambda.$$

Let $\Lambda$ be the set of $\lambda$ for which a schedule exists. For $\lambda \in \Lambda$ the project cost curve $C(\lambda)$ is the minimum of (7) over schedules for $\lambda$. Write $\lambda_c = t_n^{(0)}(d)$ (all jobs crashed) and $\lambda_N = t_n^{(0)}(D)$ (all jobs normal).

Formalization targets

Goal: the shape of the project cost curve (p. 165)

$$C \text{ is non-increasing on } \Lambda,\qquad C \text{ is piecewise linear on } \Lambda,\qquad C \text{ is convex on } \Lambda .$$

Piecewise linear means finitely many breakpoints $\beta_0<\dots<\beta_m$ with $\Lambda\subseteq[\beta_0,\infty)$, and affine pieces on $\Lambda\cap[\beta_k,\beta_{k+1}]$ and on $\Lambda\cap[\beta_m,\infty)$. The goal fixes no breakpoints or slopes. It asserts only the shape the paper claims, on the whole of $\Lambda$.

Milestones

1. Feasible range (p. 165, "until no further reduction in project completion time is possible"): $\Lambda = [\lambda_c, \infty)$.
2. Existence of optimal schedules (p. 165, the linear program (8), (9)): for every $\lambda\in\Lambda$ the minimum of (7) is attained.
3. All-normal solution (p. 165): $(D, t^{(0)}(D))$ is a minimum cost schedule for $\lambda = \lambda_N$.
4. $\lambda$ is the earliest completion time (p. 165, "within the limits of most interest"): for $\lambda_c\le\lambda\le\lambda_N$ some minimum cost schedule $(y,t)$ for $\lambda$ has $t_n^{(0)}(y)=\lambda$.

Significance

The cost curve is the output of CPM's cost analysis. Its convexity is what makes the paper's parametric procedure valid: jobs are expedited in order of increasing marginal cost, and the curve is traced from $\lambda_N$ down to $\lambda_c$ one linear piece at a time. Monotonicity justifies reading the curve as a trade-off. Piecewise linearity with finitely many pieces means the whole curve is determined by finitely many characteristic schedules, the vertices plotted in the paper's Fig. 3. The milestones identify the domain of the curve, show that it is well defined, and fix its right end at the all-normal solution.

These facts are classical: they follow from parametric linear programming, and Kelley (1961) and Fulkerson (1961) develop them in detail. No machine-checked proof of them is known. Prove2Me has a related result, LinearOptimization.lp_optimal_cost_convex_in_rhs (Bertsimas–Tsitsiklis, Theorem 5.1): convexity of the optimal cost of a standard-form LP in its right-hand side. It covers convexity only, for a different LP form, and says nothing about monotonicity or finitely many pieces. This mission adds a formal model of CPM's time–cost program and the full three-part shape theorem.

Difficulty

Convexity alone follows from the usual argument: a convex combination of optimal schedules for two durations is a schedule for the combined duration. Monotonicity needs the structure of the network: when $\lambda$ increases, only the constraints (8) on jobs ending at the terminus loosen, because no job leaves the terminus. The hard part is piecewise linearity with finitely many pieces. Convexity does not imply it, and a general result on value functions of linear programs has to be tied to this specific program, whose right-hand side depends on $\lambda$ only through $t_n = \lambda$. The domain is also unbounded, so the argument must show that the curve is eventually a single affine (in fact constant) piece. It cannot just produce finitely many pieces on a compact interval.

Formalization scope

Events are Fin (n + 1) with origin 0 and terminus Fin.last n, and 1 ≤ n. Jobs are a Finset of ordered pairs, with at most one job per ordered pair. The standing assumptions of pp. 161–162 are fields of ProjectNetwork: labels increase along jobs, and reachability via Relation.ReflTransGen from the origin and to the terminus. Times and durations are real. Job data are functions Fin (n+1) → Fin (n+1) → ℝ, constrained and read only on $P$. The hypotheses $0\le d_{ij}\le D_{ij}$, $a_{ij}\le 0$ and $b_{ij}\ge 0$ are fields of JobData. Recursion (1) is earliest, defined by well-founded recursion on the label. It uses a fallback value $0$ for an event without predecessors, which occurs only at the origin. The paper's $\lambda$ is written lam. Constraint (9) fixes $t_n = \lambda$ exactly, and the event times are otherwise unconstrained.

The goal takes $C : \mathbb{R}\to\mathbb{R}$ with the hypothesis that $C(\lambda)$ is the least element of the set of costs of schedules for $\lambda$, for every $\lambda \in \Lambda$. All three conclusions are stated on $\Lambda$ only. This rules out the trivializing formalizations:

• a junk-valued infimum off $\Lambda$ plays no role;
• $C$ is tied to the program, and the hypothesis on $C$ is satisfiable by milestone 2;
• piecewise linearity requires finitely many pieces that cover all of $\Lambda$;
• all three properties are claimed, not convexity alone.

The goal keeps $a_{ij}\le 0$, as the page does throughout §3, although monotonicity and convexity would hold without it.

Disclosed readings:

• Milestone 1 renders "until no further reduction in project completion time is possible" as $\Lambda=[\lambda_c,\infty)$.
• Milestone 4 reads "within the limits of most interest" as $\lambda_c\le\lambda\le\lambda_N$. It asserts that some optimal schedule has $t_n^{(0)}(y)=\lambda$. "Every" is false: when all $a_{ij}=0$, the all-crash durations are optimal for every $\lambda$.

A complete development needs:

• the existence of LP optima under a bounded objective, or a direct compactness argument on the feasible polyhedron;
• a parametric-LP or polyhedral argument for finitely many linear pieces;
• basic facts on the recursion (1).

The one-variable notion IsPiecewiseLinearOn and the facts on earliest event times can be reused in scheduling missions. Proofs of the milestones, of any of the three goal conjuncts separately, and general lemmas on parametric LP value functions are all welcome.

Not formalized: general piecewise linear convex job costs (deferred by the paper to its references [7], [8]), and the primal–dual procedure itself (a method, not a claim).

Selected references

• J. E. Kelley, Jr. and M. R. Walker, Critical-Path Planning and Scheduling, Proc. Eastern Joint IRE-AIEE-ACM Computer Conference, 1959, pp. 160–173. https://doi.org/10.1145/1460299.1460318
• J. E. Kelley, Jr., Critical-Path Planning and Scheduling: Mathematical Basis, Operations Research 9(3), 1961, pp. 296–320. https://doi.org/10.1287/opre.9.3.296
• D. R. Fulkerson, A Network Flow Computation for Project Cost Curves, Management Science 7(2), 1961, pp. 167–178. https://doi.org/10.1287/mnsc.7.2.167
• D. Bertsimas and J. N. Tsitsiklis, Introduction to Linear Optimization, Athena Scientific, 1997, §5.2 (the optimal cost as a function of the right-hand side).

Motivation

A system that processes a stream of tasks can often be configured in several ways, and the configuration affects both the cost of the current task and the cost of switching before the next one: paging schemes, replicated files, server placements. When the future is unknown, the natural worst-case yardstick is competitive analysis, introduced by Sleator and Tarjan for list update and paging (Sleator–Tarjan 1985): an on-line strategy is compared with the optimal strategy that knows the whole input in advance.

Borodin, Linial and Saks (J. ACM 1992; conference version STOC 1987) proposed metrical task systems as a single model containing all such problems, and determined the exact deterministic competitive ratio of every such system. Their theorem is the starting point of the on-line-algorithms literature on metrical task systems, the $k$-server problem (Manasse–McGeoch–Sleator 1990) and their randomized variants.

Timeline. 1985: Sleator and Tarjan introduce competitive analysis for paging and list update. 1987: Borodin, Linial and Saks prove $w(S,d)=2n-1$ for every $n$-state metrical task system (journal version 1992). 1990: Manasse, McGeoch and Sleator extend the task-system model to restricted task sets and pose the $k$-server conjecture. The randomized ratio of the uniform task system, bounded in the same paper between $H(n)$ and $2H(n)$, is the subject of the companion mission.

Setting

A task system $(S,d)$ has a finite set $S$ of $n$ states and a transition-cost matrix $d$ with $d(i,i)=0$, $d(i,j)>0$ for $i\neq j$, and the triangle inequality $d(i,j)+d(j,k)\ge d(i,k)$. It is metrical if also $d(i,j)=d(j,i)$.

A task $T$ is a vector of nonnegative processing costs $T(s)$, $s\in S$. Given a task sequence $\mathbf T=T^1\cdots T^m$ and an initial state $s_0$, a schedule is a map $\sigma:\{0,\dots,m\}\to S$ with $\sigma(0)=s_0$; task $T^i$ is processed in state $\sigma(i)$, and the cost is

$$c(\mathbf T;\sigma)=\sum_{i=1}^m d(\sigma(i-1),\sigma(i))+\sum_{i=1}^m T^i(\sigma(i)).$$

The off-line optimum $c_0(\mathbf T)$ is the minimum over all schedules. An on-line algorithm $A$ chooses $\sigma(i)$ knowing only $s_0$ and $T^1,\dots,T^i$; its cost is $c_A(\mathbf T)$. For $w>0$, $A$ is $w$-competitive if there is a constant $K_w$ with $c_A(\mathbf T)\le w\,c_0(\mathbf T)+K_w$ for every finite task sequence. The competitive ratio of $A$ is $w(A)=\inf\{w: A\text{ is }w\text{-competitive}\}$, and the competitive ratio of the task system is $w(S,d)=\inf_A w(A)$.

For the upper bound the paper also uses continuous-time schedules, in which task $T^i$ occupies the interval $[i,i+1)$ and the scheduler may change state at any real time, paying $\int_i^{i+1}T^i(\sigma(t))\,dt$ for processing. For a general (possibly asymmetric) matrix $d$, the cycle offset ratio $\psi(d)$ is the maximum over closed walks $s_0,\dots,s_k=s_0$ of $\sum_i d(s_{i-1},s_i)\big/\sum_i d(s_i,s_{i-1})$; it equals $1$ when $d$ is symmetric.

Formalization targets

Goal: Theorem 1.1

For every metrical task system $(S,d)$ with $n$ states,

$$w(S,d)=2n-1 .$$

The value depends on $n$ only, not on the distances.

Milestones

• Lemma 2.1. If $c_0(T^1\cdots T^m)\to\infty$ along an infinite task sequence $\mathbf T$, then $w(A)\ge w_{\mathbf T}(A)=\limsup_m c_A/c_0$.
• Theorem 2.2. Against the cruel taskmaster $M(\varepsilon)$, which charges $\varepsilon$ in the state the algorithm currently occupies,

$$w_{\mathbf T(\varepsilon)}(A)\ge\frac{2n-1}{1+\varepsilon/\min_{i\neq j}d(i,j)} .$$

• Lemma 3.1. Every on-line continuous-time algorithm is matched, on every task sequence, by an on-line discrete-time algorithm.
• Lemmas 6.3, 6.4, 6.2. Properties of the functions $f_k$ that drive the algorithm $A^*_d$: $f_k(s)-f_k(s')\le d(s',s)$; the identity $2\sum_{s\ne s_k}f_k(s)+f_k(s_k)=C_{k-1}+\sum_{i\le k}d(s_i,s_{i-1})$; and $f_k\le h_k$, the off-line cost at the $k$-th transition time.
• Theorem 6.1 (= Theorem 1.2). For every task system, symmetric or not, $A^*_d$ has competitive ratio at most $(2n-1)\psi(d)$.

Significance

The theorem settles the deterministic competitive ratio of the whole class of metrical task systems: the lower bound says that no deterministic on-line strategy can beat $2n-1$ on any metric, and the upper bound supplies one algorithm that achieves it on every metric. For asymmetric costs the same algorithm gives $(2n-1)\psi(d)$. The $2n-1$ lower bound is also the benchmark against which restricted models, such as paging and the $k$-server problem, measure their improvements, and the randomized question it leaves open drove much of the later work on metrical task systems.

The result was proved in 1987 and is standard; to the best of our knowledge no machine-checked proof exists. A formal development would provide a reusable model of deterministic on-line algorithms and competitiveness (on-line maps from task prefixes, additive competitiveness, infima over algorithms), an adversary construction by mutual recursion with an arbitrary algorithm, and an exact treatment of continuous-time schedules with piecewise-constant task costs. These pieces are reusable for other competitive-analysis results.

Difficulty

The lower bound is not a single bad input: the adversary is built from the algorithm it plays against, so the hard task sequence exists only as a recursion interleaved with the algorithm's choices, and the bound must hold for every deterministic on-line map, including ones that behave erratically. Obtaining the exact constant $2n-1$, rather than some $\Omega(n)$ bound, requires a sharp estimate of the off-line cost of that sequence.

The upper bound needs an algorithm defined in continuous time, whose transition times are determined by accumulated processing costs; the budgets can be zero, so transitions can be instantaneous, and a formal cost must remain well defined before one knows that only finitely many transitions occur. Relating the off-line cost function at those times to the recursively defined $f_k$ (Lemma 6.2) requires reasoning about all continuous-time off-line schedules. Finally, the goal combines both directions through infima over all on-line algorithms, and the discretization of Lemma 3.1 must be composed with the continuous-time algorithm.

Formalization scope

States form a finite type S (Fintype, DecidableEq, Nonempty); the goal is stated for all $n\ge1$, where $n=1$ gives $w(S,d)=1$. Task costs are finite nonnegative reals; the paper also allows $+\infty$ entries, which are excluded (this affects neither bound). A task sequence is T : Fin m → S → ℝ, with T i the paper's $T^{i+1}$, and a schedule is σ : Fin (m+1) → S. An on-line algorithm is a map sending $(s_0,[T^1,\dots,T^i])$ to $\sigma(i)$, so on-line behaviour is built into the type. Competitiveness is written additively, $c_A\le w\,c_0+K$, with $K$ independent of the task sequence and of $s_0$.

The competitive ratio competitiveRatio d is the real infimum of the set of all $w$ for which some on-line algorithm is $w$-competitive. It is not defined as an infimum of per-algorithm real infima: a non-competitive algorithm has $W_A=\emptyset$, whose real infimum is $0$, and that would drag $w(S,d)$ to $0$ for every system. Since the goal's value $2n-1$ is at least $1$ while the empty set's real infimum is $0$, the goal cannot hold vacuously.

Continuous-time algorithms are given as lists of $(\text{state},\text{length})$ pieces per unit interval; processing integrals are exact finite sums. The algorithm $A^*_d$ minimizes over states different from the current one, as its proof requires (the printed rule ranges over all states, and would stall); ties are left arbitrary. Its budgets may be $0$, its entry times are Option ℝ, and its cost is a sum in $[0,\infty]$, so that Theorem 6.1 itself asserts that only finitely many transitions occur. The ratio $\psi(d)$ excludes closed walks that never move, and Theorems 2.2 and 6.1 require $n\ge2$, where $\min_{i\ne j}d(i,j)$ and $\psi(d)$ are defined. Lemma 3.1 is stated comparing $A$ with $A'$ (the printed statement says "as well as $A$").

A complete development needs: the discrete model and off-line optimum (finite minimum over schedules), limsup arguments in EReal, continuous-time schedules with piecewise-constant costs, and the recursion defining $A^*_d$. Proofs of any milestone, including the purely combinatorial Lemmas 6.3 and 6.4, are welcome, as is a formal composition of Lemma 3.1 with Theorem 6.1.

Selected references

• A. Borodin, N. Linial, M. E. Saks, An optimal on-line algorithm for metrical task system, Journal of the ACM 39(4):745–763, 1992. https://doi.org/10.1145/146585.146588
• D. D. Sleator, R. E. Tarjan, Amortized efficiency of list update and paging rules, Communications of the ACM 28(2):202–208, 1985. https://doi.org/10.1145/2786.2793
• M. S. Manasse, L. A. McGeoch, D. D. Sleator, Competitive algorithms for server problems, Journal of Algorithms 11(2):208–230, 1990. https://doi.org/10.1016/0196-6774(90)90003-W

Motivation

Minimizing a smooth function over a closed convex set $\Omega\subseteq\mathbb R^n$ on which projection is cheap (a box, a ball, a simplex) is a routine subproblem in large-scale optimization: box-constrained minimization is the inner solver of augmented Lagrangian methods, and bound-constrained least squares, image restoration and density estimation all have this form. The classical projected gradient method of Goldstein and of Levitin and Polyak is simple and needs only gradients and projections, but with constant or Armijo-type step lengths it is slow.

Spectral projected gradient (SPG) methods, introduced by Birgin, Martínez and Raydan (paper), combine three ingredients: the projected gradient direction; the Barzilai–Borwein (spectral) step length $\alpha_{k+1}=\langle s_k,s_k\rangle/\langle s_k,y_k\rangle$, an inverse Rayleigh quotient of the average Hessian along the last step; and the nonmonotone line search of Grippo, Lampariello and Lucidi, which compares a trial value with the worst of the last $M$ objective values instead of the current one. The method is widely used in practice, and its analysis is the template for many later nonmonotone projected methods.

Timeline:

• 1964–1966: Goldstein; Levitin and Polyak introduce gradient projection.
• 1976: Bertsekas analyses the Armijo rule along the projection arc.
• 1986: Grippo, Lampariello and Lucidi introduce the nonmonotone line search for unconstrained problems.
• 1988: Barzilai and Borwein propose the two-point step size; Raydan (1993, 1997) proves convergence for quadratics and combines it with nonmonotone search in the unconstrained case.
• 2000: Birgin, Martínez and Raydan define SPG1 and SPG2 for convex constraints (SIAM J. Optim. 10(4)).
• 2003: the same authors publish the convergence proof that Theorem 2.1 refers to, in the inexact setting (IMA J. Numer. Anal. 23).

Setting

Let $\Omega\subseteq\mathbb R^n$ be nonempty, closed and convex, with the Euclidean inner product $\langle\cdot,\cdot\rangle$ and norm $\|\cdot\|$. Let $f$ have continuous partial derivatives on an open set $U\supseteq\Omega$ and write $g(x)=\nabla f(x)$. The orthogonal projection $P(z)$ is the unique point of $\Omega$ nearest to $z$. The scaled projected gradient is $g_t(x)=P(x-t\,g(x))-x$ for $x\in\Omega$, $t>0$. A point $\bar x$ is a constrained stationary point if $\langle g(\bar x),x-\bar x\rangle\ge0$ for all $x\in\Omega$.

The parameters are an integer $M\ge1$, reals $0<\alpha_{\min}<\alpha_{\max}$, a sufficient-decrease constant $\gamma\in(0,1)$ and safeguards $0<\sigma_1<\sigma_2<1$. Algorithm SPG2 starts from $x_0\in\Omega$ and $\alpha_0\in[\alpha_{\min},\alpha_{\max}]$ and at iteration $k=0,1,\dots$:

1. Stop test. If $\|P(x_k-g(x_k))-x_k\|=0$, stop: $x_k$ is stationary.
2. Backtracking. Set $d_k=P(x_k-\alpha_k g(x_k))-x_k$ and $\lambda=1$. While

$$f(x_k+\lambda d_k)\le\max_{0\le j\le\min\{k,M-1\}}f(x_{k-j})+\gamma\lambda\langle d_k,g(x_k)\rangle\qquad(3)$$

fails, replace $\lambda$ by any $\lambda_{\rm new}\in[\sigma_1\lambda,\sigma_2\lambda]$. When (3) holds, $\lambda_k=\lambda$ and $x_{k+1}=x_k+\lambda_kd_k$. 3. Spectral step. With $s_k=x_{k+1}-x_k$, $y_k=g(x_{k+1})-g(x_k)$, $b_k=\langle s_k,y_k\rangle$: $\alpha_{k+1}=\alpha_{\max}$ if $b_k\le0$, else $\alpha_{k+1}=\min\{\alpha_{\max},\max\{\alpha_{\min},\langle s_k,s_k\rangle/b_k\}\}$.

In Lean the projection is a function P with the predicate IsProjOnto Ω P, $g_t$ is scaledProjGrad P f t, stationarity is IsConstrainedStationary Ω f, the maximum in (3) is nonmonotoneRef f x M k, and an infinite run is IsSPG2Run Ω f P M αmin αmax γ σ₁ σ₂ x α.

Formalization targets

Goal: Theorem 2.1, accumulation points are stationary

For every infinite run $(x_k,\alpha_k)$ of SPG2 and every accumulation point $\bar x$ of $(x_k)$,

$$\langle g(\bar x),x-\bar x\rangle\ge0\qquad\text{for all }x\in\Omega.$$

The statement fixes no parameter values and assumes neither convexity of $f$ nor a bounded level set.

Milestones

• Lemma 2.1 (ii). For $\bar x\in\Omega$ and $t\in(0,\alpha_{\max}]$: $g_t(\bar x)=0$ iff $\bar x$ is a constrained stationary point.
• Lemma 2.1 (i). For $x\in\Omega$ and $t\in(0,\alpha_{\max}]$:

$$\langle g(x),g_t(x)\rangle\le-\tfrac1t\|g_t(x)\|_2^2\le-\tfrac1{\alpha_{\max}}\|g_t(x)\|_2^2.$$

• Theorem 2.1, first clause (SPG2 is well defined). At a point where Step 1 does not stop, every admissible backtracking sequence reaches a step satisfying (3). The step is stated for an arbitrary reference value $R\ge f(x)$, which covers the maximum in (3).
• Section 2, p. 4. The iterates remain in $\Omega_0=\{x\in\Omega:f(x)\le f(x_0)\}$.

Significance

Theorem 2.1 is the global convergence guarantee for SPG2. It holds without monotone decrease of $f$ and without any restriction on the spectral step beyond the safeguards. These are the two features that make the method fast in practice, and together they mean that no classical monotone projected-gradient argument applies directly. The same statement underlies the convergence claims of the SPG software (ACM TOMS Algorithm 813) and of the many methods that reuse the nonmonotone spectral framework: inexact SPG, augmented Lagrangian inner solvers, and projected BB methods for machine learning.

Status: the theorem is proved in the literature. This paper's proof reads "See [7]", a pointer to Birgin, Martínez and Raydan (2003). No Lean formalization of this theorem, of the nonmonotone Armijo analysis, or of the projected-gradient stationarity lemma is known. The mission produces a formal proof and a reusable Lean interface for projection-based first-order methods on convex sets.

Difficulty

The obvious argument for monotone descent methods is to show that $f(x_k)$ decreases, so that the total decrease is finite and the per-iteration decrease $\gamma\lambda_k|\langle d_k,g(x_k)\rangle|$ tends to zero. Here $f(x_k)$ need not decrease. Only the reference value $\max_{0\le j\le\min\{k,M-1\}}f(x_{k-j})$ is nonincreasing, and a small decrease of this maximum along the whole sequence does not by itself give a small decrease at the iterates that approach a given accumulation point $\bar x$. A second difficulty is that the accepted step lengths $\lambda_k$ may tend to zero along the subsequence, while $f$ is $C^1$ only on a neighbourhood of $\Omega$ and no Lipschitz constant for $g$ is available, so no uniform sufficient-decrease estimate holds. The spectral steps $\alpha_k$ vary within $[\alpha_{\min},\alpha_{\max}]$, so the directions $d_k$ are not a fixed function of $x_k$.

Formalization scope

• Space and data. The space is EuclideanSpace ℝ (Fin n) with inner ℝ and the 2-norm. $f$ is a total function EuclideanSpace ℝ (Fin n) → ℝ with ContDiffOn ℝ 1 f U on an open U ⊇ Ω, and $g$ is Mathlib's gradient f. The algorithm evaluates $f$ and $g$ only at points of $\Omega$.
• Iteration and trials. Iterations are indexed from $0$. The backtracking choice (2) is universally quantified: a run carries, at each iteration, a finite trial list $\lambda^{(0)}=1$, $\lambda^{(i+1)}\in[\sigma_1\lambda^{(i)},\sigma_2\lambda^{(i)}]$, in which test (3) fails at every trial but the last and holds at the last.
• Step size. $\alpha_{k+1}$ is given by Step 3 exactly.
• Accumulation point. An accumulation point is MapClusterPt x̄ atTop x.
• Excluded simplifications. A run predicate that accepts any positive step, or lets $\alpha_{k+1}$ range freely over $[\alpha_{\min},\alpha_{\max}]$, is not SPG2. Nor is a goal stating $g_t(\bar x)=0$ instead of the variational inequality, or one that adds convexity of $f$, a Lipschitz gradient or a bounded level set.
• Non-vacuity. The hypotheses of the goal are satisfiable: for $f(x)=\|x\|^2$, $\Omega=\mathbb R^n$, $M=1$, $\alpha_{\min}=1/8$, $\alpha_{\max}=1/4$, $\gamma=1/2$ and $v\ne0$, the iterates $x_k=2^{-k}v$ with $\alpha_k=1/4$ form an infinite run with accumulation point $0$.
• Infrastructure. A complete development needs the variational characterization of the projection (Mathlib has it for the iInf form: norm_eq_iInf_iff_real_inner_le_zero), continuity properties of the projection, a mean-value estimate for $C^1$ functions on segments in $\Omega$, and the nonmonotone reference-value bookkeeping. The projection lemmas and the nonmonotone bookkeeping are reusable beyond this mission, in particular for the companion mission on SPG1, and contributions of them as separate lemmas are welcome.

Selected references

• E. G. Birgin, J. M. Martínez, M. Raydan, Nonmonotone spectral projected gradient methods on convex sets, SIAM J. Optim. 10(4) (2000) 1196–1211; authors' updated version, July 2004. https://doi.org/10.1137/S1052623497330963, https://www.ime.unicamp.br/~martinez/bmr.pdf
• E. G. Birgin, J. M. Martínez, M. Raydan, Inexact spectral projected gradient methods on convex sets, IMA J. Numer. Anal. 23 (2003) 539–559. https://doi.org/10.1093/imanum/23.4.539
• J. Barzilai, J. M. Borwein, Two-point step size gradient methods, IMA J. Numer. Anal. 8 (1988) 141–148. https://doi.org/10.1093/imanum/8.1.141
• L. Grippo, F. Lampariello, S. Lucidi, A nonmonotone line search technique for Newton's method, SIAM J. Numer. Anal. 23 (1986) 707–716. https://doi.org/10.1137/0723046
• M. Raydan, The Barzilai and Borwein gradient method for the large scale unconstrained minimization problem, SIAM J. Optim. 7 (1997) 26–33. https://doi.org/10.1137/S1052623494266365
• D. P. Bertsekas, On the Goldstein–Levitin–Polyak gradient projection method, IEEE Trans. Automat. Control 21 (1976) 174–184. https://doi.org/10.1109/TAC.1976.1101194