proposition: Theorem sqrt_error_N_FLX x (Fx : format x) : (Rabs (round beta (FLX_exp prec) (Znearest choice) (sqrt x) - sqrt x) <= (1 - 1 / sqrt (1 + 2 * u_ro)) * Rabs (sqrt x))%R.
proof: Proof. assert (Peps := u_ro_pos beta prec). assert (Peps' : (0 < u_ro)%R). { unfold u_ro; apply Rmult_lt_0_compat; [lra|apply bpow_gt_0]. } assert (Pb := om1ds1p2u_ro_pos). assert (Pb' := s1p2u_rom1_pos). destruct (Rle_or_lt x 0) as [Nx|Px]. { rewrite (sqrt_neg _ Nx), round_0, Rabs_R0, Rmult_0_r; [|apply valid_rnd_N]. now unfold Rminus; rewrite Rplus_0_l, Rabs_Ropp, Rabs_R0; right. } destruct (sqrt_error_N_FLX_aux1 _ Fx Px) as (mu, (e, (Fmu, (Hmu, (HmuGe1, HmuLtsqradix))))). pose (t := sqrt x). set (rt := round _ _ _ _). assert (Ht : (t = sqrt mu * bpow e)%R). { unfold t; rewrite Hmu, sqrt_mult_alt; [|now apply (Rle_trans _ _ _ Rle_0_1)]. now rewrite sqrt_bpow. } destruct (sqrt_error_N_FLX_aux2 _ Fmu HmuGe1) as [Hmu'|[Hmu'|Hmu']]. { unfold rt; fold t; rewrite Ht, Hmu', sqrt_1, Rmult_1_l. rewrite round_generic; [|now apply valid_rnd_N|]. { rewrite Rminus_diag_eq, Rabs_R0; [|now simpl]. now apply Rmult_le_pos; [|apply Rabs_pos]. } apply generic_format_bpow'; [now apply FLX_exp_valid|]. unfold FLX_exp; lia. } { assert (Hsqrtmu : (1 <= sqrt mu < 1 + u_ro)%R); [rewrite Hmu'; split|]. { rewrite <- sqrt_1 at 1; apply sqrt_le_1_alt; lra. } { rewrite <- sqrt_square; [|lra]; apply sqrt_lt_1_alt; split; [lra|]. ring_simplify; assert (0 < u_ro ^ 2)%R; [apply pow_lt|]; lra. } assert (Fbpowe : generic_format beta (FLX_exp prec) (bpow e)). { apply generic_format_bpow; unfold FLX_exp; lia. } assert (Hrt : rt = bpow e :> R). { unfold rt; fold t; rewrite Ht; simpl; apply Rle_antisym. { apply round_N_le_midp; [now apply FLX_exp_valid|exact Fbpowe|]. apply (Rlt_le_trans _ ((1 + u_ro) * bpow e)). { now apply Rmult_lt_compat_r; [apply bpow_gt_0|]. } unfold succ; rewrite Rle_bool_true; [|now apply bpow_ge_0]. rewrite ulp_bpow; unfold FLX_exp. unfold Z.sub, u_ro; rewrite !bpow_plus; right; field. } apply round_ge_generic; [now apply FLX_exp_valid|now apply valid_rnd_N|exact Fbpowe|]. rewrite <- (Rmult_1_l (bpow _)) at 1. now apply Rmult_le_compat_r; [apply bpow_ge_0|]. } fold t; rewrite Hrt, Ht, Hmu', <-(Rabs_pos_eq _ Pb), <-Rabs_mult. rewrite Rabs_minus_sym; right; f_equal; field; lra. } assert (Hsqrtmu : (1 + u_ro < sqrt mu)%R). { apply (Rlt_le_trans _ (sqrt (1 + 4 * u_ro))); [|now apply sqrt_le_1_alt]. assert (P1peps : (0 <= 1 + u_ro)%R) by now apply Rplus_le_le_0_compat; [lra|apply Peps]. rewrite <- (sqrt_square (1 + u_ro)); [|lra]. apply sqrt_lt_1_alt; split; [now apply Rmult_le_pos|]. apply (Rplus_lt_reg_r (-1 - 2 * u_ro)); ring_simplify; simpl. rewrite Rmult_1_r; apply Rmult_lt_compat_r; [apply Peps'|]. now apply (Rlt_le_trans _ 1); [apply u_ro_lt_1|lra]. } assert (Hulpt : (ulp beta (FLX_exp prec) t = 2 * u_ro * bpow e)%R). { unfold ulp; rewrite Req_bool_false; [|apply Rgt_not_eq, Rlt_gt]. { unfold u_ro; rewrite <-Rmult_assoc, Rinv_r, Rmult_1_l, <-bpow_plus; [|lra]. f_equal; unfold cexp, FLX_exp. assert (Hmagt : (mag beta t = 1 + e :> Z)%Z). { apply mag_unique. unfold t; rewrite (Rabs_pos_eq _ (Rlt_le _ _ (sqrt_lt_R0 _ Px))). fold t; split. { rewrite Ht; replace (_ - _)%Z with e by ring. rewrite <- (Rmult_1_l (bpow _)) at 1; apply Rmult_le_compat_r. { apply bpow_ge_0. } now rewrite <- sqrt_1; apply sqrt_le_1_alt. } rewrite bpow_plus, bpow_1, Ht; simpl. apply Rmult_lt_compat_r; [now apply bpow_gt_0|]. rewrite <- sqrt_square. { apply sqrt_lt_1_alt; split; [lra|]. apply (Rlt_le_trans _ _ _ HmuLtsqradix); right. now unfold bpow, Z.pow_pos; simpl; rewrite Zmult_1_r, mult_IZR. } apply IZR_le, (Z.le_trans _ 2), Zle_bool_imp_le, radix_prop; lia. } rewrite Hmagt; ring. } rewrite Ht; apply Rmult_lt_0_compat; [|now apply bpow_gt_0]. now apply (Rlt_le_trans _ 1); [lra|rewrite <- sqrt_1; apply sqrt_le_1_alt]. } assert (Pt : (0 < t)%R). { rewrite Ht; apply Rmult_lt_0_compat; [lra|apply bpow_gt_0]. } assert (H : (Rabs ((rt - sqrt x) / sqrt x) <= 1 - 1 / sqrt (1 + 2 * u_ro))%R). { unfold Rdiv; rewrite Rabs_mult, (Rabs_pos_eq (/ _)); [|now left; apply Rinv_0_lt_compat]. apply (Rle_trans _ ((u_ro * bpow e) / t)). { unfold Rdiv; apply Rmult_le_compat_r; [now left; apply Rinv_0_lt_compat|]. apply (Rle_trans _ _ _ (error_le_half_ulp _ _ _ _)). fold t; rewrite Hulpt; right; field. } apply (Rle_trans _ (u_ro / sqrt (1 + 4 * u_ro))). { apply (Rle_trans _ (u_ro * bpow e / (sqrt (1 + 4 * u_ro) * bpow e))). { unfold Rdiv; apply Rmult_le_compat_l; [now apply Rmult_le_pos; [apply Peps|apply bpow_ge_0]|]. apply Rinv_le. { apply Rmult_lt_0_compat; [apply sqrt_lt_R0; lra|apply bpow_gt_0]. } now rewrite Ht; apply Rmult_le_compat_r; [apply bpow_ge_0|apply sqrt_le_1_alt]. } right; field; split; apply Rgt_not_eq, Rlt_gt; [apply sqrt_lt_R0; lra|apply bpow_gt_0]. } apply sqrt_error_N_FLX_aux3. } revert H; unfold Rdiv; rewrite Rabs_mult, Rabs_Rinv; [|fold t; lra]; intro H. apply (Rmult_le_reg_r (/ Rabs (sqrt x))); [apply Rinv_0_lt_compat, Rabs_pos_lt; fold t; lra|]. apply (Rle_trans _ _ _ H); right; field; split; [apply Rabs_no_R0;fold t|]; lra. Qed.
imports: Coq(ZArith Reals Psatz) Core Operations Relative Sterbenz Mult_error