Linear algebra

Motivation

Nonlinear programming algorithms are routinely advertised as finding a local minimum. Most of them only certify first-order conditions (a KKT point), and the natural next question is whether a given feasible point is in fact a local minimum. For smooth problems where the Hessian is nonsingular this is a second-order test, but in the degenerate case the test itself becomes a combinatorial question about a quadratic form restricted to a cone.

K. G. Murty and S. N. Kabadi (Math. Programming 39 (1987) 117–129) showed that this question is intractable already in its simplest instance: deciding whether $x = 0$ is a local minimum of a quadratic form $x^{\mathsf T}Dx$ on the nonnegative orthant is NP-complete, and so is deciding whether a square integer matrix is copositive. The paper is the standard reference for the hardness of copositivity testing, which matters for copositive programming, for second-order optimality checks in constrained optimization, and for the complexity of local search in nonconvex optimization.

Setting

Let $D$ be a real square matrix of order $n$ and $Q(x) = x^{\mathsf T}Dx$. The matrix $D$ is copositive if $Q(x) \ge 0$ for every $x \ge 0$ (coordinatewise). The paper considers the quadratic program

$$\text{(7)}\qquad \text{minimize } Q(x) \quad \text{subject to } x \ge 0,$$

and the following questions, each phrased so that "yes" is the interesting answer:

• Problem 1. Is $x = 0$ not a local minimum of (7)?
• Problem 2. Is $Q$ not bounded below on $\{x \ge 0\}$?
• Problem 3. Is there an $x \ge 0$ with $Q(x) < 0$ (is $D$ not copositive)?
• Problem 4. Given $a_0 > 0$, is there an $x \ge 0$ with $e^{\mathsf T}x = a_0$ and $Q(x) < 0$? Here $e$ is the all-ones vector.
• Problems 11, 12. With $h(u) = (u_1^2, \dots, u_n^2)\,D\,(u_1^2, \dots, u_n^2)^{\mathsf T}$, the objective of the unconstrained problem (15): is $u = 0$ not a local minimum of $h$ on $\mathbb R^n$, and is $h$ not bounded below?

The source problem is subset sum (Problem 5): given positive integers $d_0; d_1, \dots, d_n$, is there $y \in \{0,1\}^n$ with $\sum_j d_j y_j = d_0$? Let $l$ be the total number of digits in the data. The paper fixes an integer $\delta > 4\big(d_0\sum_j d_j\big)^2 n^3$ and a rational $\varepsilon$ with $0 < \varepsilon < 2^{-nl^2}$, and defines functions of $2n$ nonnegative variables $(y, s)$:

$$f_1(y,s) = \Big(\sum_j d_j y_j - d_0\Big)^2 + \delta\sum_j (y_j + s_j - 1)^2 + \sum_j y_j s_j,$$

$f_2 = f_1 + 2d_0\sum_j d_j y_j(1 - y_j)$, a homogeneous quadratic $f_4$ that agrees with $f_2$ on the set

$$P = \Big\{(y,s) : y \ge 0,\ s \ge 0,\ \sum_j (y_j + s_j) = n\Big\},$$

and $f_5 = f_4 - (\varepsilon/n^2)\big(\sum_j (y_j + s_j)\big)^2$. The function $f_5$ is a quadratic form $x^{\mathsf T}Mx$ in $x = (y, s) \in \mathbb R^{2n}$; the symmetric matrix $M$, with entries computed explicitly from $d_0, d, \delta, \varepsilon$, is the output of the reduction.

Formalization targets

Goal: the reduction is correct

For positive integer data $d_0; d_1, \dots, d_n$ and $\delta, \varepsilon$ as above, with $M$ the matrix of $f_5$, the following are equivalent:

$$\text{subset sum is solvable} \iff \text{P1}(M) \iff \text{P2}(M) \iff \text{P3}(M) \iff M \text{ not copositive} \iff \text{P4}(M, n) \iff \text{P11}(M) \iff \text{P12}(M).$$

This is the mathematical content of Theorems 1–3 and §4: a polynomially computable map from subset sum instances to matrices under which every one of these questions has the answer of the subset sum instance.

Milestones along the paper's chain

1. Problems 5 and 6 are equivalent: subset sum is solvable iff some $(y,s) \in P$ has $f_1 \le 0$.
2. Problems 6 and 7 are equivalent ($f_1$ vs. $f_2$ on $P$).
3. Problems 7 and 8 are equivalent ($f_2$ vs. $f_4$ on $P$).
4. Lemma 2: for an integer symmetric $D$ of size $L$, the minimum of $Q$ over $[0,1]^n$ is $0$ or at most $-2^{-L}$.
5. Problems 8 and 9 are equivalent: $\exists\,(y,s)\in P$ with $f_4 \le 0$ iff $\exists\,(y,s)\in P$ with $f_5 < 0$.
6. Problem 9 is a special case of Problem 4: $f_5(y,s) = x^{\mathsf T}Mx$, and Problem 9 is Problem 4 for $(M, a_0 = n)$.
7. Problems 3 and 4 are equivalent, for any $D$ and $a_0 > 0$.
8. Problems 1 and 2 are equivalent to Problem 3, for any $D$.
9. Problems 11 and 12 are equivalent to Problems 1 and 2, for any $D$.

Significance

The result. The equivalence shows that checking local optimality of a feasible point, checking boundedness of a quadratic objective on a cone, and checking copositivity are all at least as hard as subset sum, hence NP-hard. Through (15) the same holds for local minimality and boundedness of a quartic polynomial with no constraints at all. These facts are the standard justification for why nonconvex solvers settle for KKT points, and the copositivity part underlies the hardness of copositive programming.

The formalization. The paper's claims are proved, and have been cited for decades, but the proof as printed is a sketch: several steps are stated as "clearly" or "it can be verified", and one step of the proof of Theorem 1 (p. 125) is false as written. The inequality $(\delta/2)(y_j + s_j - 1)^2 + 2d_0 d_j y_j(1 - y_j) \ge 0$ for $y_j > 1$ fails for $y_j$ slightly above $1$, and the pointwise implication "$f_2 \le 0 \Rightarrow f_1 \le 0$ on $P$" has an explicit counterexample. The equivalence of Problems 6 and 7 itself survives numerical checks. A machine-checked proof of the goal settles the correctness of the reduction with the paper's own constants. No machine-checked proof of these statements exists on the platform.

Difficulty

The combinatorial direction (a subset sum solution gives a point with $f_5 < 0$) is a computation. The converse direction carries the content, in two places.

First, passing from $f_1$ to $f_2$ trades the linear penalty for a quadratic one. This is harmless on the box $0 \le y \le 1$ but not for $y_j > 1$, where $2d_0d_jy_j(1-y_j)$ is negative. The printed argument handles this coordinate by coordinate and is wrong there; a correct argument has to show that no point of $P$ with $f_2 \le 0$ exists unless a point with $f_1 \le 0$ does, which requires a global use of the size of $\delta$.

Second, passing from $f_4 \le 0$ to $f_5 < 0$ requires a quantitative gap: if $f_4 > 0$ on $P$ then $\min_P f_4 \ge \varepsilon$ with $\varepsilon$ of only polynomially many bits. Lemma 2 provides such a gap for the unit box and integer matrices, but $P$ is not the box and $f_4$ has rational coefficients, so the lemma does not apply verbatim.

The remaining equivalences (Problems 1, 2, 3, 4, 11, 12 for a fixed matrix) follow from homogeneity of the quadratic form and the substitution $x_j = u_j^2$, and are routine.

Formalization scope

The data $d_0, d_j, \delta$ are natural numbers and $\varepsilon$ is rational; all are cast to $\mathbb R$ in $f_1, \dots, f_5$ and $M$. The index $j = 1, \dots, n$ is Fin n, and the $2n$ variables of $M$ are indexed by Fin n ⊕ Fin n with $x =$ Sum.elim y s. Standing hypotheses: $d_j > 0$ and $d_0 > 0$ (the paper's "all positive integers"); $\delta > 4(d_0\sum_j d_j)^2n^3$ in $\mathbb N$; $\varepsilon > 0$ and $\varepsilon \cdot 2^{nl^2} < 1$ in $\mathbb Q$. The size $l$ counts decimal digits. Problem 1 is local minimality relative to the orthant, Problem 11 is unconstrained local minimality, both in the Euclidean topology. Lemma 2 is stated for integer symmetric $D$ (as §4 of the paper says, "as before … symmetric"), with $L$ Schrijver's encoding size, since the paper does not define "the size of $D$"; the "optimum is $0$ or $\le -2^{-L}$" is stated as a disjunction without an infimum. The milestones on $f_2, f_4, f_5$ over $P$ assume $n \ge 1$, which the paper assumes tacitly; the goal needs no such hypothesis.

Not formalized: membership in NP (Lemma 1), the polynomial-time computability of $M$ and its encoding size, the NP-completeness of subset sum (cited by the paper from Garey–Johnson), Theorem 4, and any of the words "NP-complete" or "NP-hard". The platform's complexity layer (Turing machines over bitstrings) has no subset sum problem and no encoding of rational matrices. The matrix $M$ has rational entries; a positive integer multiple of it is the integer matrix of Theorem 3, has the same answer to every question, and the rescaling is not formalized. The §3 standing assumption "D is not PSD" is not imposed; the constructed matrix can be PSD and the equivalence holds regardless.

The goal is about the explicit matrix $M$, whose entries are given in the definitions; it is not about "some matrix whose quadratic form is $f_5$", and the constants $\delta$ and $\varepsilon$ are the paper's explicit bounds, not "sufficiently large" or "for some $\varepsilon$". A formalization that quantifies existentially over the matrix or the precision would be trivially true and is ruled out.

Contributions welcome: proofs of the matrix identity and the homogeneity arguments (milestones 6–9), a proof of Lemma 2 (which needs a Cramer/Hadamard bound on basic solutions of the linear complementarity system (9)), and a correct proof of the $f_1 \leftrightarrow f_2$ step. The copositivity definitions and milestones 7–9 are reusable for any later work on copositive programming.

Selected references

• K. G. Murty and S. N. Kabadi, Some NP-complete problems in quadratic and nonlinear programming, Mathematical Programming 39 (1987) 117–129. https://doi.org/10.1007/BF02592948
• M. R. Garey and D. S. Johnson, Computers and Intractability: A Guide to the Theory of NP-Completeness, W. H. Freeman, 1979 (subset sum, problem [SP13]).
• A. Schrijver, Theory of Linear and Integer Programming, Wiley, 1986 (encoding sizes, §2.1).

Motivation

Linear systems of ordinary differential equations $\dot x = A(t)x$ are the first case of the general theory in which solutions can be described structurally rather than merely shown to exist. They arise directly (circuits, mechanical vibrations, Hill's equation for the stability of periodic motions, the Mathieu equation for parametric resonance) and indirectly, as the linearization of a nonlinear system along an equilibrium or a periodic orbit. In the second role they decide stability: the linearized stability theorems of later chapters, the Poincaré map of a periodic orbit and the stable-manifold theory all reduce questions about nonlinear flows to questions about linear systems with constant or periodic coefficients.

For periodic coefficients the central structural result is due to Gaston Floquet (Sur les équations différentielles linéaires à coefficients périodiques, Annales scientifiques de l'École Normale Supérieure 12 (1883), 47–88, doi:10.24033/asens.220), preceded by G. W. Hill's 1877 work on the lunar perigee, which studied the same equation in the scalar second-order case. Liouville's formula for the Wronski determinant goes back to Abel (1829) and Liouville (1838). This mission follows Chapter 3 of G. Teschl, Ordinary Differential Equations and Dynamical Systems (AMS Graduate Studies in Mathematics 140, 2012), in the author's preliminary version.

Setting

Fix $n \in \mathbb{N}$, an interval $I \subseteq \mathbb{R}$ and a continuous matrix function $A : I \to \mathbb{R}^{n\times n}$. A solution of the linear first-order system

$$\dot x(t) = A(t)\,x(t) \qquad (3.79)$$

on $I$ is a function $x$ with values in $\mathbb{R}^n$ that at every $t \in I$ has derivative $A(t)x(t)$ (one-sided at an endpoint of $I$). In Lean this is TeschlODE.Linear.IsSolution A I x.

The principal matrix solution $\Pi(t,t_0)$ is, for each $t_0 \in I$, the matrix function that solves the matrix initial value problem

$$\frac{d}{dt}\Pi(t,t_0) = A(t)\,\Pi(t,t_0), \qquad \Pi(t_0,t_0) = \mathbb{I}. \qquad (3.83)$$

Its columns are the solutions of (3.79) starting at the canonical basis vectors, and every solution is $x(t) = \Pi(t,t_0)x(t_0)$. In Lean the property is TeschlODE.Linear.IsPrincipalMatrixSolution A I Φ (the Lean name of $\Pi$ is Φ). The Wronski determinant of $n$ solutions $\varphi_1,\dots,\varphi_n$ is $W(t) = \det(\varphi_1(t),\dots,\varphi_n(t))$.

A system is periodic if $I = \mathbb{R}$ and $A(t+T) = A(t)$ for all $t$, for some $T > 0$ (3.117). Its monodromy matrix is $M(t_0) = \Pi(t_0+T, t_0)$ (3.119). A logarithm of a square matrix $M$ is any matrix $B$ with $\exp(B) = M$ (3.200), where $\exp$ is the matrix exponential.

Formalization targets

Goal: Floquet's theorem (Theorem 3.15)

Let $A \in C(\mathbb{R}, \mathbb{R}^{n\times n})$ be $T$-periodic with $T > 0$ and $\Pi$ its principal matrix solution. For every $t_0 \in \mathbb{R}$ there exist $Q(t_0) \in \mathbb{C}^{n\times n}$ and $P(\cdot,t_0) : \mathbb{R} \to \mathbb{C}^{n\times n}$ with

$$\Pi(t,t_0) = P(t,t_0)\,\exp\big((t-t_0)\,Q(t_0)\big), \qquad P(t+T,t_0) = P(t,t_0), \qquad P(t_0,t_0) = \mathbb{I}.$$

Milestones

In attack order:

1. Theorem 3.9 — for continuous $A$ on an interval $I$ the initial value problem $\dot x = A(t)x$, $x(t_0) = x_0$ has a unique solution, defined on all of $I$.
2. Theorem 3.10 — the solutions form an $n$-dimensional vector space, and a principal matrix solution exists with $x(t) = \Pi(t,t_0)x_0$.
3. Lemma 3.14 — for $T$-periodic $A$, $\Pi(t+T, t_0+T) = \Pi(t,t_0)$.
4. Lemma 3.11 (Liouville's formula) — $W(t) = W(t_0)\exp\big(\int_{t_0}^t \operatorname{tr} A(s)\,ds\big)$.
5. Lemma 3.34 — a complex matrix has a logarithm iff its determinant is nonzero; a real matrix whose real eigenvalues are all positive has a real logarithm; the square of an invertible real matrix has a real logarithm.
6. Theorem 3.12 (variation of constants) — the solution of $\dot x = A(t)x + g(t)$, $x(t_0) = x_0$, is $x(t) = \Pi(t,t_0)x_0 + \int_{t_0}^t \Pi(t,s)g(s)\,ds$.

Milestone 6 is not needed for the goal; it completes the linear theory of Section 3.4 on the same definitions and is the input to the perturbation results of Section 3.7.

Significance

Floquet's theorem reduces a periodic linear system to one with constant coefficients by a periodic change of variables. Its consequences in the book are the real version with doubled period (Corollary 3.16), the stability criterion in terms of Floquet multipliers (the eigenvalues of $M(t_0)$, Corollary 3.17), the reduction $y = P^{-1}x$ to $\dot y = Q y$ (Corollary 3.18), and the stability analysis of Hill's equation. Later in the book, the stability of a periodic orbit of a nonlinear system is read off from the Floquet multipliers of its linearization, via the Poincaré map.

All results of this mission have been proved for more than a century. What the mission adds is a machine-checked development: Mathlib has the matrix exponential (NormedSpace.exp, with Matrix.exp_add_of_commute), Picard–Lindelöf and Gronwall-type uniqueness for ODEs, but no principal matrix solution, no Liouville formula, no matrix logarithm and no Floquet theory. The two platform nodes named after Floquet's and Liouville's theorems are retired placeholders whose statement is True; no faithful formalization exists on the platform.

Difficulty

The obvious guess for the solution, $\exp\big(\int_{t_0}^t A(s)\,ds\big)x_0$, is wrong as soon as the values $A(t)$ do not commute, so neither Floquet's theorem nor the variation-of-constants formula can be obtained from the constant-coefficient theory by substitution. Floquet's theorem as stated is equivalent to the existence of a logarithm of the monodromy matrix: once $M(t_0) = \exp(TQ)$ is available, the periodicity of $P$ follows from Lemma 3.14. The existence of a matrix logarithm for every invertible complex matrix is the main missing piece of infrastructure; the logarithm is not unique and not continuous in the matrix, and a real invertible matrix need not have a real logarithm (a negative eigenvalue with a single Jordan block prevents it), which is why $Q(t_0)$ is complex. Liouville's formula, which gives $\det M(t_0) = \exp\big(\int_0^T \operatorname{tr} A\big) \neq 0$, needs the derivative of a determinant along a matrix solution.

Formalization scope

• Vectors are Fin n → ℝ, matrices Matrix (Fin n) (Fin n) ℝ; $A(t)x$ is Matrix.mulVec. The goal and Lemma 3.34 use Matrix (Fin n) (Fin n) ℂ for $P$, $Q$ and the complex logarithm; $\Pi(t,t_0)$ is compared with $P\exp(\cdot)$ after coercing its entries to $\mathbb{C}$.
• An interval is a set $I$ with Set.OrdConnected I; continuity is ContinuousOn A I. Derivatives are HasDerivWithinAt … I t at every $t \in I$ (one-sided at endpoints in $I$); the matrix derivative in (3.83) is taken entrywise. In Section 3.6 (Lemma 3.14 and the goal) $I = \mathbb{R}$ and $A$ is Continuous.
• The principal matrix solution enters as a hypothesis IsPrincipalMatrixSolution A Set.univ Φ. Theorem 3.10 shows such a $\Phi$ exists and Theorem 3.9 shows it is unique on $I\times I$, so this hypothesis names the principal matrix solution, not an arbitrary function.
• Periodicity is 0 < T and ∀ t, A (t + T) = A t. The periodicity of $P$ is required with the same $T$ and for all $t$, together with $P(t_0,t_0) = \mathbb{I}$. Without the periodicity of $P$ the goal would be trivial ($Q = 0$, $P = \Pi$); with it, the goal carries the full content of the theorem.
• The matrix exponential is Mathlib's NormedSpace.exp. "Logarithm" means any $B$ with NormedSpace.exp B = M; no branch is fixed. In the third claim of Lemma 3.34 the hypothesis $\det A \ne 0$ is explicit: the book states it in the paragraph containing the lemma, and without it the claim fails at $A = 0$.
• "The solutions form an $n$-dimensional vector space" is expressed through restrictions to $I$: closure under linear combinations, $n$ solutions linearly independent on $I$, and every solution a combination of them on $I$.
• Needed infrastructure, reusable beyond this mission: existence and uniqueness for linear systems on arbitrary intervals, the principal matrix solution and its cocycle property $\Pi(t,t_1)\Pi(t_1,t_0) = \Pi(t,t_0)$, the derivative of $\det$ along a matrix solution, and the matrix logarithm (via the Jordan form or via the holomorphic functional calculus). Contributions of any of these as standalone lemmas are welcome.

Selected references

• G. Teschl, Ordinary Differential Equations and Dynamical Systems, Graduate Studies in Mathematics 140, American Mathematical Society, 2012; author's preliminary version, Chapter 3. https://www.mat.univie.ac.at/~gerald/ftp/book-ode/ode.pdf (published version: doi:10.1090/gsm/140).
• G. Floquet, Sur les équations différentielles linéaires à coefficients périodiques, Annales scientifiques de l'École Normale Supérieure 12 (1883), 47–88. doi:10.24033/asens.220
• W. J. Culver, On the existence and uniqueness of the real logarithm of a matrix, Proceedings of the AMS 17 (1966), 1146–1151. doi:10.1090/S0002-9939-1966-0202740-6

Motivation

A semidefinite program (SDP) minimizes a linear objective $c^Tx$ subject to a linear matrix inequality $F(x) = F_0 + \sum_i x_i F_i \succeq 0$. When the data $F_i$ are uncertain, El Ghaoui, Oustry and Lebret (SIAM J. Optim. 9(1), 1998) proposed to optimize against the worst case over a norm-bounded family of perturbations: the robust SDP. Their Theorem 3.1 shows that, for unstructured ("full") perturbations, the robust SDP is itself an SDP in the enlarged variable $(x,\tau)$. Section 4 of the paper then asks what robustification does to the solution. Nominal SDPs are often ill-posed: the optimal set can be a whole face, and optimal points can jump under small data changes. Section 4 shows that, under explicit hypotheses, the robust problem has a unique solution with quadratic growth, which is the sense in which the paper describes robustness as a regularization of SDPs. This mission formalizes that result, Theorem 4.2.

Setting

Fix natural numbers $m, n, p, q$, matrices $F_0, \dots, F_m \in \mathbb{R}^{n\times n}$ (symmetric), $R_0, \dots, R_m \in \mathbb{R}^{q\times n}$, $L \in \mathbb{R}^{n\times p}$, and an objective vector $c \in \mathbb{R}^m$, $c \neq 0$. Write $F(x) = F_0 + \sum_{i=1}^m x_iF_i$ and $R(x) = R_0 + \sum_{i=1}^m x_iR_i$.

With full perturbations, uncertainty level $\rho = 1$ and $D = 0$ (the standing choices of §4), the robust SDP is the SDP

$$\text{minimize } c^Tx \quad\text{subject to}\quad \mathcal{F}(x,\tau) = \begin{bmatrix} F(x) - \tau LL^T & R(x)^T \\ R(x) & \tau I\end{bmatrix} \succeq 0 \tag{15}$$

in the variables $y = (x,\tau) \in \mathbb{R}^m \times \mathbb{R}$. A point is feasible if $\mathcal{F}(x,\tau) \succeq 0$ (symmetric positive semidefinite) and optimal if it is feasible and minimizes $c^Tx$ over all feasible $(x',\tau')$. The solution is the pair $(x,\tau)$.

The paper's hypotheses (§4.1):

• H1 (Slater): $\mathcal{F}(x,\tau) \succ 0$ for some $(x,\tau)$.
• H2 (inf-compactness): every sublevel set $\{(x,\tau)\ \text{feasible} : c^Tx \le M\}$ is bounded.
• H3(a): the nullspace of the pencil $\lambda R_0 + \sum_i x_iR_i$ is one and the same proper subspace $N \subsetneq \mathbb{R}^n$ for every $(\lambda,x) \neq (0,0)$.
• H3(b): for every $x$ the stacked matrix $\begin{bmatrix} L^T \\ R(x)\end{bmatrix}$ has full column rank.

For $\tau > 0$ put $G(x,\tau) = F(x) - \tau LL^T - \frac{1}{\tau}R(x)^TR(x)$, the Schur complement of the block $\tau I$ in $\mathcal{F}(x,\tau)$.

The quadratic growth condition (QGC) holds at an optimal point $y^\star = (x^\star,\tau^\star)$ if there are $\alpha, \varepsilon > 0$ with

$$c^Tx \ \ge\ c^Tx^\star + \alpha\,\|y - y^\star\|^2 \qquad \text{for every feasible } y = (x,\tau),\ \|y - y^\star\| < \varepsilon .$$

Formalization targets

Goal: Theorem 4.2 (p. 39)

Under $c \neq 0$, symmetry of the $F_i$, H1, H2, H3(a) and H3(b):

$$\bigl(\forall\, y^\star \text{ optimal for (15)}:\ \text{QGC holds at } y^\star\bigr)\quad\text{and}\quad \exists!\, y = (x,\tau) \text{ optimal for (15)} .$$

Both halves are stated; uniqueness is of the pair $(x,\tau)$, and existence is part of the claim.

Milestones, in the order the paper's proof uses them

1. §4.1 (p. 38): H3(a) implies $R(x) \neq 0$ for every $x$.
2. §4.2 (p. 39): under H3(a), every feasible $\tau$ is positive; in particular $\tau_{\mathrm{opt}} > 0$.
3. §4.2, Eq. (16): for $\tau > 0$, $\mathcal{F}(x,\tau) \succeq 0 \iff G(x,\tau) \succeq 0$.
4. Appendix A (p. 49): at every optimal $(x,\tau)$ there is a dual matrix $Z \succeq 0$, $Z \neq 0$, with $\operatorname{Tr} ZG(x,\tau) = 0$, $\operatorname{Tr} Z\,\partial G/\partial x_i = c_i$ and $\tau^2\operatorname{Tr}LL^TZ = \operatorname{Tr}R(x)^TR(x)Z$.
5. Appendix A (p. 49): H3(b) rules out $\operatorname{Tr}LL^TZ = \operatorname{Tr}R(x)^TR(x)Z = 0$ for $Z \succeq 0$, $Z \neq 0$, hence $\operatorname{Tr}R(x)^TR(x)Z > 0$.
6. Appendix A (pp. 49–50): under H3(a), with $\tau > 0$, $Z \succeq 0$ and $\operatorname{Tr}R(x)^TR(x)Z > 0$, the Hessian of the Lagrangian $c^Tx - \operatorname{Tr} Z\,G(x,\tau)$ is positive definite.

Significance

The result. Theorem 4.2 turns the robust SDP into a well-posed problem: a unique solution with quadratic growth. Quadratic growth is the property from which the paper's Hölder-stability results (Theorem 4.3, Corollaries 4.1–4.2) follow through the perturbation theory of Bonnans, Cominetti and Shapiro, and it is what justifies using the robust SDP as a regularization of ill-conditioned SDPs (§5.4). The remark after the theorem notes a geometric reading: the growth holds for every objective, so the boundary of the robust feasible set contains no facets.

Formalizing it. The theorem has a published proof (Appendix A), which relies on a second-order sufficient condition for nonlinear SDPs cited from Bonnans, Cominetti and Shapiro. There is no machine-checked proof of it or of any second-order optimality result for SDPs that we know of. A formalization provides a complete account of the dual attainment, complementarity and second-order steps for this concrete problem class, and it checks the paper's computations; one of them, the intermediate display for the second derivative in Appendix A, has a factor error in its cross term that does not affect the conclusion.

Difficulty

The feasible set of (15) is a spectrahedron, and linear objectives over spectrahedra do not in general have unique minimizers, since optimal faces can be flat. Uniqueness therefore cannot come from convexity alone. It has to come from curvature of the boundary at the optimum, and that curvature is carried only by the nonlinear term $\frac{1}{\tau}R(x)^TR(x)$ of the Schur complement, which is degenerate along some directions. Positive definiteness of the Hessian must be recovered from the structural hypotheses H3(a) and H3(b), which interact with a dual matrix $Z$ that is known only to exist. The natural first attempt is to use $\tau > 0$ and the positive semidefiniteness of $Z$ directly. That attempt fails: the second derivative is $\operatorname{Tr} Z\,\mathcal{R}^T\mathcal{R}$ for a direction-dependent matrix $\mathcal{R}$, which vanishes on the kernel of $Z$, so it is not positive without H3(a) relating the kernels of all members of the pencil. Dual attainment and complementarity for (15) also have to be established, and the local second-order bound then has to be converted into a statement about every nearby feasible point.

Formalization scope

• Representation. Data are bundled in RobustSDP.Uniqueness.SDPData m n p q; decision points are pairs y : (Fin m → ℝ) × ℝ; the coefficient Fs i, i : Fin m, is the paper's $F_{i+1}$. $\succeq 0$ and $\succ 0$ are Mathlib's Matrix.PosSemidef and Matrix.PosDef, which include symmetry, as the paper's notation does.
• Conventions fixed. §4's standing choices $D = 0$ and $\rho = 1$ are built into (15). The standing assumptions $c \neq 0$ and symmetric $F_i$ (p. 33) are explicit hypotheses. H2 is read as bounded sublevel sets of the feasible set in $(x,\tau)$; the paper's wording ("any unbounded sequence of feasible points produces an unbounded sequence of objectives") is meant in this sense, as its claim that H1 and H2 give existence of optimal points shows. H3(b)'s full column rank is injectivity of $\xi \mapsto (L^T\xi, R(x)\xi)$. The QGC uses the Euclidean norm on $\mathbb{R}^{m+1}$ in its local form, which is equivalent to the paper's $o(\|y - y_{\mathrm{opt}}\|^2)$ form. It is stated for (15) rather than for the paper's reformulation (16), with which (15) coincides near the optimum because $\tau_{\mathrm{opt}} > 0$. The auxiliary constraint $\tau \ge 0.99\,\tau_{\mathrm{opt}}$ of (16) is not formalized. $G$ uses Lean's τ⁻¹, which is $0$ at $\tau = 0$, so every statement about $G$ assumes $\tau > 0$ or $\tau \ne 0$.
• No trivializing reading. The goal cannot be satisfied by stating only uniqueness of $x$, by reading H2 as "the objective is bounded below", or by reading H3(a) as "$R(x) \neq 0$". The statement quantifies over the pair $(x,\tau)$, and both the quadratic growth and the existence and uniqueness halves are required. The hypotheses are jointly satisfiable: for example $m = 1$, $n = p = 2$, $q = 4$, $F(x) = \mathrm{diag}(3+x, 3-x)$, $L = I_2$, $R(x) = [1; x]\otimes I_2$ and $c = 1$.
• Infrastructure. A complete development needs Schur complements for positive semidefinite block matrices (available in Mathlib), strong duality with dual attainment for inequality-form SDPs under Slater's condition (ConvexOptimization.sdp_strong_duality on the platform, in another Mathlib environment), existence of minimizers on closed bounded sets, second derivatives of matrix-valued maps, and a local second-order argument for convex problems. The duality and second-order parts can be reused beyond this mission. Contributions to any milestone, or alternative proofs that avoid the general Bonnans–Cominetti–Shapiro theory, are welcome.

Selected references

• L. El Ghaoui, F. Oustry and H. Lebret, Robust Solutions to Uncertain Semidefinite Programs, SIAM J. Optim. 9(1), 33–52, 1998. https://doi.org/10.1137/S1052623496305717
• J. F. Bonnans, R. Cominetti and A. Shapiro, Sensitivity analysis of optimization problems under second order regular constraints, Math. Oper. Res. 23(4), 806–831, 1998 (the paper's reference [10]). https://doi.org/10.1287/moor.23.4.806
• A. Shapiro, First and second order analysis of nonlinear semidefinite programs, Math. Programming Ser. B 77, 301–320, 1997. https://doi.org/10.1007/BF02614439
• R. T. Rockafellar, Convex Analysis, Princeton University Press, 1970. https://doi.org/10.1515/9781400873173

Motivation

Uzan's review Varying Constants, Gravitation and Cosmology (Living Rev. Relativity 14 (2011) 2) is about testing whether the fundamental constants of physics change in space or time. Before any experiment can be read, §2.1 of the review settles which constants it even makes sense to ask about. Two claims from that section carry the whole programme:

• p. 14 (§2.1.2): "only the variation of dimensionless constants can be measured and in case such a variation is detected, it is impossible to determine, which dimensional constant is varying";
• pp. 15–17 (§2.1.1): three independent constants (e.g. Stoney's $G,e,c$ or Planck's $c,G,\hbar$) can be used to define the three mechanical units, after which "all other constants are dimensionless quantities" whose values do not depend on the units, and "any variation of constants that will leave these numbers unaffected is actually just a redefinition of units".

Both claims are informal statements of a classical result of dimensional analysis, the Buckingham $\pi$ theorem (E. Buckingham, Phys. Rev. 4 (1914) 345). This mission states them precisely and asks for formal proofs.

Setting

Fix $d$ base units $U_1,\dots,U_d$ (for mechanics $d=3$: $L,M,T$) and $n$ dimensional constants. The dimension matrix $D\in\mathbb R^{n\times d}$ says that constant $i$ has dimension $\prod_j U_j^{D_{ij}}$; real exponents are allowed, because Gaussian units give the charge the dimension $M^{1/2}L^{3/2}T^{-1}$ (Uzan, p. 15). A configuration $x\in\mathbb R^n_{>0}$ lists the numerical values of the constants in some system of units.

• A change of units is a vector $s\in\mathbb R^d_{>0}$. It acts by $(s\cdot x)_i=x_i\prod_j s_j^{D_{ij}}$.
• An observable $f:\mathbb R^n\to\mathbb R$ is unit invariant if $f(s\cdot x)=f(x)$ for all positive $x,s$.
• The dimensionless exponents are $\mathcal Z_D=\{a\in\mathbb R^n:\sum_i a_iD_{ij}=0\ \forall j\}$, and $\pi_a(x)=\prod_i x_i^{a_i}$ is the corresponding power product.
• Given $d$ chosen constants $e_1,\dots,e_d$, a system of natural units for $x$ is a positive $s$ with $(s\cdot x)_{e_j}=1$ for all $j$.
• The Planck dimension matrix $D_{\rm P}$ has rows $[c]=LT^{-1}$, $[G]=L^3M^{-1}T^{-2}$ and $[\hbar]=L^2MT^{-1}$.

Target

Goal (Buckingham $\pi$). Put $k=n-\operatorname{rank}D$. There are linearly independent $a^{(1)},\dots,a^{(k)}\in\mathcal Z_D$ such that every unit-invariant $f$ factors as

$$f(x)=F\big(\pi_{a^{(1)}}(x),\dots,\pi_{a^{(k)}}(x)\big)\qquad (x\in\mathbb R^n_{>0}).$$

Milestones, in order:

1. Dimensionless power products are unit invariant: $a\in\mathcal Z_D\Rightarrow \pi_a(s\cdot x)=\pi_a(x)$.
2. $\dim\mathcal Z_D=n-\operatorname{rank}D$.
3. If $\det(D_{e_j,k})\ne0$, every positive $x$ has exactly one system of natural units.
4. Values in natural units do not depend on the system of units you started from: $s'\cdot(r\cdot x)=s\cdot x$.
5. Planck units: $c=G=\hbar=1$ holds exactly when $s=(1/\ell_P,1/m_P,1/t_P)$, where $\ell_P=\sqrt{G\hbar/c^3}$, $m_P=\sqrt{\hbar c/G}$ and $t_P=\sqrt{G\hbar/c^5}$ (Uzan, p. 16).
6. Two positive configurations agree on every $\pi_a$ with $a\in\mathcal Z_D$ $\iff$ they differ by a change of units (Uzan, p. 17).
7. Only dimensionless variations are measurable: a unit-invariant $f$ takes equal values on configurations that agree on all dimensionless $\pi_a$ (Uzan, p. 14).

Significance

The result. Milestones 6–7 justify the review's working convention of quoting all bounds as bounds on dimensionless numbers ($\alpha_{\rm EM}$, $\mu=m_p/m_e$, $\alpha_G=Gm_p^2/\hbar c$, …). Milestone 2 counts how many independent numbers there are. Milestones 3–5 say that natural units are well defined.

Formalizing it. The mathematics is classical and the proofs are known. What remains is to formalize it: a reusable, general dimensional-analysis layer (unit actions, dimensionless monomials, natural units) that later missions on Uzan's review can build on.

Difficulty

Taking logarithms turns the multiplicative problem into linear algebra: changes of units act on $\log x$ by translations through the column space of $D$, and $\mathcal Z_D$ is the annihilator of that space. Two steps carry the content. The first is the duality step (annihilator of the annihilator) behind milestone 6 and the goal. The second is the bookkeeping with real powers (Real.rpow), which only behaves well on positive bases. The goal also needs a choice of $F$ with no regularity, defined through the quotient $\mathbb R^n/\operatorname{col}(D)$.

Formalization scope

• Lean namespace VaryingConstants; all definitions are in one definition file. Indices are Fin n and Fin d, and exponents are real.
• Positivity is an explicit hypothesis (IsPositive), not a subtype. Unit invariance constrains $f$ only on positive configurations.
• $\mathcal Z_D$ is the kernel of $a\mapsto a^{\mathsf T}D$, and the rank is Matrix.rank over $\mathbb R$.
• Natural units are indexed by a map e : Fin d → Fin n. The independence hypothesis is $\det D_e\ne0$.
• There is no trivializing reading: the goal fixes the number of $\pi$-groups to $n-\operatorname{rank}D$ and requires them to be linearly independent and dimensionless, so the identity map does not qualify.

Contributions such as general lemmas on Real.rpow products and on annihilators of matrix column spaces are reusable beyond this mission.

Selected references

• J.-P. Uzan, Varying Constants, Gravitation and Cosmology, Living Rev. Relativity 14 (2011) 2. http://www.livingreviews.org/lrr-2011-2 (§2.1, pp. 9–17).
• E. Buckingham, On physically similar systems; illustrations of the use of dimensional equations, Phys. Rev. 4 (1914) 345. https://doi.org/10.1103/PhysRev.4.345

Motivation

Linear recurrences with constant coefficients are among the first tools of enumerative combinatorics and the analysis of algorithms: the Fibonacci numbers, the number of binary strings avoiding a pattern, the running time of a divide-and-conquer loop and the number of tilings of a strip all satisfy relations of the form $c_0 a_{n+k} + c_1 a_{n+k-1} + \cdots + c_k a_n = 0$. Chapter 9 of Keller and Trotter's Applied Combinatorics (appliedcombinatorics.org) treats such relations as linear equations in an operator, the advancement operator, in direct analogy with linear differential equations with constant coefficients. The chapter's Principal Theorem says that the solution set of such an equation is a vector space whose dimension equals the order of the recurrence; everything else in the chapter (general solutions for distinct and repeated roots, and the reduction of nonhomogeneous equations to homogeneous ones) is organised around it.

This mission formalizes that theorem, in the book's model of functions on all of $\mathbb{Z}$, together with the numbered lemmas of Section 9.5 on which its outline rests.

Setting

Let $V$ be the real vector space of all functions $f : \mathbb{Z} \to \mathbb{R}$, with pointwise addition and scalar multiplication. The advancement operator $A : V \to V$ is

$$A f(n) = f(n+1) \qquad (n \in \mathbb{Z}),$$

a linear operator with $A^p f(n) = f(n+p)$ for $p \ge 0$. For a nonnegative integer $k$ and real constants $c_0, c_1, \dots, c_k$, the advancement operator polynomial is

$$p(A) = c_0 A^k + c_1 A^{k-1} + c_2 A^{k-2} + \cdots + c_k,$$

so that $p(A) f(n) = c_0 f(n+k) + c_1 f(n+k-1) + \cdots + c_k f(n)$. The solution space of the homogeneous equation $p(A) f = 0$ is

$$W = \{ f \in V : p(A) f = 0 \},$$

the kernel of $p(A)$, a linear subspace of $V$. In Lean, $A$ is AppliedComb.Recurrence.advance : Module.End ℝ (ℤ → ℝ), $p(A)$ is opPoly k c for a coefficient vector c : Fin (k + 1) → ℝ with c i multiplying $A^{k-i}$, and $W$ is solutionSpace k c := LinearMap.ker (opPoly k c). A factor $A - r$ is advance - r • 1.

Formalization targets

Goal: Theorem 9.18 (the Principal Theorem)

For a positive integer $k$ and real constants $c_0, \dots, c_k$ with $c_0 \neq 0$ and $c_k \neq 0$,

$$\dim_{\mathbb{R}} \{ f : \mathbb{Z} \to \mathbb{R} \;:\; (c_0 A^k + c_1 A^{k-1} + \cdots + c_k) f = 0 \} = k.$$

Milestones

1. Lemma 9.19. If $r \neq 0$ and $(A - r) f = 0$, then $f(n) = f(0)\, r^n$ for every $n \in \mathbb{Z}$.
2. Lemma 9.20. If $c_0, c_k \neq 0$, $g \in V$ is arbitrary and $p(A) f_0 = g$, then every solution of $p(A) f = g$ is $f = f_0 + f_1$ with $f_1 \in W$.
3. Theorem 9.21. If $r_1, \dots, r_k$ are distinct nonzero reals, every solution of $(A - r_1)(A - r_2)\cdots(A - r_k) f = 0$ has the form

$$f(n) = c_1 r_1^n + c_2 r_2^n + \cdots + c_k r_k^n \qquad (n \in \mathbb{Z}).$$

4. Lemma 9.22. If $k \ge 1$ and $r \neq 0$, then $(A - r)^k f = 0$ holds exactly when

$$f(n) = c_1 r^n + c_2 n r^n + c_3 n^2 r^n + \cdots + c_k n^{k-1} r^n \qquad (n \in \mathbb{Z})$$

for some real constants $c_1, \dots, c_k$.

Significance

The result itself. Theorem 9.18 is what justifies the standard recipe for solving a recurrence: once $k$ linearly independent solutions are found, every solution is a linear combination of them, and $k$ initial values pin down a unique solution. Theorems 9.21 and 9.22 name those $k$ solutions when the characteristic polynomial splits over $\mathbb{R}$, and Lemma 9.20 reduces nonhomogeneous recurrences, the ones arising from counting problems with a forcing term, to one particular solution plus the homogeneous space.

Formalizing it. The book states Theorem 9.18 without a full proof ("we won't prove the full result") and leaves Lemma 9.22 as an exercise, so a formalization supplies arguments the source omits. Mathlib's LinearRecurrence develops the $\mathbb{N}$-indexed, monic version (LinearRecurrence.solSpace_rank); the bi-infinite, two-sided setting of the book, in which the constant term $c_k$ must be nonzero, is not in Mathlib or on the platform as of this mission.

Difficulty

The subspace part of Theorem 9.18 is immediate; the content is the dimension count, in both directions. On $\mathbb{Z}$ a solution must extend to negative indices as well as positive ones, so an argument that works for sequences indexed by $\mathbb{N}$ does not transfer: there the constant term plays no role and the dimension is $k$ whatever the constant term, while on $\mathbb{Z}$ a vanishing $c_k$ changes the answer. The outline in the book proceeds through factorisations of $p(A)$, which over $\mathbb{R}$ need not exist (complex roots), so the goal cannot be obtained by combining Theorems 9.21 and 9.22 alone. Lemma 9.22 requires showing that the functions $n^j r^n$ solve $(A - r)^k f = 0$ and that they exhaust its solutions, for $n$ ranging over all integers.

Formalization scope

• Functions are $\mathbb{Z} \to \mathbb{R}$ (the real space $V$ of Section 9.5, p. 199), not $\mathbb{N} \to \mathbb{R}$ and not $\mathbb{Z} \to \mathbb{C}$. Powers $r^n$ with $n$ negative are integer powers (zpow), always of a nonzero base.
• Dimension in the goal is Module.rank ℝ (solutionSpace k c) = k, a cardinal equality, so it also asserts finite-dimensionality.
• The coefficient vector is c : Fin (k + 1) → ℝ; c 0 is the leading coefficient $c_0$ and c (Fin.last k) the constant term $c_k$.
• Theorem 9.21 asserts one inclusion, as on the page ("every solution has the form"); Lemma 9.22 asserts both ("the general solution"), as an iff.
• Lemma 9.22 carries the hypothesis $r \neq 0$, the standing assumption of Section 9.5.2 (p. 200), which the lemma's own sentence does not repeat.
• No O(·), "≈" or unspecified threshold occurs in these statements, so no explicit constant is instantiated.
• Ruling out the trivial reading: $W$ is defined as the kernel of the operator polynomial $p(A)$, never as the span of $k$ chosen functions, so the goal is not a statement about the rank of a spanning family.

Needed infrastructure: linear operators on the function space ℤ → ℝ, powers and products in Module.End, Module.rank/Module.finrank, and finite sums; Mathlib's LinearRecurrence may be adapted for the forward half. The definitions advance, opPoly and solutionSpace are reusable for any later mission on recurrences or on generating functions for bi-infinite sequences. Contributions of the milestone lemmas, of a proof of the goal that avoids factorisation over $\mathbb{C}$, and of a complex-valued variant are welcome.

Selected references

• Mitchel T. Keller and William T. Trotter, Applied Combinatorics, 2017 Edition, Chapter 9 (Recurrence Equations), CC BY-SA 4.0. appliedcombinatorics.org
• Mathlib, Mathlib/Algebra/LinearRecurrence.lean (ℕ-indexed linear recurrences, LinearRecurrence.solSpace_rank). Mathlib docs
• Ronald L. Graham, Donald E. Knuth and Oren Patashnik, Concrete Mathematics, 2nd ed., Addison-Wesley, 1994, Section 7.3 (solving recurrences). ISBN 978-0-201-55802-9.

Motivation

Many computations need the trace of a matrix $A$ that is never formed explicitly and can only be applied to vectors: the trace of a matrix function such as $\mathrm{trace}(A^{-1})$ or $\log\det A = \mathrm{trace}(\log A)$ in statistics and lattice QCD, the Frobenius norm $\|B\|_F^2 = \mathrm{trace}(B^TB)$ of an operator, or the number of triangles of a graph. The standard tool is Monte-Carlo estimation, introduced by M. F. Hutchinson (Hutchinson 1989): average $M$ quadratic forms $z_i^TAz_i$ over random sign vectors $z_i$. Each sample costs one matrix–vector product, uses one random bit per entry, and needs only additions and subtractions.

Before Avron and Toledo (2011), only the variance of such estimators had been analysed. A small variance does not say how many samples guarantee a given relative error with a given probability. Avron and Toledo gave the first bounds of this kind for several estimators. This mission covers the bound for Hutchinson's estimator, their Theorem 7.1.

Setting

A Rademacher random variable takes the values $+1$ and $-1$, each with probability $1/2$. Let $A \in \mathbb{R}^{n\times n}$ be symmetric positive semi-definite. Draw $M \ge 1$ independent random vectors $z_1, \ldots, z_M \in \mathbb{R}^n$ whose $Mn$ entries are independent Rademacher variables. Hutchinson's trace estimator is

$$H_M = \frac{1}{M}\sum_{i=1}^{M} z_i^TAz_i .$$

A single sample $z^TAz$ is an unbiased estimator of $\mathrm{trace}(A)$ (Lemma 2.1 of the paper, due to Hutchinson). For symmetric $A$ its variance is $2\bigl(\|A\|_F^2 - \sum_i A_{ii}^2\bigr)$, twice the squared Frobenius mass of $A$ off the diagonal.

Given $\epsilon > 0$ and $\delta \in (0,1)$, a random estimator $T$ is an $(\epsilon,\delta)$-approximator of $\mathrm{trace}(A)$ if

$$\Pr\bigl(|T - \mathrm{trace}(A)| \le \epsilon\,\mathrm{trace}(A)\bigr) \ge 1 - \delta .$$

The rank $\mathrm{rank}(A)$ is the number of nonzero eigenvalues $\lambda_j$ of $A$, counted with multiplicity.

Formalization targets

Goal: Theorem 7.1, sample bound for $H_M$

For every symmetric positive semi-definite $A$, every $0 < \epsilon \le 1/2$ and every $0 < \delta < 1$,

$$M \ge 6\epsilon^{-2}\ln\!\left(\frac{2\,\mathrm{rank}(A)}{\delta}\right) \quad\Longrightarrow\quad H_M \text{ is an } (\epsilon,\delta)\text{-approximator of } \mathrm{trace}(A).$$

The bound depends on $A$ only through its rank. It does not depend on the dimension $n$, on the condition number, or on how the trace is spread over the diagonal.

Milestones

1. Lemma 7.2 (Achlioptas 2001, Lemma 5). For a unit vector $\alpha \in \mathbb{R}^n$ and $S = \frac1M\sum_{i=1}^M(\alpha^Tz_i)^2$, for every $\epsilon > 0$,

$$\Pr(|S - 1| \ge \epsilon) \le 2\exp\!\left(-\frac{M}{2}\left(\frac{\epsilon^2}{2} - \frac{\epsilon^3}{3}\right)\right).$$

2. Per-direction bound (proof of Theorem 7.1, p. 8:11). Let $r \ge 1$ and $0 < \epsilon \le 1/2$. If $M \ge 6\epsilon^{-2}\ln(2r/\delta)$, then $\Pr(|S - 1| \ge \epsilon) \le \delta/r$.
3. From directions to the trace (proof of Theorem 7.1, p. 8:11). Write $A = U\Lambda U^T$ and $y_i = U^Tz_i$. If $|\frac1M\sum_i y_{ij}^2 - 1| \le \epsilon$ for every $j$ with $\lambda_j \ne 0$, then $|H_M - \mathrm{trace}(A)| \le \epsilon\,\mathrm{trace}(A)$. This step is deterministic.
4. Lemma 2.1 (Hutchinson). $\mathrm{E}(z^TAz) = \mathrm{trace}(A)$, and for symmetric $A$, $\mathrm{Var}(z^TAz) = 2(\|A\|_F^2 - \sum_i A_{ii}^2)$.

Significance

Theorem 7.1 gives a practitioner an explicit number of matrix–vector products after which Hutchinson's method is guaranteed to reach relative accuracy $\epsilon$ with confidence $1-\delta$. No bound was available before, although the method had been in wide use for two decades. The bound exceeds the one the same paper proves for Gaussian test vectors ($20\epsilon^{-2}\ln(2/\delta)$) by a $\ln(\mathrm{rank}(A))$ factor. The authors conjecture that this factor is not needed. Later work removed it: Roosta-Khorasani and Ascher (2015) proved a rank-free bound for the Rademacher case, and Cortinovis and Kressner (2022) extended sample bounds to indefinite matrices. Sample bounds of this type underlie variance-reduced estimators such as Hutch++ (Meyer et al. 2021).

The theorem is proved in the literature. To the platform's knowledge it has no machine-checked proof. Formalizing it produces a checked Rademacher concentration inequality for averages of squared linear forms (Achlioptas' lemma, which is also a core lemma of database-friendly Johnson–Lindenstrauss projections), and a checked spectral reduction from a quadratic-form estimator to its eigen-directions. Neither is currently in Mathlib.

Difficulty

The estimator is not an average of independent copies of a bounded variable with a small range: a single sample $z^TAz$ can deviate from $\mathrm{trace}(A)$ by an amount comparable to $\|A\|_F\sqrt n$. Chebyshev's inequality with the variance of Lemma 2.1 only gives a polynomial dependence on $1/\delta$. Hoeffding's inequality applied to $z^TAz$ directly gives a dependence on $n$ and on the size of the entries of $A$. Unlike the Gaussian case, the estimator cannot be written as a weighted sum of independent chi-squared variables, because rotating a Rademacher vector does not give another Rademacher vector. The central difficulty is Lemma 7.2: a tail bound for $(\alpha^Tz)^2$ that holds uniformly over every unit direction $\alpha$, including directions in which $\alpha^Tz$ is far from Gaussian (for $\alpha = e_1$ it is a constant).

Formalization scope

Matrices are Matrix (Fin n) (Fin n) ℝ, and positive semi-definiteness is Matrix.PosSemidef. The Rademacher law is $\tfrac12(\delta_1 + \delta_{-1})$ on ℝ. The sample space is Fin M → Fin n → ℝ with the product of $Mn$ copies of this law (hutchinsonSampleMeasure), so the law of the estimator is constructed, not assumed. $H_M(\omega) = (M:\mathbb{R})^{-1}\sum_i \omega_i\cdot(A\omega_i)$. Probabilities are Measure.real, and the $(\epsilon,\delta)$-approximator is Definition 4.1 verbatim. $\mathrm{rank}(A)$ is Matrix.rank. In Lemma 7.2 the i.i.d. copies $Q_i$ are the functions $(\alpha\cdot z_i)^2$ on this product space.

Corrections of printed statements, each recorded in the item's Formalization Note:

• Theorem 7.1 is printed without a range on $\epsilon$. Its proof needs $\frac{M}{2}(\frac{\epsilon^2}{2} - \frac{\epsilon^3}{3}) \ge \frac{M\epsilon^2}{6}$, which holds exactly when $\epsilon \le 1/2$. Without a range the printed statement is false: take $A = \frac1n\mathbf 1\mathbf 1^T$, $n = 10^4$, $\epsilon = 10$ and $\delta = 10^{-4}$. The condition admits $M = 1$, while $\Pr(H_1 > 11) \approx 9\cdot10^{-4} > \delta$. The goal and milestone 2 are therefore stated for $0 < \epsilon \le 1/2$.
• Lemma 2.1 is printed for an arbitrary $n\times n$ matrix. The variance formula fails for non-symmetric $A$: for $A = \begin{pmatrix}0&1\\0&0\end{pmatrix}$ the variance is $1$, not $2$. Symmetry is assumed for the variance part only.
• Proof of Theorem 7.1. The proof writes $\Lambda = UAU^T$ together with $y_i = U^Tz_i$. These fit together only for $A = U\Lambda U^T$, which is the convention of milestone 3.

For $A = 0$ the threshold involves $\ln 0$. Lean's Real.log 0 = 0 turns the condition into $M \ge 0$, which agrees with the paper's reading $\ln 0 = -\infty$. The conclusion then holds because $H_M = \mathrm{trace}(A) = 0$, so no extra hypothesis is added. A formalization that took the law of the samples as a hypothesis could make that hypothesis unsatisfiable and the theorem vacuous; the constructed product space rules this out.

A complete development needs:

• the product Rademacher measure and moment generating functions of Rademacher sums;
• a Chernoff bound for averages of i.i.d. bounded variables on a product space;
• the spectral theorem for real symmetric matrices, with rank equal to the number of nonzero eigenvalues;
• a finite union bound.

The Rademacher concentration results are reusable beyond this mission. Proofs of any milestone are welcome, as are alternative proofs of Lemma 7.2.

Selected references

• H. Avron and S. Toledo, Randomized algorithms for estimating the trace of an implicit symmetric positive semi-definite matrix, Journal of the ACM 58(2), Article 8, 2011. https://doi.org/10.1145/1944345.1944349
• M. F. Hutchinson, A stochastic estimator of the trace of the influence matrix for Laplacian smoothing splines, Communications in Statistics – Simulation and Computation 18(3), 1059–1076, 1989. https://doi.org/10.1080/03610918908812806
• D. Achlioptas, Database-friendly random projections, Proceedings of PODS 2001, 274–281. https://doi.org/10.1145/375551.375608
• F. Roosta-Khorasani and U. Ascher, Improved bounds on sample size for implicit matrix trace estimators, Foundations of Computational Mathematics 15, 1187–1212, 2015. https://doi.org/10.1007/s10208-014-9220-1
• A. Cortinovis and D. Kressner, On randomized trace estimates for indefinite matrices with an application to determinants, Foundations of Computational Mathematics 22, 875–903, 2022. https://doi.org/10.1007/s10208-021-09525-9
• R. A. Meyer, C. Musco, C. Musco and D. P. Woodruff, Hutch++: Optimal stochastic trace estimation, SOSA 2021, 142–155. https://doi.org/10.1137/1.9781611976496.16