LAProof.accuracy_proofs.common

Common Definitions and Lemmas for Floating-Point Accuracy Proofs

This file provides foundational definitions and lemmas used throughout the LAProof library. It establishes the core vocabulary for reasoning about floating-point rounding, error bounds, and accumulation of rounding errors in numerical computations.

The main concepts defined here are:

Floating-point type: We use VCFloat's concept of a floating-point type, that contains a precision (mantissa size), max-exponent, and some properties of them; all based on Flocq's underlying constructions. Examples are Tdouble (64-bit IEEE double precision) and Tsingle (32-bit), and any other (legal) combination of fprec and femax that one might need.
Floating-point rounding: The function rounded captures round-to-nearest-even (RNE) rounding in radix-2 floating point, parameterized by a floating-point type t that fixes the precision fprec t and exponent range femax t.
Special floating-point values: neg_zero and pos_zero represent the IEEE 754 signed zero values for a given type t , useful when reasoning about sign-sensitive floating-point operations.
Zero testing: iszero is a boolean predicate on floating-point values that returns true exactly when the value is an IEEE 754 zero (of either sign). The lemma iszeroR_iszeroF connects this structural test to the real-number interpretation FT2R x = 0.
Counting nonzeros: nnzF (resp. nnzR ) counts the number of nonzero elements in a list of floating-point (resp. real) values. Several lemmas formalize the behavior in the case when the nonzero count is zero.
Default relative error bound default_rel:

default_rel = (1/2) * 2^(-(fprec t) + 1)

This is the unit roundoff (machine epsilon) for type t . It satisfies default_rel > 0 and numerous inequalities involving 1 + default_rel and its powers that are needed in error analysis.
Default absolute error bound default_abs:

default_abs = (1/2) * 2^(3 - femax t - fprec t)

This bounds the absolute error introduced when rounding a subnormal result. It satisfies default_abs <= default_rel .
Relative error accumulation factor g n:

g n = (1 + default_rel)^n - 1

This bounds the accumulated relative rounding error after n floating-point operations. Key properties include g_pos, le_g_Sn (monotonicity), and the recurrence one_plus_d_mul_g, which expresses how one additional rounding step advances the bound. See also: Kellison et al., "LAProof: A Library of Formal Proofs of Accuracy and Correctness of Linear Algebra Programs", 2023, equation (4), where it is called "h".
Mixed absolute/relative error accumulation factor g1 n1 n2:

g1 n1 n2 = INR n1 * default_abs * (1 + g n2)

This bounds the accumulated absolute rounding error when n1 subnormal rounding errors, each of size default_abs, are each amplified by up to (1 + default_rel)^n2 subsequent multiplications. Numerous lemmas establish how g1 grows as its arguments increase, supporting inductive error analyses. See also: Kellison et al., "LAProof:...", equation (5).

Hint database: All positivity, monotonicity, and ordering lemmas for default_rel, default_abs, g, and g1 are registered in the commonDB hint database for use with auto and eauto.

From LAProof.accuracy_proofs Require Import preamble real_lemmas.

Import Zorder.

Global Definitions and Setup

The definitions below are independent of any particular floating-point type and are available without opening Section WithType.

Definition rounded t r :=
  Generic_fmt.round Zaux.radix2 (SpecFloat.fexp (fprec t) (femax t))
    (BinarySingleNaN.round_mode BinarySingleNaN.mode_NE) r.

Definition neg_zero {t : type} := Binary.B754_zero (fprec t) (femax t) true.
Definition pos_zero {t : type} := Binary.B754_zero (fprec t) (femax t) false.

Definition Beq_dec_t {t : type} (x y : ftype t) : {x = y } + {x ≠ y } :=
  Beq_dec (fprec t) (femax t) x y.

Create HintDb commonDB discriminated.

Global Hint Resolve
  bpow_gt_0 bpow_ge_0 pos_INR lt_0_INR pow_le : commonDB.

Delimit Scope R_scope with Re.
Open Scope R_scope.

Lemma rev_list_rev : @rev = @List.rev.

Lemma size_not_empty_nat {A} (l : seq A) : l ≠ [] → Nat.le 1 (size l).

Parameterization by Nans and type

Any IEEE floating-point implementation must instantiate,

precision (mantissa size) and exponent size
propagation rules for Not-a-Number

The Rocq types for these are, respectively, type and FPCore.Nans. LAProof's accuracy proofs are, in general, parameterized to work with any instantiation of these. We express that by Rocq's Section and Context commands.

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

Zero Predicates

Definition iszero {t} (x : ftype t) : bool :=
  match x with
  | Binary.B754_zero _ _ _ ⇒ true
  | _ ⇒ false
  end.

Lemma iszeroR_iszeroF : ∀ x : ftype t,
    Binary.is_finite x → FT2R x = R0 → iszero x.

Number of Nonzeros

Definition nnzF : seq (ftype t) → nat :=
count (fun x ⇒ negb (iszero x)).

Definition nnzR : seq R → nat :=
count (fun x ⇒ negb (0 == x)).

Lemma nnzF_zero l : (nnzF l == 0%nat) = (size l == count iszero l ).

Lemma nnzR_zero l : (nnzR l == 0%nat) = (size l == count (eq_op 0) l ).

Lemma nnzF_lemma l : (nnzF l == 0%nat) = all iszero l.

Lemma nnzR_lemma l : (nnzR l == 0%nat) = (all (eq_op R0) l ).

Lemma nnzF_is_zero_cons a l : nnzF (a :: l) == 0%nat → nnzF l == 0%nat.

Lemma nnzR_is_zero_cons a l : nnzR (a :: l) == 0%nat → nnzR l == 0%nat.

Lemma nnzR_cons l : nnzR (0%Re :: l) == nnzR l.

Error Bound Constants

Fundamental error parameters for a floating-point type t . All ordering and positivity lemmas are collected in commonDB.

Definition default_rel : R :=
/ 2 × Raux.bpow Zaux.radix2 (- fprec t + 1).

Definition default_abs : R :=
/ 2 × Raux.bpow Zaux.radix2 (3 - femax t - fprec t).

Lemma default_rel_sep_0 : default_rel ≠ R0.
Hint Resolve default_rel_sep_0 : commonDB.

Lemma default_rel_gt_0 : 0 < default_rel.
Hint Resolve default_rel_gt_0 : commonDB.

Lemma default_rel_ge_0 : 0 ≤ default_rel.
Hint Resolve default_rel_ge_0 : commonDB.

Lemma default_rel_plus_1_ge_1 : 1 ≤ 1 + default_rel.
Hint Resolve default_rel_plus_1_ge_1 : commonDB.

Lemma default_rel_plus_0_ge_1 : 0 ≤ 1 + default_rel.
Hint Resolve default_rel_plus_0_ge_1 : commonDB.

Lemma default_rel_plus_1_gt_1 : 1 < 1 + default_rel.
Hint Resolve default_rel_plus_1_gt_1 : commonDB.

Lemma default_rel_plus_1_gt_0 : 0 < 1 + default_rel.
Hint Resolve default_rel_plus_1_gt_0 : commonDB.

Lemma default_rel_plus_1_ge_1' n : 1 ≤ (1 + default_rel) ^ n.
Hint Resolve default_rel_plus_1_ge_1' : commonDB.

Lemma default_abs_gt_0 : 0 < default_abs.
Hint Resolve default_abs_gt_0 : commonDB.

Lemma default_abs_ge_0 : 0 ≤ default_abs.
Hint Resolve default_abs_ge_0 : commonDB.

Lemma abs_le_rel : default_abs ≤ default_rel.

fmax is the largest finite value representable in type t .

Definition fmax := bpow Zaux.radix2 (femax t).

Lemma bpow_femax_lb : (1 < femax t)%Z.

Lemma bpow_fmax_lb_4 :
4 ≤ fmax.

Lemma bpow_fprec_lb_2 :
2 ≤ bpow Zaux.radix2 (fprec t).

Upper bounds on default_abs and default_rel

The default absolute rounding error default_abs is at most fmax.

Lemma default_abs_le_fmax :
default_abs ≤ fmax.

default_abs t is at most 1.

Lemma default_abs_ub :
default_abs ≤ 1.

default_rel t is at most 1.

Lemma default_rel_ub :
  default_rel ≤ 1.

End WithType.

Global Hint Resolve
  default_rel_sep_0
  default_rel_gt_0
  default_rel_ge_0
  default_rel_plus_1_ge_1
  default_rel_plus_1_gt_0
  default_rel_plus_1_ge_1'
  default_abs_ge_0
  default_abs_gt_0
  default_rel_plus_1_ge_1
  abs_le_rel
  default_rel_plus_0_ge_1
: commonDB.

Error Accumulation Factors

g n and g1 n1 n2 bound accumulated rounding errors over sequences of floating-point operations. The commonDB hint database is populated with the lemmas below to support automated error bound proofs.

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

Notation D := (@default_rel t).
Notation E := (@default_abs t).

Relative Error Factor g

Definition g (n : nat) : R := (1 + D) ^ n - 1.

Lemma g_pos n : 0 ≤ g n.
Hint Resolve g_pos : commonDB.

Lemma le_g_Sn n : g n ≤ g (S n).
Hint Resolve le_g_Sn : commonDB.

Lemma d_le_g n : D ≤ g (n + 1).
Hint Resolve d_le_g : commonDB.

Lemma d_le_g_1 n : (1 ≤ n)%nat → D ≤ g n.
Hint Resolve d_le_g_1 : commonDB.

Lemma one_plus_d_mul_g a n :
(1 + D) × g n × a + D × a = g (n + 1) × a.
Hint Resolve one_plus_d_mul_g : commonDB.

Mixed Error Factor g1

Definition g1 (n1 n2 : nat) : R :=
  INR n1 × E × (1 + g n2).

Lemma g1_pos n m : 0 ≤ g1 n m.
Hint Resolve g1_pos : commonDB.

Lemma one_plus_d_mul_g1 n :
  (1 ≤ n)%nat →
  g1 n (n - 1) × (1 + D) = g1 n n.
Hint Resolve g1_pos : commonDB.

Lemma one_plus_d_mul_g1' n m :
  g1 n m × (1 + D) = g1 n (S m).
Hint Resolve g1_pos : commonDB.

Hint Resolve fprec_lt_femax : commonDB.

Lemma e_le_g1 n : (1 ≤ n)%nat → E ≤ g1 n n.
Hint Resolve e_le_g1 : commonDB.

Lemma plus_d_e_g1_le' n m :
  (1 ≤ n)%nat → (1 ≤ m)%nat →
  g1 n m + (1 + D) × E ≤ g1 (S n) m.
Hint Resolve plus_d_e_g1_le' : commonDB.

Lemma mult_d_e_g1_le' n m :
  (1 ≤ n)%nat → (1 ≤ m)%nat →
  g1 n m × (1 + D) + E ≤ g1 (S n) (S m).
Hint Resolve mult_d_e_g1_le' : commonDB.

Lemma plus_d_e_g1_le n :
  (1 ≤ n)%nat →
  g1 n n + (1 + D) × E ≤ g1 (S n) n.
Hint Resolve plus_d_e_g1_le : commonDB.

Lemma plus_e_g1_le n : g1 n n + E ≤ g1 (S n) n.
Hint Resolve plus_e_g1_le : commonDB.

Lemma g1n_le_g1Sn n :
  (1 ≤ n)%nat →
  g1 n (n - 1) ≤ g1 (S n) (S (n - 1)).
Hint Resolve g1n_le_g1Sn : commonDB.

Lemma g1n_le_g1Sn' n : g1 n n ≤ g1 (S n) (S n).
Hint Resolve g1n_le_g1Sn' : commonDB.

Lemma g1n_lt_g1Sn n :
  (1 ≤ n)%nat →
  g1 n (n - 1) < g1 (S n) (S (n - 1)).

Floating-Point Operation Lemmas

Structural identities: rounding is idempotent on exact FP values (round_FT2R); signed-zero behavior of BPLUS and BMINUS; commutativity of BPLUS; and the equivalence BMINUS x y = BPLUS x (BOPP y).

Lemma round_FT2R a :
  Generic_fmt.round Zaux.radix2 (SpecFloat.fexp (fprec t) (femax t))
    (BinarySingleNaN.round_mode BinarySingleNaN.mode_NE) (FT2R a) = @FT2R t a.

Lemma BMINUS_neg_zero : ∀ c : ftype t, feq (BMINUS neg_zero (BOPP c)) c.

Lemma foldl_congr :
  ∀ (op : ftype t → ftype t → ftype t)
         (Hop : ∀ x y, feq x y → ∀ x' y', feq x' y' →
                   feq (op x x') (op y y'))
         (u v : ftype t) al bl,
    feq u v → Forall2 feq al bl → feq (foldl op u al) (foldl op v bl).

Lemma BPLUS_neg_zero : ∀ c : ftype t, feq (BPLUS c neg_zero) c.

Lemma BPLUS_comm : ∀ x y : ftype t, feq (BPLUS x y) (BPLUS y x).

Lemma MINUS_PLUS_BOPP : ∀ x y : ftype t, feq (BMINUS x y) (BPLUS x (BOPP y)).

End WithType.

Global Hint Resolve
  g_pos
  le_g_Sn
  d_le_g
  d_le_g_1
  g1_pos
  plus_d_e_g1_le'
  one_plus_d_mul_g1
  e_le_g1
  mult_d_e_g1_le'
  plus_d_e_g1_le
  plus_e_g1_le
  g1n_le_g1Sn
  Rplus_le_lt_compat
  Rmult_le_lt_compat
  g1n_lt_g1Sn
: commonDB.

Automation

field_simplify_Rabs reduces a goal of the form Rabs e ≤ _ by simplifying the expression e and splitting denominator non-zero side conditions.

Ltac field_simplify_Rabs :=
  match goal with
  | |- Rabs ?a ≤ _ ⇒
      field_simplify a;
      repeat split;
      try match goal with
          | |- ?z ≠ 0 ⇒ field_simplify z
          end
  end.