From c5a426f410f95d6e410da5b2afd5a9225e902c3e Mon Sep 17 00:00:00 2001
From: David Monniaux <David.Monniaux@univ-grenoble-alpes.fr>
Date: Mon, 24 Jan 2022 19:21:56 +0100
Subject: twostep_div_longu_smallb_correct

---
 kvx/FPDivision64.v | 82 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 82 insertions(+)

(limited to 'kvx/FPDivision64.v')

diff --git a/kvx/FPDivision64.v b/kvx/FPDivision64.v
index a53169b8..758061b5 100644
--- a/kvx/FPDivision64.v
+++ b/kvx/FPDivision64.v
@@ -1385,3 +1385,85 @@ Proof.
   rewrite zlt_false by lia.
   reflexivity.
 Qed.
+
+Definition twostep_div_longu a b :=
+  let q1 := step1_div_longu a b in
+  let q2 := step2_div_long (Val.subl a (Val.mull b q1)) b in
+  Val.addl q1 q2.
+
+Lemma unsigned_repr_sub :
+  forall a b,
+    Int64.repr (a - b) = Int64.repr (a - Int64.unsigned (Int64.repr b)).
+Proof.
+  intros.
+  apply Int64.eqm_samerepr.
+  apply Int64.eqm_sub.
+  - apply Int64.eqm_refl.
+  - apply Int64.eqm_unsigned_repr.
+Qed.
+
+Lemma unsigned_repr_add :
+  forall a b,
+    Int64.repr (a + b) = Int64.repr (a + Int64.unsigned (Int64.repr b)).
+Proof.
+  intros.
+  apply Int64.eqm_samerepr.
+  apply Int64.eqm_add.
+  - apply Int64.eqm_refl.
+  - apply Int64.eqm_unsigned_repr.
+Qed.
+    
+Lemma twostep_div_longu_smallb_correct :
+    forall a b
+      (b_RANGE : (1 < Int64.unsigned b <= smallb_thresh)%Z),
+          (twostep_div_longu (Vlong a) (Vlong b)) =
+            (Val.maketotal (Val.divlu (Vlong a) (Vlong b))).
+Proof.
+  intros.
+  unfold twostep_div_longu.
+  destruct (step1_div_longu_correct a b b_RANGE) as (q1 & C1R & C1E).
+  rewrite C1R.
+  set (q1' := Int64.unsigned q1) in *.
+  set (b' := Int64.unsigned b) in *.
+  set (a' := Int64.unsigned a) in *.
+  assert ( Z.abs (Int64.signed (Int64.sub a (Int64.mul b q1))) <= smalla_thresh)%Z as r1_SMALL.
+  { unfold smalla_thresh, smallb_approx_range in *.
+    unfold Int64.sub, Int64.mul.
+    fold q1' b' a'.
+    rewrite <- unsigned_repr_sub.
+    rewrite Int64.signed_repr ; cycle 1.
+    { change Int64.min_signed with (-9223372036854775808)%Z.
+      change Int64.max_signed with (9223372036854775807)%Z.
+      lia.
+    }
+    lia.
+  }
+  assert (0 < b')%Z as b_NOT0 by lia.
+  assert (b' <= Int64.max_signed)%Z as b_NOTBIG.
+  { change Int64.max_signed with (9223372036854775807)%Z.
+    unfold smallb_thresh in b_RANGE.
+    lia.
+  }
+  cbn.
+  rewrite (step2_div_long_smalla_correct (Int64.sub a (Int64.mul b q1)) b r1_SMALL b_NOT0 b_NOTBIG).
+  unfold Int64.add, Int64.sub, Int64.mul, Int64.divu.
+  fold q1' b' a'.
+  rewrite <- unsigned_repr_sub.
+  rewrite <- unsigned_repr_add.
+  rewrite Int64.signed_repr ; cycle 1.
+  {
+    change Int64.min_signed with (-9223372036854775808)%Z.
+    change Int64.max_signed with (9223372036854775807)%Z.
+    unfold smallb_approx_range in *.
+    lia.
+  }
+  rewrite Z.add_comm.
+  rewrite <- Z.div_add by lia.
+  replace (a' - b' * q1' + q1' * b')%Z with a' by ring.
+  rewrite Int64.eq_false ; cycle 1.
+  { intro Z. unfold b' in b_NOT0. rewrite Z in b_NOT0.
+    rewrite Int64.unsigned_zero in b_NOT0.
+    lia.
+  }
+  reflexivity.
+Qed.
-- 
cgit