1 files changed, 85 insertions, 0 deletions

diff --git a/theories/nix/floats.v b/theories/nix/floats.v
new file mode 100644
index 0000000..246e0c3
--- /dev/null
+++ b/theories/nix/floats.v
@@ -0,0 +1,85 @@
+From stdpp Require Import prelude ssreflect.
+From Flocq.IEEE754 Require Import
+  Binary BinarySingleNaN (mode_NE, mode_DN, mode_UP) Bits.
+From stdpp Require Import options.
+
+Global Arguments B754_zero {_ _}.
+Global Arguments B754_infinity {_ _}.
+Global Arguments B754_nan {_ _}.
+Global Arguments B754_finite {_ _}.
+
+(** The bit representation of floats is not observable in Nix, and it appears
+that only negative NaNs are ever produced. So we setup the Flocq floats in
+the way that it always produces [-NaN], i.e., [indef_nan]. *)
+Definition float := binary64.
+
+Variant round_mode := Floor | Ceil | NearestEven.
+Global Instance round_mode_eq_dec : EqDecision round_mode.
+Proof. solve_decision. Defined.
+
+Module Float.
+  Definition prec : Z := 53.
+  Definition emax : Z := 1024.
+
+  Lemma Hprec : (0 < 53)%Z.
+  Proof. done. Qed.
+  Lemma Hprec_emax : (53 < 1024)%Z.
+  Proof. done. Qed.
+
+  Global Instance inhabited : Inhabited float := populate (B754_zero false).
+
+  Global Instance eq_dec : EqDecision float.
+  Proof.
+    refine (λ f1 f2,
+      match f1, f2 with
+      | B754_zero s1, B754_zero s2 => cast_if (decide (s1 = s2))
+      | B754_infinity s1, B754_infinity s2 => cast_if (decide (s1 = s2))
+      | B754_nan s1 pl1 _, B754_nan s2 pl2 _ =>
+          cast_if_and (decide (s1 = s2)) (decide (pl1 = pl2))
+      | B754_finite s1 m1 e1 _, B754_finite s2 m2 e2 _ =>
+          cast_if_and3 (decide (s1 = s2)) (decide (m1 = m2)) (decide (e1 = e2))
+      | _, _ => right _
+      end); abstract naive_solver (f_equal; apply (proof_irrel _)).
+  Defined.
+
+  Definition of_Z (x : Z) : float :=
+    binary_normalize prec emax (refl_equal _) (refl_equal _) mode_NE x 0 false.
+
+  Definition to_Z (f : float) : option Z :=
+    match f with
+    | B754_zero _ => Some 0
+    | B754_finite s m e _ => Some (Zaux.cond_Zopp s (Zpos m) ≪ e)
+    | _ => None
+    end%Z.
+
+  (** QNaN Floating-Point Indefinite; see Table 4-3. Floating-Point Number and
+  NaN Encodings. *)
+  Definition indef_nan : { f | is_nan prec emax f = true } :=
+    @B754_nan prec emax true (2^51) (refl_equal _) ↾ eq_refl _.
+
+  Definition to_flocq_round_mode (m : round_mode) : BinarySingleNaN.mode :=
+    match m with Floor => mode_DN | Ceil => mode_UP | NearestEven => mode_NE end.
+
+  Definition round (m : round_mode) : float → float :=
+    Bnearbyint prec emax (refl_equal _) (λ _, indef_nan) (to_flocq_round_mode m).
+
+  (* For add: not [mode_DN]; otherwise [+0.0 + -0.0] would yield [-0.0], but
+  [inf / (+0.0 + -0.0)] yields [inf] in C++, so this cannot be the case. *)
+  (* C++ compiles floating point addition to the x86 ADDSD instruction. Looking
+  at the Intel x86 Software Developer's Manual, it seems that the default
+  rounding mode on x86 is Round to Nearest (even); see table 4-8. (In §4.8.4.) *)
+  Definition add : float → float → float :=
+    Bplus _ _ Hprec Hprec_emax (λ _ _, indef_nan) mode_NE.
+  Definition sub : float → float → float :=
+    Bminus _ _ Hprec Hprec_emax (λ _ _, indef_nan) mode_NE.
+  Definition mul : float → float → float :=
+    Bmult _ _ Hprec Hprec_emax (λ _ _, indef_nan) mode_NE.
+  Definition div : float → float → float :=
+    Bdiv _ _ Hprec Hprec_emax (λ _ _, indef_nan) mode_NE.
+
+  Definition eqb (f1 f2 : float) : bool :=
+    bool_decide (b64_compare f1 f2 = Some Eq).
+
+  Definition ltb (f1 f2 : float) : bool :=
+    bool_decide (b64_compare f1 f2 = Some Lt).
+End Float.