proposition: Lemma midpoint_correct : forall x y, sensible_format = true -> real x = true -> real y = true -> (toR x <= toR y)%R -> real (midpoint x y) = true /\ (toR x <= toR (midpoint x y) <= toR y)%R.
proof: Proof. intros x y _. rewrite !real_correct. unfold toR, toX, toF. rewrite <-!B2SF_Prim2B. set (b_x := Prim2B x). set (b_y := Prim2B y). intros Hx Hy Hxy. unfold midpoint. replace (Prim2B (if is_infinity _ then _ else _)) with (if is_infinity ((x + y) / 2) then Prim2B (x / 2 + y / 2) else Prim2B ((x + y) / 2)). 2:{ now case is_infinity. } rewrite is_infinity_equiv. rewrite add_equiv, !div_equiv, add_equiv. fold b_x; fold b_y. set (b2 := Prim2B 2). assert (Nz2 : B2R b2 <> 0%R). { compute; lra. } revert Hx Hxy. set (bplus := Bplus _). set (bdiv := Bdiv _). case b_x; [intros sx..| |intros sx mx ex Hmex]; [ |intro H; discriminate H..| ]; intros _. { revert Hy. case b_y; [intros sy..| |intros sy my ey Hmey]; [ |intro H; discriminate H..| ]; intros _. { now case sx, sy. } case sy; [intro Hy; simpl in Hy|intros _]. { generalize (Generic_proof.FtoR_Rneg radix2 my ey); lra. } change (bplus (B754_zero sx) _) with (B754_finite false my ey Hmey). set (by2 := bdiv (B754_finite false my ey Hmey) b2). elim (Bdiv2_correct (B754_finite false my ey Hmey) eq_refl). fold bdiv; fold b2; fold by2. intros _ [Fy2 [Sy2 Hy2']]; revert Fy2 Sy2 Hy2'. case by2 => [sy2|sy2| |sy2 my2 ey2 Hmey2]; [ |intro H; discriminate H..| ]; intros _; simpl. { intros _ _. split; [reflexivity|split; [now right| ]]. apply Rlt_le, Generic_proof.FtoR_Rpos. } intros ->. change (Z.pos my) with (cond_Zopp false (Z.pos my)). rewrite <-!FtoR_split, !Generic_proof.FtoR_abs. intro H; split; [reflexivity|split; [ |exact H]]. apply Rlt_le, Generic_proof.FtoR_Rpos. } revert Hy. case b_y; [intros sy..| |intros sy my ey Hmey]; [ |intro H; discriminate H..| ]; intros _. { case sx; [intros _|intros Hx; simpl in Hx]. 2:{ generalize (Generic_proof.FtoR_Rpos radix2 mx ex); lra. } change (bplus _ (B754_zero sy)) with (B754_finite true mx ex Hmex). set (bx2 := bdiv (B754_finite true mx ex Hmex) b2). elim (Bdiv2_correct (B754_finite true mx ex Hmex) eq_refl). fold bdiv; fold b2; fold bx2. intros _ [Fx2 [Sx2 Hx2]]; revert Fx2 Sx2 Hx2. case bx2 => [sx2|sx2| |sx2 mx2 ex2 Hmex2]; [ |intro H; discriminate H..| ]; intros _; simpl. { intros _ _. split; [reflexivity|split; [ |now right]]. apply Rlt_le, Generic_proof.FtoR_Rneg. } intros ->. change (Z.neg mx) with (cond_Zopp true (Z.pos mx)). rewrite <-!FtoR_split, !Generic_proof.FtoR_abs. intro H; split; [reflexivity|split]. 2:{ apply Rlt_le, Generic_proof.FtoR_Rneg. } change true with (negb false). rewrite <-!Generic_proof.FtoR_neg. now apply Ropp_le_contravar. } clear x y b_x b_y. set (b_x := B754_finite sx mx ex Hmex). set (b_y := B754_finite sy my ey Hmey). intros Hxy; simpl in Hxy. generalize (Bplus_correct _ _ Hprec Hmax mode_NE b_x b_y eq_refl eq_refl). fold bplus. case Rlt_bool_spec => Hb. { intros [Rxpy [Fxpy Sxpy]]. elim (Bdiv2_correct (bplus b_x b_y) Fxpy). fold bdiv; fold b2. intros _ [Fxpy2 _]. replace (match bdiv _ _ with B754_infinity _ => true | _ => _ end) with false; [ |now revert Fxpy2; case bdiv]. split; [now revert Fxpy2; case bdiv| ]. elim (Bdiv2_correct _ Fxpy); fold bdiv b2. intros Rxpy2 _. simpl. set (rx := FtoR radix2 sx mx ex). set (ry := FtoR radix2 sy my ey). revert Rxpy Rxpy2. set (fexp := FLT.FLT_exp _ _). set (m := round_mode _). intros Rxpy Rxpy2. rewrite <-(Generic_fmt.round_generic radix2 fexp m rx). 2:{ unfold rx; rewrite FtoR_split; change (Defs.F2R _) with (B2R b_x). apply generic_format_B2R. } rewrite <-(Generic_fmt.round_generic radix2 fexp m ry). 2:{ unfold ry; rewrite FtoR_split; change (Defs.F2R _) with (B2R b_y). apply generic_format_B2R. } replace rx with ((rx + rx) / 2)%R; [ |lra]. replace ry with ((ry + ry) / 2)%R; [ |lra]. replace (proj_val _) with (B2R (bdiv (bplus b_x b_y) b2)). 2:{ change (binary_normalize _ _ _ _ _ _ _ _) with (bplus b_x b_y). case bdiv => [s|s| |sb mb eb Hmeb]; [reflexivity..| ]. now unfold B2R; rewrite <-FtoR_split. } rewrite Rxpy2, Rxpy. split; (apply Generic_fmt.round_le; [now apply FLT.FLT_exp_valid|now apply Generic_fmt.valid_rnd_N| ]; unfold Rdiv; apply Rmult_le_compat_r; [lra| ]). { rewrite <-(Generic_fmt.round_generic radix2 fexp m (rx + rx)). { apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } unfold B2R, b_x, b_y; rewrite <-!FtoR_split. now apply Rplus_le_compat_l. } replace (rx + rx)%R with (rx * bpow radix2 1)%R; [ |simpl; lra]. apply mult_bpow_pos_exact_FLT; [ |lia]. unfold rx; rewrite FtoR_split; change (Defs.F2R _) with (B2R b_x). apply generic_format_B2R. } rewrite <-(Generic_fmt.round_generic radix2 fexp m (ry + ry)). { apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } unfold B2R, b_x, b_y; rewrite <-!FtoR_split. now apply Rplus_le_compat_r. } replace (ry + ry)%R with (ry * bpow radix2 1)%R; [ |simpl; lra]. apply mult_bpow_pos_exact_FLT; [ |lia]. unfold ry; rewrite FtoR_split; change (Defs.F2R _) with (B2R b_y). apply generic_format_B2R. } change (binary_overflow _ _ _ _) with (@B2SF FloatOps.prec emax (B754_infinity sx)). intros [H H']; revert H'; rewrite (B2SF_inj _ _ _ _ H); clear H. intro Hsxy; simpl in Hsxy. change (match bdiv _ _ with B754_infinity _ => true | _ => _ end) with true. revert Hb. change (SpecFloat.fexp _ _) with (FLT.FLT_exp (3 - emax - FloatOps.prec) FloatOps.prec). set (fexp := FLT.FLT_exp _ _). set (m := round_mode _). elim (Plus_error.FLT_plus_error_N_ex _ _ _ (fun x : Z => negb (Z.even x)) _ _ (generic_format_B2R _ _ b_x) (generic_format_B2R _ _ b_y)). change (Generic_fmt.Znearest _) with (round_mode mode_NE). unfold emin. fold fexp m. intros eps [Heps ->]. rewrite Rabs_mult. intro Hb. assert (R1peps : (0 < Rabs (1 + eps))%R). { apply Rabs_gt; right. generalize (Rle_trans _ _ _ Heps (Relative.u_rod1pu_ro_le_u_ro _ _)). intro H; generalize (Rabs_le_inv _ _ H); compute; lra. } generalize (Rmult_le_compat_r _ _ _ (Rlt_le _ _ (Rinv_0_lt_compat _ R1peps)) Hb). rewrite Rmult_assoc, Rinv_r, ?Rmult_1_r; [ |lra]. clear Hb; intro Hb. generalize (Rle_trans _ _ _ Hb (Rabs_triang _ _)). clear Hb; intro Hb. assert (Hb' : (1 / 256 <= bpow radix2 emax * / Rabs (1 + eps) - (bpow radix2 emax - bpow radix2 (emax - FloatOps.prec)))%R). { rewrite Rcomplements.Rle_minus_r. apply (Rmult_le_reg_r _ _ _ R1peps). rewrite Rmult_assoc, Rinv_l, ?Rmult_1_r; [ |lra]. refine (Rle_trans _ _ _ (Rmult_le_compat_l _ _ _ _ (Rle_trans _ _ _ (Rabs_triang _ _) (Rplus_le_compat_l _ _ _ Heps))) _). { apply Rplus_le_le_0_compat; [lra| ]. now apply Rle_0_minus, bpow_le; compute. } rewrite (Rabs_pos_eq _ Rle_0_1). compute; lra. } assert (Hx2h : (1 / 256 <= Rabs (toR (B2Prim b_x)))%R). { unfold toR, toX, toF; rewrite Prim2SF_B2Prim; unfold b_x; simpl. apply (Rle_trans _ _ _ Hb'). apply (Rle_trans _ _ _ (Rplus_le_compat_r _ _ _ Hb)). rewrite FtoR_split; change (Defs.F2R _) with (B2R b_x). apply (Rplus_le_reg_r (- Rabs (B2R b_y))). ring_simplify. unfold Rminus; rewrite Rplus_assoc. apply Rplus_le_compat_l. generalize (abs_B2R_le_emax_minus_prec _ emax Hprec b_y). lra. } assert (Hy2h : (1 / 256 <= Rabs (toR (B2Prim b_y)))%R). { unfold toR, toX, toF; rewrite Prim2SF_B2Prim; unfold b_y; simpl. apply (Rle_trans _ _ _ Hb'). apply (Rle_trans _ _ _ (Rplus_le_compat_r _ _ _ Hb)). rewrite FtoR_split; change (Defs.F2R _) with (B2R b_y). apply (Rplus_le_reg_r (- Rabs (B2R b_x))). ring_simplify. unfold Rminus; rewrite Rplus_assoc, Rplus_comm. apply Rplus_le_compat_r. generalize (abs_B2R_le_emax_minus_prec _ emax Hprec b_x). lra. } generalize (div2_correct _ (refl_equal _) Hy2h). generalize (div2_correct _ (refl_equal _) Hx2h). intros Hx2 Hy2. assert (Fx2 : is_finite (bdiv b_x b2) = true). { revert Hx2; unfold toX, toF, div2. rewrite <-B2SF_Prim2B, div_equiv, Prim2B_B2Prim; fold bdiv b2. rewrite Prim2SF_B2Prim. unfold Xdiv', Xbind2; rewrite is_zero_false; [ |lra]. now case bdiv => [s|s| |s m' e Hme]. } assert (Fy2 : is_finite (bdiv b_y b2) = true). { revert Hy2; unfold toX, toF, div2. rewrite <-B2SF_Prim2B, div_equiv, Prim2B_B2Prim; fold bdiv b2. rewrite Prim2SF_B2Prim. unfold Xdiv', Xbind2; rewrite is_zero_false; [ |lra]. now case bdiv => [s|s| |s m' e Hme]. } generalize (Bplus_correct _ _ Hprec Hmax mode_NE _ _ Fx2 Fy2). fold bplus fexp m. replace (B2R (bdiv b_x b2)) with (B2R b_x / 2)%R. 2:{ revert Hx2; unfold toX, toF, div2. rewrite <-B2SF_Prim2B, div_equiv, Prim2B_B2Prim; fold bdiv b2. rewrite Prim2SF_B2Prim. unfold Xdiv', Xbind2; rewrite is_zero_false; [ |lra]. case bdiv => [s|s| |s m' e Hme]; [ |intro H; discriminate H..| ]. { now intro H; inversion H as (H'); simpl; rewrite H', FtoR_split. } intro H; inversion H as (H'); revert H'; simpl. now rewrite !FtoR_split => ->. } replace (B2R (bdiv b_y b2)) with (B2R b_y / 2)%R. 2:{ revert Hy2; unfold toX, toF, div2. rewrite <-B2SF_Prim2B, div_equiv, Prim2B_B2Prim; fold bdiv b2. rewrite Prim2SF_B2Prim. unfold Xdiv', Xbind2; rewrite is_zero_false; [ |lra]. case bdiv => [s|s| |s m' e Hme]; [ |intro H; discriminate H..| ]. { now intro H; inversion H as (H'); simpl; rewrite H', FtoR_split. } intro H; inversion H as (H'); revert H'; simpl. now rewrite !FtoR_split => ->. } rewrite Rlt_bool_true. 2:{ unfold b_x, b_y; rewrite <-Hsxy. case_eq sx => Hsx. { apply (Rle_lt_trans _ (Rabs (B2R b_x))). 2:{ apply abs_B2R_lt_emax. } rewrite Rabs_left1. 2:{ rewrite <-(Generic_fmt.round_0 radix2 fexp m). apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } simpl. change (Z.neg mx) with (cond_Zopp true (Z.pos mx)). change (Z.neg my) with (cond_Zopp true (Z.pos my)). rewrite <-!FtoR_split. generalize (Generic_proof.FtoR_Rneg radix2 mx ex). generalize (Generic_proof.FtoR_Rneg radix2 my ey). lra. } rewrite Rabs_left1. 2:{ simpl. rewrite <-FtoR_split, Hsx. generalize (Generic_proof.FtoR_Rneg radix2 mx ex). lra. } apply Ropp_le_contravar. rewrite <-(Generic_fmt.round_generic radix2 fexp m (B2R b_x)). { apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } replace (B2R b_x) with (B2R b_x / 2 + B2R b_x / 2)%R by field. rewrite <-Hsx; apply Rplus_le_compat_l. apply Rmult_le_compat_r; [lra| ]. now revert Hxy; rewrite !FtoR_split, <-Hsxy. } apply generic_format_B2R. } apply (Rle_lt_trans _ (Rabs (B2R b_y))). 2:{ apply abs_B2R_lt_emax. } rewrite Rabs_pos_eq. 2:{ rewrite <-(Generic_fmt.round_0 radix2 fexp m). apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } simpl. change (Z.pos mx) with (cond_Zopp false (Z.pos mx)). change (Z.pos my) with (cond_Zopp false (Z.pos my)). rewrite <-!FtoR_split. generalize (Generic_proof.FtoR_Rpos radix2 mx ex). generalize (Generic_proof.FtoR_Rpos radix2 my ey). lra. } rewrite Rabs_pos_eq. 2:{ simpl. rewrite <-FtoR_split, <-Hsxy, Hsx. generalize (Generic_proof.FtoR_Rpos radix2 my ey). lra. } rewrite <-(Generic_fmt.round_generic radix2 fexp m (B2R b_y)). { apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } replace (B2R b_y) with (B2R b_y / 2 + B2R b_y / 2)%R by field. rewrite <-Hsx, Hsxy; apply Rplus_le_compat_r. apply Rmult_le_compat_r; [lra| ]. now revert Hxy; rewrite !FtoR_split, Hsxy. } apply generic_format_B2R. } intros [Rx2py2 [Fx2py2 _]]. split. { revert Fx2py2; case bplus => [s|s| |s m' e Hme]; [ |intro H; discriminate H..| ]; reflexivity. } unfold proj_val at -2 3. replace (proj_val _) with (B2R (bplus (bdiv b_x b2) (bdiv b_y b2))). 2:{ now case bplus => [s|s| |s m' e Hme]; [..|simpl; rewrite <-FtoR_split]. } unfold B2SF, b_x, b_y, FtoX; fold b_x b_y. rewrite FtoR_split; change (Defs.F2R _) with (B2R b_x). rewrite FtoR_split; change (Defs.F2R _) with (B2R b_y). rewrite Rx2py2. rewrite <-(Generic_fmt.round_generic radix2 fexp m (B2R b_x)) at 1. 2:{ apply generic_format_B2R. } rewrite <-(Generic_fmt.round_generic radix2 fexp m (B2R b_y)) at 3. 2:{ apply generic_format_B2R. } split. { apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } replace (B2R b_x) with (B2R b_x / 2 + B2R b_x / 2)%R at 1 by field. apply Rplus_le_compat_l. apply Rmult_le_compat_r; [lra| ]. now revert Hxy; rewrite !FtoR_split. } apply Generic_fmt.round_le. { now apply FLT.FLT_exp_valid. } { now apply Generic_fmt.valid_rnd_N. } replace (B2R b_y) with (B2R b_y / 2 + B2R b_y / 2)%R at 2 by field. apply Rplus_le_compat_r. apply Rmult_le_compat_r; [lra| ]. now revert Hxy; rewrite !FtoR_split. Qed.
imports: Coq(ZArith Reals) Int63Compat Coq(Floats Psatz) Flocq(Zaux Raux BinarySingleNaN PrimFloat Sterbenz Mult_error) Missing.Stdlib Missing.Flocq Xreal Basic Sig