Classical and Modified Gram-Schmidt Orthogonalization

This file formalizes the algorithmic variants of Gram-Schmidt orthogonalization: the classical Gram-Schmidt (CGS) process operating on lists, and the modified Gram-Schmidt (MGS) process operating on matrices indexed by Fin m.

Mathlib already provides the abstract mathematical Gram-Schmidt process (gramSchmidt in Mathlib.Analysis.InnerProductSpace.GramSchmidtOrtho), which uses well-founded recursion over an arbitrary linearly ordered index type. The definitions here are complementary: they model the specific computational procedures as implemented in practice, and prove their correctness properties (orthogonality, span preservation, nonzero outputs) directly. A bridging theorem connecting classicalGS to Mathlib's gramSchmidt is a natural next step.

Origin

Translated from the Rocq/MathComp formalization in Gram_Schmidt/GS_functional_mathcomp.v. The original defines classical and modified GS over 'cV[R]_n.+1 (real column vectors) using MathComp's dotmul. This Lean version works over a general real inner product space E with @inner ℝ E _ as the inner product.

Main definitions

• projection u v: orthogonal projection of v onto the direction of u, equal to (⟪u, v⟫ / ⟪u, u⟫) • u.
• computeUk orthVecs vk: the CGS step vk - Σ projection(uᵢ, vk).
• classicalGSAux, classicalGS: the list-based classical Gram-Schmidt process.
• modifiedGSAux', modifiedGSAux, matrixConstruction, modifiedGS: the matrix-based modified Gram-Schmidt process.

Main results

• classicalGS_orthogonal: the classical GS output is pairwise orthogonal and preserves successive spans (gram_schmidt).
• nonzeroResultVecHelper: the main inductive lemma — first k CGS output vectors are nonzero, pairwise orthogonal, and span the same subspace as the first k inputs.
• orthoLinInd: pairwise orthogonal nonzero vectors are linearly independent.

Admitted / sorry status

The following items use sorry, matching Admitted in the original Rocq source:

• nonzeroVec: nonzero vector has nonzero squared inner product
• projCreatesOrtho: residual a - projection b a is orthogonal to b
• orthoProj: projection onto orthogonal direction is zero
• cgsAuxUk: characterization of the k-th CGS auxiliary output vector
• mgsLemmaOne: MGS matrix entry characterization
• mgsLemmaM: MGS output characterization
• gsEquiv: CGS and MGS equivalence (commented out / Admitted in Rocq)

References

• [Trefethen and Bau, Numerical Linear Algebra][trefethen1997], Chapter 8
• [Golub and Van Loan, Matrix Computations][golub2013], Section 5.2.8

Projection

noncomputable def NumLib.LinearAlgebra.GramSchmidt.projection {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (u v : E) : E

Orthogonal projection of v onto the direction of u: projection u v = (⟪u, v⟫ / ⟪u, u⟫) • u. This is the vector projection (not the subspace projection); it is the component of v in the direction of u.

Equations

• NumLib.LinearAlgebra.GramSchmidt.projection u v = (inner ℝ u v / inner ℝ u u) • u

Classical Gram-Schmidt

noncomputable def NumLib.LinearAlgebra.GramSchmidt.computeUk {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (orthVecs : List E) (vk : E) : E

Compute the k-th orthogonal vector: vk - Σ projection(uᵢ, vk). This is one step of the classical Gram-Schmidt process.

Equations

• NumLib.LinearAlgebra.GramSchmidt.computeUk orthVecs vk = vk - (List.map (fun (u : E) => NumLib.LinearAlgebra.GramSchmidt.projection u vk) orthVecs).sum

def NumLib.LinearAlgebra.GramSchmidt.classicalGSAux {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] : List E → List E → List E

Auxiliary classical GS recursion. Processes remaining input vectors, appending orthogonalized results to the accumulator.

Equations

• NumLib.LinearAlgebra.GramSchmidt.classicalGSAux [] x✝ = x✝
• NumLib.LinearAlgebra.GramSchmidt.classicalGSAux (v :: vs) x✝ = NumLib.LinearAlgebra.GramSchmidt.classicalGSAux vs (x✝ ++ [NumLib.LinearAlgebra.GramSchmidt.computeUk x✝ v])

noncomputable def NumLib.LinearAlgebra.GramSchmidt.classicalGS {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVecs : List E) : List E

The classical Gram-Schmidt process: given a list of vectors, produces a list of pairwise orthogonal vectors with the same successive spans.

Equations

• NumLib.LinearAlgebra.GramSchmidt.classicalGS initVecs = NumLib.LinearAlgebra.GramSchmidt.classicalGSAux initVecs []

Modified Gram-Schmidt

The modified Gram-Schmidt process uses matrices of vectors (Fin m → Fin m → E). It is numerically more stable than the classical variant in floating-point arithmetic, though mathematically equivalent in exact arithmetic.

def NumLib.LinearAlgebra.GramSchmidt.replaceItemMx {E : Type u_1} {m : ℕ} (A : Fin m → Fin m → E) (item : E) (row col : Fin m) : Fin m → Fin m → E

Replace a single entry in a matrix of vectors.

Equations

• NumLib.LinearAlgebra.GramSchmidt.replaceItemMx A item row col i j = if i = row ∧ j = col then item else A i j

def NumLib.LinearAlgebra.GramSchmidt.diagElements {E : Type u_1} {m : ℕ} (A : Fin m → Fin m → E) : Fin m → E

Extract diagonal elements of a matrix as a function on Fin m.

Equations

• NumLib.LinearAlgebra.GramSchmidt.diagElements A i = A i i

def NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux' {E : Type u_1} {m : ℕ} (mx : Fin m → Fin m → E) : ℕ → (row : Fin m) → Fin m → Fin m → E

Inner loop of modified GS: given a matrix, a row index, and an iteration count, update entries by subtracting projections.

Equations

• NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux' mx 0 x✝ = mx
• NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux' mx j.succ x✝ = NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux' mx j x✝

def NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux {E : Type u_1} {m : ℕ} : ℕ → (Fin m → Fin m → E) → Fin m → Fin m → E

Outer loop of modified GS: iterate through rows.

Equations

• NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux 0 x✝ = x✝
• NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux i.succ x✝ = NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux i x✝

def NumLib.LinearAlgebra.GramSchmidt.matrixConstruction {E : Type u_1} [NormedAddCommGroup E] {m : ℕ} (initVec : Fin m → E) : Fin m → Fin m → E

Construct the MGS matrix from input vectors by running the full modified GS process.

Equations

• NumLib.LinearAlgebra.GramSchmidt.matrixConstruction initVec = NumLib.LinearAlgebra.GramSchmidt.modifiedGSAux m fun (i j : Fin m) => if ↑i = 0 then initVec j else 0

def NumLib.LinearAlgebra.GramSchmidt.modifiedGS {E : Type u_1} [NormedAddCommGroup E] {m : ℕ} (initVec : Fin m → E) : Fin m → E

The modified Gram-Schmidt process: returns diagonal of the constructed matrix.

Equations

• NumLib.LinearAlgebra.GramSchmidt.modifiedGS initVec = NumLib.LinearAlgebra.GramSchmidt.diagElements (NumLib.LinearAlgebra.GramSchmidt.matrixConstruction initVec)

noncomputable def NumLib.LinearAlgebra.GramSchmidt.normalizedModifiedGS {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] {m : ℕ} (initVec : Fin m → E) : Fin m → E

Normalized modified GS: normalize each output vector to unit length.

Equations

• One or more equations did not get rendered due to their size.

Helper lemmas

Algebraic lemmas about sums and projections that support the main theorems. Where Mathlib provides a direct equivalent (e.g. inner_sub_right, inner_add_right, inner_zero_right, real_inner_comm), we use it inline rather than wrapping it.

theorem NumLib.LinearAlgebra.GramSchmidt.greaterZeroSum {ι : Type u_2} [Fintype ι] (P : ι → Prop) [DecidablePred P] (f : ι → ℝ) (hf : ∀ (i : ι), P i → f i ≥ 0) : ∑ i : ι with P i, f i ≥ 0

If all summands indexed by a predicate P are ≥ 0, the sum is ≥ 0.

theorem NumLib.LinearAlgebra.GramSchmidt.zeroSum {ι : Type u_2} [Fintype ι] (P : ι → Prop) [DecidablePred P] (f : ι → ℝ) (hf : ∀ (i : ι), P i → f i = 0) : ∑ i : ι with P i, f i = 0

If all summands indexed by a predicate P are 0, the sum is 0.

theorem NumLib.LinearAlgebra.GramSchmidt.nonzeroVec {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a : E) (ha : a ≠ 0) : inner ℝ a a ≠ 0

A nonzero vector has nonzero inner product with itself. Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.sumEquiv1 {m : ℕ} (a b c : Fin m → ℝ) : ∑ i : Fin m, (a i * b i - a i * c i) = ∑ i : Fin m, a i * b i - ∑ i : Fin m, a i * c i

Σ (aᵢ * bᵢ - aᵢ * cᵢ) = Σ aᵢ * bᵢ - Σ aᵢ * cᵢ.

theorem NumLib.LinearAlgebra.GramSchmidt.sumEquiv2 {m : ℕ} (a b c : Fin m → ℝ) : ∑ i : Fin m, (a i * b i + a i * c i) = ∑ i : Fin m, a i * b i + ∑ i : Fin m, a i * c i

Σ (aᵢ * bᵢ + aᵢ * cᵢ) = Σ aᵢ * bᵢ + Σ aᵢ * cᵢ.

theorem NumLib.LinearAlgebra.GramSchmidt.projectionDistrib {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a b c : E) : projection c (a - b) = projection c a - projection c b

Projection distributes over subtraction in the second argument.

theorem NumLib.LinearAlgebra.GramSchmidt.projCreatesOrtho {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a b : E) (hb : b ≠ 0) : inner ℝ b (a - projection b a) = 0

If b ≠ 0, then a - projection b a is orthogonal to b: ⟪b, a - projection b a⟫ = 0. Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.orthoProj {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a b x : E) (hab : inner ℝ a b = 0) : projection a (projection b x) = 0

If ⟪a, b⟫ = 0, then projection a (projection b x) = 0. Admitted in the original Rocq source.

Projection orthogonality

theorem NumLib.LinearAlgebra.GramSchmidt.projsOrtho' {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a : E) (vecList : List E) (k : ℕ) (hk : k ≤ vecList.length) (hNonzero : ∀ (i : Fin k), vecList[↑i] ≠ 0) (hOrtho : ∀ (i j : Fin k), i ≠ j → inner ℝ vecList[↑i] vecList[↑j] = 0) (x : Fin k) : inner ℝ (a - ∑ i : Fin k, projection vecList[↑i] a) vecList[↑x] = 0

Given pairwise orthogonal nonzero vectors, the residual a - Σ projection(vᵢ, a) is orthogonal to each vector vₓ in the collection. This is the key step in showing that each new Gram-Schmidt output is orthogonal to all previous outputs.

Linear independence and span lemmas

theorem NumLib.LinearAlgebra.GramSchmidt.takeConcat {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a : E) (X : List E) : List.take X.length (X ++ [a]) = X

List.take X.length (X ++ [a]) = X.

theorem NumLib.LinearAlgebra.GramSchmidt.getConcat {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a : E) (X : List E) : (X ++ [a])[X.length] = a

(X ++ [a])[X.length] = a.

theorem NumLib.LinearAlgebra.GramSchmidt.takeConcatLt {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (a : E) (X : List E) (i : ℕ) (hi : i < X.length) : List.take i (X ++ [a]) = List.take i X

For i < |X|, (X ++ [a]).take i = X.take i.

theorem NumLib.LinearAlgebra.GramSchmidt.orthoLinInd {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (vecList : List E) (hOrtho : ∀ (i j : ℕ) (hi : i < vecList.length) (hj : j < vecList.length), i ≠ j → inner ℝ vecList[i] vecList[j] = 0) (hNonzero : ∀ (i : ℕ) (hi : i < vecList.length), vecList[i] ≠ 0) : LinearIndependent ℝ fun (i : Fin vecList.length) => vecList[↑i]

Pairwise orthogonal nonzero vectors are linearly independent. This is a standard result: take the inner product of any dependency relation with a single vector to isolate its coefficient.

Main CGS theorems

theorem NumLib.LinearAlgebra.GramSchmidt.cgsAuxSize {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (vs rs : List E) : (classicalGSAux vs rs).length = vs.length + rs.length

The auxiliary GS function produces output of length |vs| + |rs|.

theorem NumLib.LinearAlgebra.GramSchmidt.initResultVecSameLength {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) : (classicalGS initVec).length = initVec.length

The classical GS output has the same length as the input.

theorem NumLib.LinearAlgebra.GramSchmidt.cgsAuxUk {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) (iter : ℕ) (rs : List E) (hiter : iter = initVec.length - rs.length) (k : ℕ) (hk : k < iter) : let resultVec := classicalGSAux initVec rs; resultVec[k] = initVec[k] - ∑ i : Fin k, projection resultVec[↑i] initVec[k]

Characterization of the k-th CGS output vector as the input vector minus the sum of projections onto all previous output vectors. Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.cgsUk {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) (k : ℕ) (hk : k < initVec.length) : let resultVec := classicalGS initVec; resultVec[k] = initVec[k] - ∑ i : Fin k, projection resultVec[↑i] initVec[k]

Direct characterization for classicalGS: the k-th output vector equals the k-th input vector minus the sum of projections onto all previous output vectors.

theorem NumLib.LinearAlgebra.GramSchmidt.nonzeroResultVecHelper {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) (hfree : LinearIndependent ℝ fun (i : Fin initVec.length) => initVec[↑i]) (k : ℕ) (hk : k ≤ initVec.length) : let resultVec := classicalGS initVec; (∀ (i : Fin k), resultVec[↑i] ≠ 0) ∧ (∀ (i j : Fin k), i ≠ j → inner ℝ resultVec[↑i] resultVec[↑j] = 0) ∧ Submodule.span ℝ (Set.range fun (i : Fin k) => resultVec[↑i]) = Submodule.span ℝ (Set.range fun (i : Fin k) => initVec[↑i])

Main inductive lemma. For k ≤ |initVec|, the first k CGS output vectors are (1) nonzero, (2) pairwise orthogonal, and (3) span the same subspace as the first k input vectors. The proof follows the Rocq original by induction on k.

theorem NumLib.LinearAlgebra.GramSchmidt.nonzeroResultVec {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) (hfree : LinearIndependent ℝ fun (i : Fin initVec.length) => initVec[↑i]) (k : ℕ) (hk : k ≤ initVec.length) : let resultVec := classicalGS initVec; (∀ (i j : Fin k), i ≠ j → inner ℝ resultVec[↑i] resultVec[↑j] = 0) ∧ Submodule.span ℝ (Set.range fun (i : Fin k) => resultVec[↑i]) = Submodule.span ℝ (Set.range fun (i : Fin k) => initVec[↑i])

Corollary: orthogonality and span preservation for the first k CGS output vectors.

theorem NumLib.LinearAlgebra.GramSchmidt.classicalGS_orthogonal {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] (initVec : List E) (hfree : LinearIndependent ℝ fun (i : Fin initVec.length) => initVec[↑i]) : ∃ (resultVec : List E), resultVec.length = initVec.length ∧ (∀ (i : ℕ) (hi : i ≤ initVec.length), Submodule.span ℝ (Set.range fun (j : Fin i) => (classicalGS initVec)[↑j]) = Submodule.span ℝ (Set.range fun (j : Fin i) => initVec[↑j])) ∧ ∀ (i j : ℕ) (hi : i < (classicalGS initVec).length) (hj : j < (classicalGS initVec).length), i ≠ j → inner ℝ (classicalGS initVec)[i] (classicalGS initVec)[j] = 0

The classical Gram-Schmidt theorem. Given linearly independent input vectors, classicalGS produces a list of pairwise orthogonal vectors with the same successive spans.

Modified GS partial results

These are attempts toward MGS characterization and the CGS/MGS equivalence proof. All are Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.mgsLemmaOne {n : ℕ} {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] {m : ℕ} (initVec : Fin (m + 1) → E) (hNonzero : ∀ (i : Fin (m + 1)), initVec i ≠ 0) (hmn : m ≤ n) (k : Fin (m + 1)) (hOrtho : ∀ (i j : Fin (m + 1)), ↑i < ↑k → ↑j < ↑k → i ≠ j → inner ℝ (matrixConstruction initVec i i) (matrixConstruction initVec j j) = 0) (h : Fin (m + 1)) : matrixConstruction initVec h k = initVec k - ∑ i : Fin (m + 1) with ↑i < ↑h, projection (matrixConstruction initVec i i) (initVec k)

MGS matrix entry characterization. Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.mgsAux {n : ℕ} {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] {m : ℕ} (initVec : Fin (m + 1) → E) (hNonzero : ∀ (i : Fin (m + 1)), initVec i ≠ 0) (hmn : m ≤ n) (k : Fin (m + 1)) (hOrtho : ∀ (i j : Fin (m + 1)), ↑i < ↑k → ↑j < ↑k → i ≠ j → inner ℝ (matrixConstruction initVec i i) (matrixConstruction initVec j j) = 0) : matrixConstruction initVec k k = initVec k - ∑ i : Fin (m + 1) with ↑i < ↑k, projection (matrixConstruction initVec i i) (initVec k)

MGS diagonal characterization: the diagonal entry M[k,k] equals the input vector minus the sum of projections onto all previous diagonal entries.

theorem NumLib.LinearAlgebra.GramSchmidt.mgsLemmaM {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] {m : ℕ} (initVec : Fin (m + 1) → E) (hNonzero : ∀ (i : Fin (m + 1)), initVec i ≠ 0) (k : Fin (m + 1)) : modifiedGS initVec k = initVec k - ∑ i : Fin (m + 1) with ↑i < ↑k, projection (modifiedGS initVec i) (initVec k) ∧ ∀ (i j : Fin (m + 1)), i ≠ j → inner ℝ (modifiedGS initVec i) (modifiedGS initVec j) = 0

MGS output characterization and pairwise orthogonality. Admitted in the original Rocq source.

theorem NumLib.LinearAlgebra.GramSchmidt.gsEquiv {n : ℕ} {E : Type u_1} [NormedAddCommGroup E] [InnerProductSpace ℝ E] {m : ℕ} (initVec : Fin (m + 1) → E) (hNonzero : ∀ (i : Fin (m + 1)), initVec i ≠ 0) (hmn : m ≤ n) (k : Fin (m + 1)) : (classicalGS (List.ofFn initVec))[↑k] = modifiedGS initVec k

CGS and MGS produce the same result in exact arithmetic. Admitted in the original Rocq source (commented out).

Auxiliary nat lemmas

theorem NumLib.LinearAlgebra.GramSchmidt.ltIN (j m : ℕ) : m - j < m + 1

m - j < m + 1.

theorem NumLib.LinearAlgebra.GramSchmidt.rowGt0 (h m : ℕ) (hh : h < m + 1) : h - 1 < m + 1

h < m + 1 → h - 1 < m + 1.