diff options
Diffstat (limited to 'theories/Reals/Machin.v')
| -rw-r--r-- | theories/Reals/Machin.v | 168 |
1 files changed, 168 insertions, 0 deletions
diff --git a/theories/Reals/Machin.v b/theories/Reals/Machin.v new file mode 100644 index 00000000..6b91719d --- /dev/null +++ b/theories/Reals/Machin.v @@ -0,0 +1,168 @@ +Require Import Fourier. +Require Import Rbase. +Require Import Rtrigo1. +Require Import Ranalysis_reg. +Require Import Rfunctions. +Require Import AltSeries. +Require Import Rseries. +Require Import SeqProp. +Require Import PartSum. +Require Import Ratan. + +Local Open Scope R_scope. + +(* Proving a few formulas in the style of John Machin to compute Pi *) + +Definition atan_sub u v := (u - v)/(1 + u * v). + +Lemma atan_sub_correct : + forall u v, 1 + u * v <> 0 -> -PI/2 < atan u - atan v < PI/2 -> + -PI/2 < atan (atan_sub u v) < PI/2 -> + atan u = atan v + atan (atan_sub u v). +intros u v pn0 uvint aint. +assert (cos (atan u) <> 0). + destruct (atan_bound u); apply Rgt_not_eq, cos_gt_0; auto. + rewrite <- Ropp_div; assumption. +assert (cos (atan v) <> 0). + destruct (atan_bound v); apply Rgt_not_eq, cos_gt_0; auto. + rewrite <- Ropp_div; assumption. +assert (t : forall a b c, a - b = c -> a = b + c) by (intros; subst; field). +apply t, tan_is_inj; clear t; try assumption. +rewrite tan_minus; auto. + rewrite !atan_right_inv; reflexivity. +apply Rgt_not_eq, cos_gt_0; rewrite <- ?Ropp_div; tauto. +rewrite !atan_right_inv; assumption. +Qed. + +Lemma tech : forall x y , -1 <= x <= 1 -> -1 < y < 1 -> + -PI/2 < atan x - atan y < PI/2. +assert (ut := PI_RGT_0). +intros x y [xm1 x1] [ym1 y1]. +assert (-(PI/4) <= atan x). + destruct xm1 as [xm1 | xm1]. + rewrite <- atan_1, <- atan_opp; apply Rlt_le, atan_increasing. + assumption. + solve[rewrite <- xm1, atan_opp, atan_1; apply Rle_refl]. +assert (-(PI/4) < atan y). + rewrite <- atan_1, <- atan_opp; apply atan_increasing. + assumption. +assert (atan x <= PI/4). + destruct x1 as [x1 | x1]. + rewrite <- atan_1; apply Rlt_le, atan_increasing. + assumption. + solve[rewrite x1, atan_1; apply Rle_refl]. +assert (atan y < PI/4). + rewrite <- atan_1; apply atan_increasing. + assumption. +rewrite Ropp_div; split; fourier. +Qed. + +(* A simple formula, reasonably efficient. *) +Lemma Machin_2_3 : PI/4 = atan(/2) + atan(/3). +assert (utility : 0 < PI/2) by (apply PI2_RGT_0). +rewrite <- atan_1. +rewrite (atan_sub_correct 1 (/2)). + apply f_equal, f_equal; unfold atan_sub; field. + apply Rgt_not_eq; fourier. + apply tech; try split; try fourier. +apply atan_bound. +Qed. + +Lemma Machin_4_5_239 : PI/4 = 4 * atan (/5) - atan(/239). +rewrite <- atan_1. +rewrite (atan_sub_correct 1 (/5)); + [ | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +replace (4 * atan (/5) - atan (/239)) with + (atan (/5) + (atan (/5) + (atan (/5) + (atan (/5) + - + atan (/239))))) by ring. +apply f_equal. +replace (atan_sub 1 (/5)) with (2/3) by + (unfold atan_sub; field). +rewrite (atan_sub_correct (2/3) (/5)); + [apply f_equal | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +replace (atan_sub (2/3) (/5)) with (7/17) by + (unfold atan_sub; field). +rewrite (atan_sub_correct (7/17) (/5)); + [apply f_equal | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +replace (atan_sub (7/17) (/5)) with (9/46) by + (unfold atan_sub; field). +rewrite (atan_sub_correct (9/46) (/5)); + [apply f_equal | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +rewrite <- atan_opp; apply f_equal. +unfold atan_sub; field. +Qed. + +Lemma Machin_2_3_7 : PI/4 = 2 * atan(/3) + (atan (/7)). +rewrite <- atan_1. +rewrite (atan_sub_correct 1 (/3)); + [ | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +replace (2 * atan (/3) + atan (/7)) with + (atan (/3) + (atan (/3) + atan (/7))) by ring. +apply f_equal. +replace (atan_sub 1 (/3)) with (/2) by + (unfold atan_sub; field). +rewrite (atan_sub_correct (/2) (/3)); + [apply f_equal | apply Rgt_not_eq; fourier | apply tech; try split; fourier | + apply atan_bound ]. +apply f_equal; unfold atan_sub; field. +Qed. + +(* More efficient way to compute approximations of PI. *) + +Definition PI_2_3_7_tg n := + 2 * Ratan_seq (/3) n + Ratan_seq (/7) n. + +Lemma PI_2_3_7_ineq : + forall N : nat, + sum_f_R0 (tg_alt PI_2_3_7_tg) (S (2 * N)) <= PI / 4 <= + sum_f_R0 (tg_alt PI_2_3_7_tg) (2 * N). +Proof. +assert (dec3 : 0 <= /3 <= 1) by (split; fourier). +assert (dec7 : 0 <= /7 <= 1) by (split; fourier). +assert (decr : Un_decreasing PI_2_3_7_tg). + apply Ratan_seq_decreasing in dec3. + apply Ratan_seq_decreasing in dec7. + intros n; apply Rplus_le_compat. + apply Rmult_le_compat_l; [ fourier | exact (dec3 n)]. + exact (dec7 n). +assert (cv : Un_cv PI_2_3_7_tg 0). + apply Ratan_seq_converging in dec3. + apply Ratan_seq_converging in dec7. + intros eps ep. + assert (ep' : 0 < eps /3) by fourier. + destruct (dec3 _ ep') as [N1 Pn1]; destruct (dec7 _ ep') as [N2 Pn2]. + exists (N1 + N2)%nat; intros n Nn. + unfold PI_2_3_7_tg. + rewrite <- (Rplus_0_l 0). + apply Rle_lt_trans with + (1 := R_dist_plus (2 * Ratan_seq (/3) n) 0 (Ratan_seq (/7) n) 0). + replace eps with (2 * eps/3 + eps/3) by field. + apply Rplus_lt_compat. + unfold R_dist, Rminus, Rdiv. + rewrite <- (Rmult_0_r 2), <- Ropp_mult_distr_r_reverse. + rewrite <- Rmult_plus_distr_l, Rabs_mult, (Rabs_pos_eq 2);[|fourier]. + rewrite Rmult_assoc; apply Rmult_lt_compat_l;[fourier | ]. + apply (Pn1 n); omega. + apply (Pn2 n); omega. +rewrite Machin_2_3_7. +rewrite !atan_eq_ps_atan; try (split; fourier). +unfold ps_atan; destruct (in_int (/3)); destruct (in_int (/7)); + try match goal with id : ~ _ |- _ => case id; split; fourier end. +destruct (ps_atan_exists_1 (/3)) as [v3 Pv3]. +destruct (ps_atan_exists_1 (/7)) as [v7 Pv7]. +assert (main : Un_cv (sum_f_R0 (tg_alt PI_2_3_7_tg)) (2 * v3 + v7)). + assert (main :Un_cv (fun n => 2 * sum_f_R0 (tg_alt (Ratan_seq (/3))) n + + sum_f_R0 (tg_alt (Ratan_seq (/7))) n) (2 * v3 + v7)). + apply CV_plus;[ | assumption]. + apply CV_mult;[ | assumption]. + exists 0%nat; intros; rewrite R_dist_eq; assumption. + apply Un_cv_ext with (2 := main). + intros n; rewrite scal_sum, <- plus_sum; apply sum_eq; intros. + rewrite Rmult_comm; unfold PI_2_3_7_tg, tg_alt; field. +intros N; apply (alternated_series_ineq _ _ _ decr cv main). +Qed. |