diff options
Diffstat (limited to 'kvx')
-rw-r--r-- | kvx/FPDivision.v | 179 |
1 files changed, 90 insertions, 89 deletions
diff --git a/kvx/FPDivision.v b/kvx/FPDivision.v index 9bfacb80..f22387d6 100644 --- a/kvx/FPDivision.v +++ b/kvx/FPDivision.v @@ -500,103 +500,104 @@ Proof. assert(Rabs (round radix2 (FLT_exp (-1074) 53) ZnearestE (IZR a' * inv_b_r) - (IZR a' * inv_b_r)) <= 1/512)%R as R1 by gappa. - rewrite <- div_approx_reals_correct with (x := B2R _ _ prod). - 2: lia. + assert ( (Rabs (B2R 53 1024 prod - IZR (Int.unsigned a) / IZR (Int.unsigned b)) < 1 / 2)%R) as NEAR. { - unfold div_approx_reals. - fold a' b' prod' prod_r. - unfold Int64.mul. - rewrite Int64.unsigned_repr by (cbn; lia). - rewrite Int64.unsigned_repr by (cbn; lia). - unfold Int64.sub. - rewrite Int64.unsigned_repr by (cbn; lia). - rewrite Int64.unsigned_repr by (cbn; nia). - assert (FMA_RANGE : Int64.min_signed <= a' - prod_r * b' <= Int64.max_signed). - { cbn. - unfold prod_r. - rewrite <- C1E in R1. - assert (IZR (-9223372036854775808) <= IZR (a' - ZnearestE prod' * b') <= IZR 9223372036854775807)%R. - 2: split; apply le_IZR; lra. - rewrite minus_IZR. - rewrite mult_IZR. - replace (IZR (ZnearestE prod')) with (prod' + (IZR (ZnearestE prod') - prod'))%R by ring. - pose proof (Znearest_imp2 (fun x => negb (Z.even x)) prod') as R2. - set (delta1 := (IZR (ZnearestE prod') - prod')%R) in *. - replace prod' with ((prod' - IZR a' * inv_b_r) + IZR a' * (inv_b_r - 1 / IZR b') - + IZR a' / IZR b')%R by (field; lra). - set (delta2 := (inv_b_r - 1 / IZR b')%R) in *. - set (delta3 := (prod' - IZR a' * inv_b_r)%R) in *. - replace (IZR a' - (delta3 + IZR a' * delta2 + IZR a' / IZR b' + delta1) * IZR b')%R with - (- (delta3 + IZR a' * delta2 + delta1) * IZR b')%R by (field; lra). + unfold prod. + pose proof (Bmult_correct 53 1024 (@eq_refl _ Lt) (@eq_refl _ Lt) Float.binop_nan mode_NE (BofZ 53 1024 (@eq_refl _ Lt) (@eq_refl _ Lt) a') inv_b) as C2. + rewrite C0E in C2. + rewrite Rlt_bool_true in C2; cycle 1. + { clear C2. + apply (Rabs_relax (bpow radix2 64)). + { apply bpow_lt. reflexivity. } + cbn. + fold inv_b_r. + replace inv_b_r with (1 / IZR b' + (inv_b_r - 1 / IZR b'))%R by ring. + set (delta := (inv_b_r - 1 / IZR b')%R) in *. unfold approx_inv_thresh in *. - assert (Rabs(IZR a' * delta2) <= 1/4)%R as R4. - { apply Rabs_le. - split; - nra. } - set (delta4 := (IZR a' * delta2)%R) in *. - apply Rabs_def2b in R1. - apply Rabs_def2b in R2. - apply Rabs_def2b in R4. + gappa. + } + destruct C2 as (C2E & C2F & _). + rewrite C2E. + fold inv_b_r a' b'. + replace ((round radix2 (FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) (IZR a' * inv_b_r)) - + (IZR a' / IZR b'))%R with + (((round radix2 (FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) (IZR a' * inv_b_r)) - + (IZR a' * inv_b_r)) + + (IZR a' * (inv_b_r - 1 / IZR b')))%R by (field ; lra). + set(delta := (inv_b_r - 1 / IZR b')%R) in *. + cbn. + (* assert(Rabs + (round radix2 (FLT_exp (-1074) 53) ZnearestE (IZR a' * inv_b_r) - (IZR a' * inv_b_r)) <= 1/512)%R as R1 by gappa. *) + unfold approx_inv_thresh in *. + assert (Rabs(IZR a' * delta) <= 3/8)%R as R2. + { apply Rabs_le. split; nra. } - case Z.ltb_spec; intro CMP. - - replace (Int64.lt (Int64.repr (a' - prod_r * b')) Int64.zero) with true; cycle 1. - { unfold Int64.lt. - change (Int64.signed Int64.zero) with 0. - rewrite Int64.signed_repr by exact FMA_RANGE. - rewrite zlt_true by lia. - reflexivity. - } - cbn. - f_equal. - admit. - - replace (Int64.lt (Int64.repr (a' - prod_r * b')) Int64.zero) with false; cycle 1. - { unfold Int64.lt. - change (Int64.signed Int64.zero) with 0. - rewrite Int64.signed_repr by exact FMA_RANGE. - rewrite zlt_false by lia. - reflexivity. - } - cbn. - unfold Int64.loword. - rewrite Int64.unsigned_repr by (cbn; lia). - reflexivity. + rewrite <- C1E. + rewrite <- C1E in R1. + pose proof (Rabs_triang (prod' - IZR a' * inv_b_r) (IZR a' * delta))%R. + lra. } - unfold prod. - pose proof (Bmult_correct 53 1024 (@eq_refl _ Lt) (@eq_refl _ Lt) Float.binop_nan mode_NE (BofZ 53 1024 (@eq_refl _ Lt) (@eq_refl _ Lt) a') inv_b) as C2. - rewrite C0E in C2. - rewrite Rlt_bool_true in C2; cycle 1. - { clear C2. - apply (Rabs_relax (bpow radix2 64)). - { apply bpow_lt. reflexivity. } - cbn. - fold inv_b_r. - replace inv_b_r with (1 / IZR b' + (inv_b_r - 1 / IZR b'))%R by ring. - set (delta := (inv_b_r - 1 / IZR b')%R) in *. + rewrite <- div_approx_reals_correct with (x := B2R _ _ prod); cycle 1. lia. exact NEAR. + + unfold div_approx_reals. + fold a' b' prod' prod_r. + unfold Int64.mul. + rewrite Int64.unsigned_repr by (cbn; lia). + rewrite Int64.unsigned_repr by (cbn; lia). + unfold Int64.sub. + rewrite Int64.unsigned_repr by (cbn; lia). + rewrite Int64.unsigned_repr by (cbn; nia). + assert (FMA_RANGE : Int64.min_signed <= a' - prod_r * b' <= Int64.max_signed). + { cbn. + unfold prod_r. + rewrite <- C1E in R1. + assert (IZR (-9223372036854775808) <= IZR (a' - ZnearestE prod' * b') <= IZR 9223372036854775807)%R. + 2: split; apply le_IZR; lra. + rewrite minus_IZR. + rewrite mult_IZR. + replace (IZR (ZnearestE prod')) with (prod' + (IZR (ZnearestE prod') - prod'))%R by ring. + pose proof (Znearest_imp2 (fun x => negb (Z.even x)) prod') as R2. + set (delta1 := (IZR (ZnearestE prod') - prod')%R) in *. + replace prod' with ((prod' - IZR a' * inv_b_r) + IZR a' * (inv_b_r - 1 / IZR b') + + IZR a' / IZR b')%R by (field; lra). + set (delta2 := (inv_b_r - 1 / IZR b')%R) in *. + set (delta3 := (prod' - IZR a' * inv_b_r)%R) in *. + replace (IZR a' - (delta3 + IZR a' * delta2 + IZR a' / IZR b' + delta1) * IZR b')%R with + (- (delta3 + IZR a' * delta2 + delta1) * IZR b')%R by (field; lra). unfold approx_inv_thresh in *. - gappa. - } - destruct C2 as (C2E & C2F & _). - rewrite C2E. - fold inv_b_r a' b'. - replace ((round radix2 (FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) (IZR a' * inv_b_r)) - - (IZR a' / IZR b'))%R with - (((round radix2 (FLT_exp (3 - 1024 - 53) 53) (round_mode mode_NE) (IZR a' * inv_b_r)) - - (IZR a' * inv_b_r)) + - (IZR a' * (inv_b_r - 1 / IZR b')))%R by (field ; lra). - set(delta := (inv_b_r - 1 / IZR b')%R) in *. - cbn. - (* assert(Rabs - (round radix2 (FLT_exp (-1074) 53) ZnearestE (IZR a' * inv_b_r) - (IZR a' * inv_b_r)) <= 1/512)%R as R1 by gappa. *) - unfold approx_inv_thresh in *. - assert (Rabs(IZR a' * delta) <= 3/8)%R as R2. - { apply Rabs_le. + assert (Rabs(IZR a' * delta2) <= 1/4)%R as R4. + { apply Rabs_le. + split; + nra. } + set (delta4 := (IZR a' * delta2)%R) in *. + apply Rabs_def2b in R1. + apply Rabs_def2b in R2. + apply Rabs_def2b in R4. split; nra. } - rewrite <- C1E. - rewrite <- C1E in R1. - pose proof (Rabs_triang (prod' - IZR a' * inv_b_r) (IZR a' * delta))%R. - lra. + case Z.ltb_spec; intro CMP. + - replace (Int64.lt (Int64.repr (a' - prod_r * b')) Int64.zero) with true; cycle 1. + { unfold Int64.lt. + change (Int64.signed Int64.zero) with 0. + rewrite Int64.signed_repr by exact FMA_RANGE. + rewrite zlt_true by lia. + reflexivity. + } + cbn. + f_equal. + admit. + - replace (Int64.lt (Int64.repr (a' - prod_r * b')) Int64.zero) with false; cycle 1. + { unfold Int64.lt. + change (Int64.signed Int64.zero) with 0. + rewrite Int64.signed_repr by exact FMA_RANGE. + rewrite zlt_false by lia. + reflexivity. + } + cbn. + unfold Int64.loword. + rewrite Int64.unsigned_repr by (cbn; lia). + reflexivity. Admitted. |