LAProof.accuracy_proofs.mv_mathcomp

Matrix and Vector Operations: MathComp–List Correspondence

This file establishes the formal connection between MathComp matrix/vector operations and their list-based counterparts, with applications to floating-point arithmetic via the VCFloat framework.

Main Results

Structure

  • Utility definitions for norms (normv, normM, sum_abs) and sequence extraction (seq_of_rV).
  • Tactics ordify and case_splitP for working with ordinal indices and block-structured matrices.
  • Enumeration lemmas connecting ord_enum, index_enum, and list operations (nth_ord_enum', nth_index_enum, index_ord_enum, etc.).
  • Norm lemmas including positivity (normv_pos, normM_pos), the submultiplicative inequality (subMultNorm), and the triangle inequality (normv_triang).
  • Module F defining floating-point matrix/vector operations (F.dotprod, F.mulmx, F.FMA_mulmx, etc.) and proving structural lemmas such as row/column block decompositions.
  • Conversion lemmas between MathComp matrices and list-of-lists representations (listlist_of_mx, mx_of_listlist, list_of_cV, cV_of_list), including roundtrip identities and size invariants.

Dependencies

This file relies on:
  • preamble, common: basic setup and shared definitions
  • dotprod_model, sum_model: relational models of dot product and summation
  • dot_acc: dot product accuracy lemmas
  • float_acc_lems: elementary floating-point accuracy lemmas

From LAProof.accuracy_proofs Require Import preamble common
    dotprod_model sum_model dot_acc float_acc_lems.

From mathcomp.algebra_tactics Require Import ring.

Open Scope ring_scope.
Open Scope order_scope.

Norm and sequence definitions

sum_abs A i is the L1 row norm of row i of matrix A .

Definition sum_abs {m n} (A : 'M[R]_(m,n)) i : R :=
  \sum_j (Rabs (A i j)).

normv v is the infinity norm of a column vector v , i.e., the maximum of the absolute values of its entries.

Definition normv {m} (v : 'cV[R]_m) : R :=
  \big[maxr/0]_(i < m) Rabs (v i 0%Ri).

normM A is the infinity matrix norm of A , i.e., the maximum over rows of the L1 row norms.

Definition normM {m n} (A : 'M[R]_(m,n)) : R :=
  \big[maxr/0]_i (sum_abs A i).

seq_of_rV x converts a MathComp row vector to a plain list.

Definition seq_of_rV {T} [n] (x : 'rV[T]_n) :=
  map (x ord0) (ord_enum n).

Tactics

ordify n i replaces a variable i of type Z or nat with a corresponding ordinal i : 'I_n , introducing the coercion hypothesis Hi : i = nat_of_ord i (or its Z analogue).

Ltac ordify n i :=
  let Hi := fresh "H" i in
  let Hj := fresh "H" i in
  let j := fresh "i" in
  match type of i with ?t ⇒
    let t' := eval hnf in t in
    match t' with
    | Z ⇒
        assert (Hi : Datatypes.is_true (ssrnat.leq (S (Z.to_nat i)) n)) by lia;
        set (j := @Ordinal n (Z.to_nat i) Hi);
        assert (Hj : i = Z.of_nat (nat_of_ord j)) by (simpl; lia)
    | nat ⇒
        assert (Hi : Datatypes.is_true (ssrnat.leq (S i) n)) by lia;
        set (j := @Ordinal n i Hi);
        assert (Hj : i = nat_of_ord j) by (simpl; lia)
    end
  end;
  clearbody j; clear Hi;
  subst i;
  rename j into i.

case_splitP j destructs an ordinal j : 'I_(a + b) into its left (lshift]) and right (rshift) components, rewriting j in the goal accordingly.

Ltac case_splitP j :=
  tryif clearbody j then
    fail "case_splitP requires a variable, but got a local definition" j
  else
    tryif is_var j then idtac
    else fail "case_splitP requires a variable, but got" j;
  match type of j with 'I_(addn ?a ?b) ⇒
    let i := fresh "j" in
    let H := fresh in
    destruct (splitP j) as [i H | i H];
    [ replace j with (@lshift a b i);
      [ | apply ord_inj; simpl; lia ]
    | replace j with (@rshift a b i);
      [ | apply ord_inj; simpl; lia ] ];
    clear j H; rename i into j
  end.

Example uses of case_splitP

Proof of A *m row_mx Bl Br = row_mx (A *m Bl) (A *m Br) using case_splitP.

Local Remark mul_mx_row' [R : pzSemiRingType] m n p1 p2
    (A : 'M[R]_(m,n))
    (Bl : 'M[R]_(n,p1))
    (Br : 'M[R]_(n,p2)) :
  A ×m row_mx Bl Br = row_mx (A ×m Bl) (A ×m Br).

Alternative proof of mul_mx_row' following the MathComp style from mathcomp.algebra.matrix.

Local Remark mul_mx_row'' [R : pzSemiRingType] m n p1 p2
    (A : 'M[R]_(m, n))
    (Bl : 'M_(n, p1))
    (Br : 'M_(n, p2)) :
  A ×m row_mx Bl Br = row_mx (A ×m Bl) (A ×m Br).

Enumeration and list lemmas

seq.nth and List.nth agree for any list.

Lemma nth_List_nth : ∀ {A : Type} (d : A) (l : seq.seq A) (n : nat),
  seq.nth d l n = List.nth n l d.

The predecessor of a positive natural number is strictly smaller.

Lemma pred_lt : ∀ [n : nat], (0 < n → n.-1 < n)%nat.

The ordinal n-1 : 'I_n for a positive n .

Definition pred_ord [n : nat] (Hn : (0 < n)%nat) : 'I_n :=
  Ordinal (pred_lt Hn).

The finite enumeration of 'I_n has size n .

Lemma ordinal_enum_size : ∀ n : nat,
  size (Finite.enum (ordinal n)) = n.

ord_enum n has size n .

Lemma size_ord_enum : ∀ n, size (ord_enum n) = n.

The i -th element of index_enum (ordinal n) is i itself.

Lemma nth_index_enum : ∀ {n : nat} (x : 'I_n) y,
  seq.nth y (index_enum (ordinal n)) x = x.

The i -th element of ord_enum n is i itself.

Lemma nth_ord_enum' : ∀ n (i : 'I_n) x, seq.nth x (ord_enum n) i = i.

index_enum (ordinal n) equals ord_enum n.

Lemma index_ord_enum : ∀ n : nat,
  index_enum (ordinal n) = ord_enum n.

seq_of_rV x has size n .

Lemma size_seq_of_rV : ∀ {T} [n] x,
  size (@seq_of_rV T n x) = n.

The i -th element of seq_of_rV x is x ord0 i.

Lemma nth_seq_of_rV : ∀ {T} [n] (d : T) (x : 'rV[T]_n) (i : 'I_n),
  nth d (seq_of_rV x) i = x ord0 i.

Lemmas about the maxr operator

maxr is commutative.

Lemma maxrC : @commutative R R maxr.

maxr is associative.

Lemma maxrA : @associative R maxr.

Scalar multiplication distributes over a big op -fold when op is "linear" in the sense expressed by the hypothesis Hc .

Lemma big_mul {n : nat} (F : ordinal n → R) op a
  (Hc: ∀ i b, op (F i) b × a = op (F i × a) (b × a)):
  R0 ≤ a →
  \big[op/0]_(i0 < n) (F i0) × a = \big[op/0]_(i0 < n) (F i0 × a).

Scalar multiplication distributes over a big maxr-fold for nonnegative scalars.

Lemma big_max_mul {n : nat} (F : ordinal n → R) a :
  R0 ≤ a →
  \big[maxr/0]_(i0 < n) (F i0) × a = \big[maxr/0]_(i0 < n) (F i0 × a).

Norm lemmas

normv v is nonnegative for any vector v .

Lemma normv_pos {m} (v : 'cV[R]_m) : R0 ≤ normv v.

normM A is nonnegative for any matrix A .

Lemma normM_pos [m n] (A : 'M[R]_(m,n)) : R0 ≤ normM A.

Triangle inequality for absolute values under a finite sum.

Lemma Rabs_sum (n : nat) : ∀ (F : ordinal n → R),
  Rabs (\sum_j F j) ≤ \sum_j Rabs (F j).

Submultiplicativity: ‖A u‖_∞ ≤ ‖A‖_∞ · ‖u‖_∞.

Lemma subMultNorm m n (A : 'M[R]_(m,n)) (u : 'cV_n) :
  normv (A ×m u) ≤ normM A × normv u.

Triangle inequality for normv: ‖u + v‖_∞ ≤ ‖u‖_∞ + ‖v‖_∞.

Lemma normv_triang m (u v : 'cV_m) :
  normv (u + v) ≤ normv u + normv v.

Auxiliary definitions

An eliminator for the empty type 'I_0.

Local Definition crazy (T : Type) : 'I_0 → T.

If witnesses exist for every entry, then a matrix with those entries exists.

Lemma exists_mx : ∀ {T} [m n] (F : 'I_m → 'I_n → T → Prop),
  (∀ i j, ∃ x, F i j x) →
  ∃ A : 'M[T]_(m,n), ∀ i j, F i j (A i j).

Reversing ord_enum n yields the list of rev_ord images.

Lemma rev_ord_enum : ∀ n,
  rev (ord_enum n) = map (@rev_ord n) (ord_enum n).

Any list u is the image of nth d u composed with the ordinal enumeration of size u.
Lemma nth_ord_enum_lemma : ∀ [T] (d : T) (u : seq T),
  u = map (nth d u \o @nat_of_ord (size u)) (ord_enum (size u)).

sumR x equals the big sum over ordinals of nth R0 x i.

Lemma sumR_sum : ∀ (x : seq R),
  sumR x = \sum_(i in 'I_(size x)) nth R0 x (nat_of_ord i).

Module F: Floating-point MathComp matrix/vector operations


Module F.

This module defines floating-point analogues of the standard matrix and vector operations using MathComp's matrix type. The operations are parameterized over a floating-point type t and a NaN payload NAN.

Section WithNAN.
Context {NAN : FPCore.Nans} {t : type}.

sum x computes the floating-point sum of the values of x over all ordinal indices, accumulating in reverse order.
Definition sum [n : nat] (x : 'I_n → ftype t) : ftype t :=
  \big[BPLUS / neg_zero]_i x (rev_ord i).

dotprod x y computes the floating-point dot product of row vector x and column vector y using pairwise BMULT and BPLUS.
Definition dotprod [n : nat] (x : 'rV[ftype t]_n) (y : 'cV[ftype t]_n)
    : ftype t :=
  \big[BPLUS / pos_zero]_i (BMULT (x ord0 (rev_ord i)) (y (rev_ord i) ord0)).

FMA_dotprod x y computes the dot product of x and y using fused multiply-add (fma_dotprod) on their list representations.
Definition FMA_dotprod [n : nat] (x : 'rV[ftype t]_n) (y : 'cV[ftype t]_n)
    : ftype t :=
  fma_dotprod (seq_of_rV x) (seq_of_rV y^T).

mulmx A B is the floating-point matrix product, with each entry computed via dotprod.
Definition mulmx [m n p] (A : 'M[ftype t]_(m,n)) (B : 'M[ftype t]_(n,p)) :=
  \matrix_(i,k) dotprod (row i A) (col k B).

FMA_mulmx A B is the floating-point matrix product using FMA_dotprod for each entry.
Definition FMA_mulmx [m n p] (A : 'M[ftype t]_(m,n)) (B : 'M[ftype t]_(n,p)) :=
  \matrix_(i,k) FMA_dotprod (row i A) (col k B).

scalemx a M scales every entry of M by a using BMULT.

Definition scalemx [m n] (a : ftype t) (M : 'M[ftype t]_(m,n)) :=
  map_mx (BMULT a) M.

addmx A B adds two matrices entry-wise using BPLUS.

Definition addmx [m n] (A B : 'M[ftype t]_(m,n)) : 'M[ftype t]_(m,n) :=
  \matrix_(i,j) BPLUS (A i j) (B i j).

mulmx distributes over right block-row matrices.

Lemma mulmx_row :
  ∀ m n p1 p2
    (A : 'M[ftype t]_(m,n))
    (Bl : 'M_(n,p1))
    (Br : 'M_(n,p2)),
  mulmx A (row_mx Bl Br) = row_mx (mulmx A Bl) (mulmx A Br).

FMA_mulmx distributes over right block-row matrices.

Lemma FMA_mulmx_row :
  ∀ m n p1 p2
    (A : 'M[ftype t]_(m,n))
    (Bl : 'M_(n,p1))
    (Br : 'M_(n,p2)),
  FMA_mulmx A (row_mx Bl Br) = row_mx (FMA_mulmx A Bl) (FMA_mulmx A Br).

mulmx distributes over left block-column matrices.

Lemma mulmx_col :
  ∀ m1 m2 n p
    (Au : 'M[ftype t]_(m1,n))
    (Ad : 'M[ftype t]_(m2,n))
    (B : 'M_(n,p)),
  mulmx (col_mx Au Ad) B = col_mx (mulmx Au B) (mulmx Ad B).

FMA_mulmx distributes over left block-column matrices.

Lemma FMA_mulmx_col :
  ∀ m1 m2 n p
    (Au : 'M[ftype t]_(m1,n))
    (Ad : 'M[ftype t]_(m2,n))
    (B : 'M_(n,p)),
  FMA_mulmx (col_mx Au Ad) B = col_mx (FMA_mulmx Au B) (FMA_mulmx Ad B).

sum x equals the list-based sumF applied to the image of x over ord_enum n.
Lemma sum_sumF : ∀ [n] (x : 'I_n → ftype t),
  sum x = sumF (map x (ord_enum n)).

dotprod x y equals the list-based dotprodF applied to the list representations of x and y^T.
Lemma dotprod_dotprodF :
  ∀ [n] (x : 'rV[ftype t]_n) (y : 'cV[ftype t]_n),
  dotprod x y = dotprodF (seq_of_rV x) (seq_of_rV (trmx y)).

For 1 × n and n × 1 matrices, mulmx A B equals the constant matrix whose sole entry is dotprodF (seq_of_rV A) (seq_of_rV B^T).
Lemma mulmx_dotprodF :
  ∀ [n] (A : 'M[ftype t]_(1,n)) (B : 'M[ftype t]_(n,1)),
  mulmx A B = const_mx (dotprodF (seq_of_rV A) (seq_of_rV (trmx B))).

For 1 × n and n × 1 matrices, FMA_mulmx A B equals the constant matrix whose sole entry is fma_dotprod (seq_of_rV A) (seq_of_rV B^T).
Lemma FMA_mulmx_fma_dotprod :
  ∀ [n] (A : 'M[ftype t]_(1,n)) (B : 'M[ftype t]_(n,1)),
  FMA_mulmx A B = const_mx (fma_dotprod (seq_of_rV A) (seq_of_rV (trmx B))).

finitemx A asserts that every entry of A is a finite floating-point number.
Definition finitemx [m n] (A : 'M[ftype t]_(m,n)) : Prop :=
  ∀ i j, Binary.is_finite (A i j).

If addmx A B is finite entry-wise, then both A and B are.

Lemma finitemx_addmx_e : ∀ [m n] (A B : 'M[ftype t]_(m,n)),
  finitemx (addmx A B) → finitemx A ∧ finitemx B.

If scalemx c A is finite entry-wise, then A is.

Lemma finitemx_scalemx_e : ∀ [m n] (c : ftype t) (A : 'M[ftype t]_(m,n)),
  finitemx (scalemx c A) → finitemx A.

End WithNAN.

End F.

Conversions between MathComp matrices and list-of-lists

listlist_of_mx A converts a MathComp matrix to a list of rows, each row being a list of entries.
Definition listlist_of_mx {T} [m n : nat] (A : 'M[T]_(m,n)) : list (list T) :=
  map (fun i : 'I_m ⇒ map (A i) (ord_enum n)) (ord_enum m).

list_of_cV V converts a MathComp column vector to a plain list.

Definition list_of_cV {T} [n : nat] (V : 'cV[T]_n) : list T :=
  map (fun i ⇒ V i ord0) (ord_enum n).

mx_of_listlist rows cols mval builds a MathComp matrix from a list-of-lists mval, using d as the default element for out-of-bounds accesses.
Definition mx_of_listlist {T} {d : T} (rows cols : nat) (mval : list (list T))
    : 'M[T]_(rows, cols) :=
  \matrix_(i,j) seq.nth (d : T) (seq.nth nil mval i) j.

cV_of_list n vval builds a MathComp column vector from a list vval, using d as the default element for out-of-bounds accesses.
Definition cV_of_list {T} {d : T} (n : nat) (vval : list T) : 'cV[T]_n :=
  \matrix_(i,j) seq.nth (d : T) vval i.

matrix_cols_nat m cols asserts that every row in m has length cols.

Definition matrix_cols_nat {T} (m : list (list T)) (cols : nat) :=
  Forall (fun r ⇒ size r = cols) m.

Round-trip: converting a list-of-lists to a MathComp matrix and back recovers the original list-of-lists, provided the dimensions match.
Lemma listlist_of_mx_of_listlist :
  ∀ {t} {d} rows cols (mval : list (list (ftype t))),
  rows = Datatypes.length mval →
  matrix_cols_nat mval cols →
  listlist_of_mx (@mx_of_listlist _ d rows cols mval) = mval.

Round-trip: converting a MathComp matrix to a list-of-lists and back recovers the original matrix.
Lemma mx_of_listlist_of_mx :
  ∀ {T} {d : T} rows cols (A : 'M[T]_(rows,cols)),
  @mx_of_listlist _ d rows cols (listlist_of_mx A) = A.

Round-trip: converting a list to a column vector and back recovers the original list, provided the sizes match.
Lemma list_of_cV_of_list :
  ∀ {T} {d : T} n (vval : list T),
  size vval = n →
  list_of_cV (@cV_of_list _ d n vval) = vval.

Round-trip: converting a column vector to a list and back recovers the original vector.
Lemma cV_of_list_of_cV :
  ∀ {T} `{d : T} n (x : 'cV[T]_n),
  @cV_of_list _ d n (list_of_cV x) = x.

listlist_of_mx A has rows rows.

Lemma matrix_rows_listlist_of_mx : ∀ {T} [rows cols] (A : 'M[T]_(rows,cols)),
  size (listlist_of_mx A) = rows.

Every row of listlist_of_mx A has length cols.

Lemma matrix_cols_listlist_of_mx : ∀ {T} [rows cols] (A : 'M[T]_(rows,cols)),
  matrix_cols_nat (listlist_of_mx A) cols.

list_of_cV vval has length n.

Lemma size_list_of_cV : ∀ {T} [n] (vval : 'cV[T]_n),
  size (list_of_cV vval) = n.

The i-th element of list_of_cV vval is vval i ord0.

Lemma nth_list_of_cV :
  ∀ {T} {d : T} [n] (vval : 'cV[T]_n) (i : 'I_n),
  nth d (list_of_cV vval) (nat_of_ord i) = vval i ord0.

List-based floating-point operations

list_dotprod v1 v2 computes the dot product of v1 and v2 using fused multiply-add, accumulating from left to right with initial value 0.
Definition list_dotprod {NAN : FPCore.Nans} {t : type}
    (v1 v2 : list (ftype t)) : ftype t :=
  foldl (fun s x12 ⇒ BFMA (fst x12) (snd x12) s) (Zconst t 0) (zip v1 v2).

matrix_vector_mult m v applies the matrix m (given as a list of rows) to the vector v (a list), computing each entry via list_dotprod.
Definition matrix_vector_mult {NAN : FPCore.Nans} {t : type}
    (m : list (list (ftype t))) (v : list (ftype t)) : list (ftype t) :=
  map (fun row ⇒ list_dotprod row v) m.

list_of_cV commutes with col_mx: stacking two column vectors and converting to a list is the same as concatenating their list representations.
Lemma list_of_cV_col_mx : ∀ {T} n1 n2 (x : 'cV[T]_n1) (y : 'cV[T]_n2),
  list_of_cV (col_mx x y) = list_of_cV x ++ list_of_cV y.

Mapping a constant function yields a repeat.

Lemma map_const_len : ∀ {A B} (c : B) (al : list A),
  map (fun _ ⇒ c) al = repeat c (length al).

listlist_of_mx commutes with col_mx: stacking two matrices and converting to a list-of-lists is the same as appending their list-of-lists representations.
Lemma listlist_of_mx_col_mx :
  ∀ {T} n1 n2 m (A : 'M[T]_(n1,m)) (B : 'M[T]_(n2,m)),
  listlist_of_mx (col_mx A B) = listlist_of_mx A ++ listlist_of_mx B.

listlist_of_mx is injective.

Lemma listlist_of_mx_inj : ∀ {T} [m n] (A B : 'M[T]_(m,n)),
  listlist_of_mx A = listlist_of_mx B → A = B.

Main theorem: floating-point matrix–vector multiplication

Fmulmx_matrix_vector_mult is the central result connecting the list-based matrix_vector_mult (using list_dotprod) to the MathComp-based F.FMA_mulmx. Given:
  • mval: a list-of-lists of floating-point values with rows rows and cols columns,
  • vval: a list of cols floating-point values,
the list matrix_vector_mult mval vval equals the column vector F.FMA_mulmx (mx_of_listlist mval) (cV_of_list vval) converted back to a list.
Lemma Fmulmx_matrix_vector_mult :
  ∀ {NAN : FPCore.Nans} {t} rows cols
    (mval : list (list (ftype t)))
    (vval : list (ftype t)),
  rows = size mval →
  cols = size vval →
  matrix_cols_nat mval cols →
  matrix_vector_mult mval vval =
    list_of_cV (F.FMA_mulmx
      (@mx_of_listlist _ (Zconst t 0) rows cols mval)
      (@cV_of_list _ (Zconst t 0) cols vval)).