Require Import Rbase. Require Import Rfunctions. Require Import SeqSeries. Require Import Rtrigo_def. Require Import Fourier. Require Import Binomial. Open Local Scope R_scope. (* This should in PartSum *) Lemma sum_f_R0_reverse : forall a N, sum_f_R0 a N = sum_f_R0 (fun i => a (N - i)%nat) N. Proof. intros a; induction N;[simpl; reflexivity | ]. rewrite tech5, (fun b => decomp_sum b (S N)), IHN;[simpl; ring | omega]. Qed. (* This one may be considered specific to this development. *) Lemma rearrange : forall a b N, sum_f_R0 (fun i => sum_f_R0 (fun k => a k * b (i - k)%nat) i) N = sum_f_R0 (fun i => (sum_f_R0 b i) * a (N - i)%nat) N. Proof. induction N; [simpl; ring | ]. rewrite tech5, IHN, tech5, (fun b => decomp_sum b (S N));[ | omega]. rewrite <- minus_n_n, <- minus_n_O. simpl (sum_f_R0 b 0); rewrite <- Rplus_assoc, <- (Rplus_comm (a (S N) * b 0%nat)). rewrite (Rmult_comm (a (S N))); apply f_equal. match goal with |- ?a + _ = ?b => set (u := a); set (v:=b) end. rewrite sum_f_R0_reverse. rewrite <- (sum_eq (fun i => b (S i) * a (N - i)%nat)). unfold u, v; simpl; rewrite <- sum_plus. apply sum_eq; intros i ile; rewrite Rmult_plus_distr_r; tauto. intros i ile; replace (S N) with ((N - i) + S i)%nat. rewrite minus_plus; ring. rewrite <- plus_n_Sm, plus_comm, le_plus_minus_r; tauto. Qed. (* This one should be in Rseries. *) Lemma series_cv_to_term_cv : forall a l, Un_cv (sum_f_R0 a) l -> Un_cv a 0. Proof. intros a l U eps ep. destruct (cv_cauchy_1 a (exist _ l U) eps ep) as [N Hn]. exists (S N); intros [ | n] nN;[inversion nN | apply le_S_n in nN ]. assert (mN : (N <= (S n))%nat) by (apply le_S, nN). assert (t := Hn _ _ nN mN); simpl in t; generalize t; unfold R_dist. unfold Rminus; rewrite Ropp_plus_distr, <- Rplus_assoc, Rplus_opp_r, Rplus_0_l. rewrite Rabs_Ropp, Ropp_0, Rplus_0_r; tauto. Qed. (* This one should be in PartSum. *) Lemma sum_Rlt : forall An Bn N, (forall n, (n <= N)%nat -> An n < Bn n) -> sum_f_R0 An N < sum_f_R0 Bn N. Proof. intros An Bn N D; induction N. solve[simpl; apply D; apply le_n]. simpl; apply Rplus_lt_compat. apply IHN; intros n Hn; apply D; apply le_S; assumption. apply D; apply le_n. Qed. Lemma Mertens : forall a b l1 l2 l3, (Un_cv (sum_f_R0 (fun x => Rabs (a x))) l3) -> (Un_cv (sum_f_R0 a) l1) -> (Un_cv (sum_f_R0 b) l2) -> (Un_cv (sum_f_R0 (fun k => sum_f_R0 (fun p => a p * b (k - p)%nat) k)) (l1 * l2)). Proof. intros a b l1 l2 l3 abscv acv bcv eps ep. assert (G : Un_growing (sum_f_R0 (fun x => Rabs (a x)))). clear; intros k; simpl; assert (t := Rabs_pos (a (S k))); fourier. assert (0 < 4) by fourier; assert (4 <> 0) by (apply not_eq_sym, Rlt_not_eq;auto). assert (l3p : 0 <= l3). solve[apply Rle_trans with (2 := growing_ineq _ _ G abscv 0), Rabs_pos]. assert (e'p : eps / (4 * (l3 + 1)) > 0). apply Rmult_lt_0_compat;[assumption | ]. apply Rinv_0_lt_compat, Rmult_lt_0_compat;[ | apply Rle_lt_trans with l3];fourier. destruct (bcv (eps / (4 * (l3 + 1))) e'p) as [n1 Hn1]. assert (t := maj_by_pos (fun i => Rabs(sum_f_R0 b i - l2))). assert (t' : Un_cv (fun i => Rabs (sum_f_R0 b i - l2)) 0). clear -bcv; intros eps ep;destruct (bcv eps ep) as [N Hn]; exists N; intros n nN. unfold R_dist; rewrite Rminus_0_r, Rabs_Rabsolu; exact (Hn n nN). destruct (t (exist _ 0 t')) as [B [Bp BB]]; clear t t'. assert (e2p : 0 < eps / ((4 * (INR (S n1) * B)))). apply Rmult_lt_0_compat;[assumption | rewrite S_INR]. apply Rinv_0_lt_compat, Rmult_lt_0_compat;[ |apply Rmult_lt_0_compat];auto. assert (t := pos_INR n1); fourier. destruct (series_cv_to_term_cv a l1 acv _ e2p) as [M HM]. assert (e3p : 0 < eps / (2 * (Rabs l2 + 1))). apply Rmult_lt_0_compat;[assumption | ]. assert (t := Rabs_pos l2); apply Rinv_0_lt_compat, Rmult_lt_0_compat; fourier. destruct (acv _ e3p) as [L HL]. set (N := MinMax.max (M + n1 + 1) L). assert (NM := MinMax.MMP.max_lub_l _ _ N (le_n N)). assert (NL := MinMax.MMP.max_lub_r _ _ N (le_n N)). exists N; intros n nN; assert (n0 : (0 < n)%nat) by omega. rewrite (rearrange a b n). rewrite <- (sum_eq (fun i => (sum_f_R0 b i - l2) * a (n - i)%nat + a (n - i)%nat * l2)); [ | intros i _; ring]. rewrite sum_plus, <- scal_sum, (Rmult_comm l2); unfold R_dist. unfold Rminus; rewrite <-Ropp_mult_distr_l_reverse, Rplus_assoc. apply Rle_lt_trans with (1 := Rabs_triang _ _). replace eps with (eps/4 + eps/4 + eps/2) by field; apply Rplus_lt_compat. rewrite (tech2 _ n1 n);[ | omega]. apply Rle_lt_trans with (1 := Rabs_triang _ _); apply Rplus_lt_compat. apply Rle_lt_trans with (1 := Rsum_abs _ _). assert (tmp' : forall i, (i <= n1)%nat -> Rabs ((sum_f_R0 b i - l2) * a (n - i)%nat) < B * (eps / (4 * (INR (S n1) * B)))). intros i in1; rewrite Rabs_mult. apply Rle_lt_trans with (B * Rabs (a (n - i)%nat)). solve[apply Rmult_le_compat_r; [apply Rabs_pos |rewrite <- Rabs_Rabsolu;apply BB]]. apply Rmult_lt_compat_l;[assumption | ]. rewrite <- (Rminus_0_r (a (n - i)%nat)); apply HM. apply le_trans with ((n1 + (M + 1)) - i)%nat. replace n1 with (i + (n1 - i))%nat by (rewrite <- le_plus_minus; omega). rewrite <- plus_assoc, minus_plus; omega. apply minus_le_compat_r; omega. generalize (sum_Rlt _ _ n1 tmp'); clear tmp'. assert (t'' :forall x y z, x = y -> z < x -> z < y) by (intros; fourier). apply t''. rewrite sum_cte, Rmult_assoc, Rmult_comm, Rmult_assoc; unfold Rdiv. field; split; apply not_eq_sym, Rlt_not_eq;[assumption | apply lt_0_INR; omega]. apply Rle_lt_trans with (1 := Rsum_abs _ _). assert (tmp' : forall i, (i <= n - S n1)%nat -> Rabs((sum_f_R0 b (S n1 + i) - l2) * a (n - (S n1 + i))%nat) <= Rabs(a (n - (S n1 + i))%nat) * (eps/(4 * (l3 + 1)))). intros; rewrite Rabs_mult, Rmult_comm. apply Rmult_le_compat_l;[solve [apply Rabs_pos] | apply Rlt_le, Hn1; omega]. apply Rle_lt_trans with (1 := sum_Rle _ _ _ tmp'); rewrite <- scal_sum. apply Rle_lt_trans with (eps/ (4 * (l3 + 1)) * sum_f_R0 (fun i => Rabs(a (n - i)%nat)) n). rewrite (tech2 _ n1 n), (Rmult_plus_distr_l (eps/(4 * (l3 + 1)))). assert (t' : forall x y, 0 <= x -> y <= x + y) by (intros; fourier); apply t'. solve[apply Rmult_le_pos;[apply Rlt_le, e'p | apply Rle_trans with (2:= Rsum_abs _ _), Rabs_pos]]. omega. assert (l3 + 1 <> 0) by (apply not_eq_sym, Rlt_not_eq; fourier). apply Rle_lt_trans with (eps / (4 * (l3 + 1)) * l3). apply Rmult_le_compat_l;[apply Rlt_le, e'p | ]. rewrite <- (sum_f_R0_reverse (fun i => Rabs (a i))). solve[eapply growing_ineq; eauto]. unfold Rdiv; rewrite Rmult_assoc; apply Rmult_lt_compat_l;[exact ep |]. rewrite <- (Rmult_1_r (/ 4)). rewrite Rinv_mult_distr, Rmult_assoc;[ | assumption | assumption]. apply Rmult_lt_compat_l;[apply Rinv_0_lt_compat; fourier | ]. apply Rmult_lt_reg_l with (l3 + 1);[fourier | ]. rewrite <- Rmult_assoc, Rinv_r, Rmult_1_l, Rmult_1_r;[fourier | assumption]. rewrite <- Rmult_plus_distr_r, Rabs_mult. match goal with |- ?a * ?b < _ => apply Rle_lt_trans with (a * (b + 1)) end. apply Rmult_le_compat_l;[apply Rabs_pos | ]; fourier. assert (l2pos' : 0 < Rabs l2 + 1) by (assert (t := Rabs_pos l2); fourier). replace (eps / 2) with ((eps / (2 * (Rabs l2 + 1))) * (Rabs l2 + 1)); [ | field; apply not_eq_sym, Rlt_not_eq; assumption]. apply Rmult_lt_compat_r;[assumption | ]. rewrite <- sum_f_R0_reverse; apply HL; omega. Qed. (* This is general, but especially useful for circle or hyperbolic trigonometry. *) Lemma infinite_sum_decomposition : forall a l1, infinite_sum a l1 -> infinite_sum (fun i => a (2 * i)%nat + a (2 * i + 1)%nat) l1. Proof. intros a l1 cv eps ep; destruct (cv _ ep) as [N Hn]; exists N; intros n nN. destruct n as [ | n]. apply Rle_lt_trans with (2 := Hn _ (le_trans _ _ _ nN (le_S _ _ (le_n 0)))); solve[apply Rle_refl]. rewrite sum_plus; set (x := sum_f_R0 (fun i => a (2 * i)%nat)(S n)). rewrite <- (sum_eq (fun l => a (S (2 * l)))), tech5, <- Rplus_assoc; unfold x. rewrite (sum_decomposition a), <- tech5; apply Hn; omega. intros i _; rewrite plus_comm; reflexivity. Qed. (* This should be in Rseries. *) Lemma Un_cv_ext : forall un vn, (forall n, un n = vn n) -> forall l, Un_cv un l -> Un_cv vn l. Proof. intros un vn quv l P eps ep; destruct (P eps ep) as [N Pn]; exists N. intro n; rewrite <- quv; apply Pn. Qed. (* This should be in Rseries. *) Lemma CV_cte : forall x, Un_cv (fun i => x) x. Proof. intros x e ep; exists 0%nat; unfold R_dist, Rminus; rewrite Rplus_opp_r, Rabs_R0. intros; exact ep. Qed. (* This is specific to trigonometry. *) Lemma ch_cv : forall x, infinite_sum (fun i => Rabs (cos_n i * Rsqr x ^ i)) ((exp x + exp (-x))/2). Proof. intros x; unfold exp; destruct (exist_exp x) as [vex pex]. destruct (exist_exp (-x)) as [vemx pemx]. unfold exp_in in pex, pemx; simpl proj1_sig. apply (Un_cv_ext (sum_f_R0 (fun i => (/INR (fact (2 * i)) * x ^ (2 * i) + /INR (fact (2 * i)) * (- x) ^ (2 * i))/2 + (/INR (fact (2 * i + 1)) * x ^ (2 * i + 1) + /INR (fact (2 * i + 1)) * (- x) ^ (2 * i + 1))/2))). intros i; apply sum_eq; intros n _. rewrite plus_comm; change (x ^ (1 + 2 * n)) with (x * (x ^ (2 * n))). change ((-x) ^ (1 + 2 * n)) with ((-x) * ((-x) ^ (2 * n))). rewrite !pow_Rsqr, <- Rsqr_neg; unfold cos_n, Rdiv. rewrite !Rabs_mult, <- RPow_abs, Rabs_Ropp, !Rabs_pos_eq, pow1, Rmult_1_l. field; solve[split; apply INR_fact_neq_0]. solve[apply pow_le, Rle_0_sqr]. solve[apply Rlt_le, Rinv_0_lt_compat, INR_fact_lt_0]. solve[apply Rle_0_1]. apply (infinite_sum_decomposition (fun i => (/INR (fact i) * x ^ i + /INR (fact i) * (-x) ^ i)/2)). apply (Un_cv_ext (fun i => sum_f_R0 (fun i => /INR (fact i) * x ^ i + /INR (fact i) * (-x) ^ i) i /2)). intros n; unfold Rdiv; rewrite <- scal_sum; ring. apply CV_mult. apply (Un_cv_ext (fun n => sum_f_R0 (fun i => /INR (fact i) * x ^ i) n + sum_f_R0 (fun i => /INR (fact i) * (-x) ^ i) n)). intros n; rewrite <- sum_plus; ring. apply CV_plus; assumption. apply CV_cte. Qed. (* This should be in Rseries. *) Lemma Un_cv_unique : forall a l1 l2, Un_cv a l1 -> Un_cv a l2 -> l1 = l2. Proof. assert (Main : forall a l1 l2, Un_cv a l1 -> Un_cv a l2 -> ~ l1 < l2). intros a l1 l2 u1 u2 d. assert (ep3 : 0 < (l2 - l1) / 3) by (apply Rmult_lt_0_compat; fourier). destruct (u1 _ ep3) as [n1 Hn1]; destruct (u2 _ ep3) as [n2 Hn2]. assert (t : (l2 - l1) < (l2 - l1)/3 + (l2 - l1)/3). set (e := l2 - l1) at 2 3; replace (l2 - l1) with ((a (n1 + n2)%nat - l1) - (a (n1 + n2)%nat - l2)) by ring. apply Rle_lt_trans with (1 := Rle_abs _). apply Rle_lt_trans with (1 := Rabs_triang _ _). apply Rplus_lt_compat; rewrite ?Rabs_Ropp;[apply Hn1 | apply Hn2]; omega. assert (0 < 3) by fourier; apply (Rmult_lt_compat_r 3) in t;[ | assumption]. revert t; unfold Rdiv; rewrite Rmult_plus_distr_r, !Rmult_assoc, Rinv_l, Rmult_1_r; [ | apply not_eq_sym, Rlt_not_eq; assumption]. replace (l2 - l1 + (l2 - l1)) with (2 * l2 - 2 *l1) by ring. replace ((l2 - l1) * 3) with (3 * l2 - 3 * l1) by ring. intros; fourier. intros a l1 l2 u1 u2; destruct (Rtotal_order l1 l2) as [a1 | [g | a2]]. case (Main _ _ _ u1 u2 a1). assumption. case (Main _ _ _ u2 u1 a2). Qed. (* This should be in Binomial. *) Lemma C_0 : forall n, C n 0 = 1. Proof. intros n; unfold C; simpl (INR (fact 0)); rewrite <- minus_n_O, Rmult_1_l. apply Rinv_r, INR_fact_neq_0. Qed. (* This should be in sin_plus. *) Lemma sin_plus : forall x y, sin (x + y) = cos x * sin y + cos y * sin x. Proof. intros x y; unfold cos, sin; destruct (exist_cos (Rsqr x)) as [cx scx]. destruct (exist_sin (Rsqr y)) as [vy svy]. destruct (exist_cos (Rsqr y)) as [cy scy]. destruct (exist_sin (Rsqr x)) as [vx svx]. destruct (exist_sin (Rsqr (x + y))) as [vs svs]. assert (t1 := Mertens _ _ cx vy _ (ch_cv x) scx svy). assert (t2 := Mertens _ _ cy vx _ (ch_cv y) scy svx). assert (t4 := CV_mult _ _ _ _ t1 (CV_cte y)). assert (t5 := CV_mult _ _ _ _ t2 (CV_cte x)); assert (sum := CV_plus _ _ _ _ t4 t5). assert (t8 := CV_mult _ _ _ _ (CV_cte (x + y)) svs); clear -t8 sum. assert (t9 := fun h => Un_cv_ext _ (fun i => (x + y) * sum_f_R0 (fun i => sin_n i * Rsqr (x + y) ^ i) i) h _ sum). rewrite (Rmult_comm y), (Rmult_comm x), <- !Rmult_assoc. apply (fun h => Un_cv_unique _ _ _ t8 (t9 h)); clear; intros n. rewrite <- (Rmult_comm y), <- (Rmult_comm x), !scal_sum, <- sum_plus. apply sum_eq; intros m _; rewrite Rmult_plus_distr_l, <- !pow_Rsqr. rewrite binomial, <- !(Rmult_comm x), <- !(Rmult_comm y), !scal_sum. (* The first term in the left hand side contains the terms that have an even degree in x, the second term contains the others. We nee to make it appear that they are instances of the same function, whic we call f. *) destruct m;[solve[unfold cos_n, sin_n, C; simpl; field] | ]. set (f := fun j => (-1) ^S m / INR(fact (S (2 * S m) - j)) * y ^ (S (2 * S m) - j) * / INR (fact j) * x ^ j). match goal with |- ?a + ?b = ?c + ?d => set (u := a); set (v := b); set (w := c); set (z := d) end. assert (t : forall a b c d e, a*(b*(c*(d*e))) = a * d * (b*c*e)) by (intros; ring). assert (v1 : v = sum_f_R0 (fun i => f (S (2 * i))) (S m)). rewrite sum_f_R0_reverse; apply sum_eq. intros i ile; unfold f, cos_n, sin_n. replace (S (2 * S m) - (S (2 * (S m - i))))%nat with (2 * i)%nat by omega. replace (2 * (S m - i) + 1)%nat with (S (2 * (S m - i))) by omega. unfold Rdiv; rewrite !Rmult_assoc, t, <- pow_add, le_plus_minus_r, <- !pow_Rsqr, <- !tech_pow_Rmult; [ring | omega]. assert (u1 : u = sum_f_R0 (fun i => f (2 * i)%nat) (S m)). unfold u, f, cos_n, sin_n; apply sum_eq; intros i ile. replace (2 * (S m - i) + 1)%nat with (S (2 * (S m) - 2 * i)) by omega. replace (S (2 * S m) - 2 * i)%nat with (S (2 * S m - 2 * i)) by omega. unfold Rdiv; rewrite <- !pow_Rsqr, mult_minus_distr_l, !Rmult_assoc, t, <- pow_add, le_plus_minus_r, <- !tech_pow_Rmult;[ring | omega]. clear t. rewrite v1;rewrite sum_N_predN, <- Rplus_assoc, u1;[simpl (pred (S m)) | omega]. unfold f at 3, Rdiv; rewrite sum_decomposition, <- minus_n_n. simpl (y ^ 0); simpl (fact 0); simpl (INR 1); rewrite Rinv_1, !Rmult_1_r. (* Here we make z follow the same shape as w and we extract the first term. *) set (g := fun i => x ^ i * y ^ (2 * S m - i) * sin_n (S m) * x). assert (plus_compat : forall a b c d, a = c -> b = d -> a + b = c + d) by (intros; subst; ring). assert (z1 : z = sum_f_R0 (fun i => C (2 * S m) (S i) * g i) (S (2 * m)) + (-1) ^ S m / INR (fact (S (2 * S m))) * y ^ (2 * S m) * y). unfold g, z; rewrite decomp_sum, C_0, <- minus_n_O; simpl (x ^ 0);[ | omega]. rewrite !Rmult_1_l, (Rplus_comm (y ^ (2 * S m) * sin_n (S m) * y)). apply plus_compat. replace (pred (2 * S m)) with (S (2 * m)) by omega. apply sum_eq; intros i ile; simpl (x ^ S i). replace (2 * S m - i)%nat with (S (2 * S m - S i))%nat;[simpl; ring | omega]. replace (S (2 * (S m))) with (2 * S m + 1)%nat by ring. rewrite <- (Rmult_comm (sin_n (S m))); reflexivity. assert (w1 : w = sum_f_R0 (fun i => C (2 * S m) i * g i) (S (2 * m)) + (-1)^ S m * /INR (fact (S (2 * S m))) * x ^ S (2 * (S m))). unfold w, g; replace (2 * S m)%nat with (S (S (2 * m))) by ring. rewrite sum_N_predN;[change (pred (S (S (2 * m)))) with (S (2 * m))|omega]. apply plus_compat; [apply sum_eq;intros i ile; ring | ]. rewrite <- !(tech_pow_Rmult x), pascal_step1, <- minus_n_n, C_0;[ | omega]. simpl (y ^ 0); unfold sin_n, Rdiv. replace (2 * S m + 1)%nat with (S(S(S (2 * m)))); ring. rewrite w1, z1; clear u1 v1 u v z1 z w1 w. assert (t : forall a b c d, b = c + d -> b + a = c + a + d); [clear; intros a b c d h1 ; rewrite h1; ring | apply t; unfold f, g]. rewrite <- Rplus_assoc, <- sum_plus, decomp_sum, Rplus_comm;[ | omega]. apply plus_compat. replace (pred (2 * S m)) with (S (2 * m)) by (simpl; ring). apply sum_eq; intros i ile. replace (S (2 * S m) - S i)%nat with (2 * S m - i)%nat by omega. rewrite <- !(Rmult_comm (sin_n (S m))). rewrite !Rmult_assoc, <- !(Rmult_assoc (sin_n (S m))), <- Rmult_plus_distr_r. transitivity (((-1) ^ S m * / INR (fact (2 * S m - i)) * / INR (fact (S i))) * (x ^ i * (y ^ (2 * S m - i) * x))). rewrite <- (tech_pow_Rmult x); unfold Rdiv; ring. unfold sin_n; rewrite pascal;[unfold C | omega]. replace (2 * S m + 1)%nat with (S (2 * (S m))) by ring. replace (S (2 * S m) - S i)%nat with (2 * S m - i)%nat by omega. field; repeat (apply conj); apply INR_fact_neq_0. simpl (INR (fact 0)); simpl (x ^ 0). rewrite <- !minus_n_O, Rinv_1, !Rmult_1_r, <- (tech_pow_Rmult y); ring. Qed. Locate pascal.