Library Plouffe
Require Import ssreflect NPeano Psatz Reals Coquelicot.
(******************************************************************************)
(* *)
(* PROOF OF THE PLOUFFE FORMULA *)
(* *)
(******************************************************************************)
Notation "a %:R " := (INR a) (at level 10).
Open Scope R_scope.
Ltac tlia := try lia.
Ltac tlra := try lra.
Ltac Rcotpatch :=
change (NormedModule.sort R_AbsRing R_NormedModule) with R; repeat
match goal with
|- context[scal ?X ?Y] => change (scal X Y) with (X * Y)
end.
Lemma pow_lt_1 x n: Rabs x < 1 -> Rabs (x ^ (S n)) < 1.
Section Hole.
Definition hole (k : nat) (a : nat -> R) (n : nat) :=
if Nat.eqb (n mod k) 0 then a (n / k)%nat else 0.
Lemma sum_f_R0_holes_small a g k n
(trm := fun i => hole k a i * g i) l :
(S l < k)%nat
-> sum_f_R0 (fun i => trm (S(k * n) + i)%nat) l = 0.
Lemma sumf_f_R0_add_r f m n : (n <> 0)%nat ->
sum_f_R0 f (m + n) =
sum_f_R0 f m + sum_f_R0 (fun i => f (S m + i)%nat) (pred n).
Lemma sum_f_R0_holes k a g n :
(k <> 0)%nat -> sum_f_R0 (fun i => hole k a i * g i) (k * n) =
sum_f_R0 (fun i => a i * g (k * i)%nat) n.
Lemma fill_holes : forall k a x, (k <> 0)%nat ->
ex_pseries a (x ^ k)->
PSeries (hole k a) x = Series (fun n => a n * x ^ (k * n)).
End Hole.
Section Radius.
Variable a : nat -> R.
Hypothesis aP : forall n, 0 <= a n <= 1.
Lemma PS_cv: forall x, 0 <= x -> x < 1 ->
Rbar.Rbar_lt (Rbar.Finite (Rabs x)) (CV_radius a) .
End Radius.
Section Plouffe.
Lemma PSerie_pow_add x k f a:
Series (fun i => (a i ) * x ^ (k + f i)) =
x ^ k * Series (fun i => (a i) * x ^ f i).
Lemma PSeries_hole x a d k : 0 <= x < 1 ->
Series (fun i : nat => (a * x ^ (d + (S k) * i))) =
PSeries (PS_incr_n (hole (S k) (fun _ => a)) d) x.
Lemma RInt_Series_hole n r: 0 <= r -> r < 1 -> (0 < n)%nat ->
RInt (fun x => Series (fun i => x ^ (n - 1 + 8 * i))) 0 r =
PSeries (PS_Int (PS_incr_n (hole 8 (fun _ => 1)) (n - 1))) r.
Lemma is_RInt_pow n b: is_RInt (fun x => (x ^ n)) 0 b (b ^ (S n)/ (S n)%:R).
Lemma RInt_pow n b : RInt (fun x => x ^ n) 0 b = b ^ (S n)/ (S n)%:R.
Lemma Series_16_Rint n (f := fun i => / (16 ^ i * (8 * i + n)%:R)) :
(0 < n)%nat ->
Series f =
sqrt 2 ^ n *
RInt (fun x => x ^ (n - 1) / (1 - x ^ 8)) 0 (/(sqrt 2)).
(* have hr2:= sqrt2_neq_0.*)
Lemma Pos_x2 x : 0 < x ^ 2 - 2 * x + 2.
Lemma Derive_x2 x :
Derive (fun x => -2 * ln (x ^ 2 - 2 * x + 2)) x =
(4 - 4 * x) / (x ^ 2 - 2 * x + 2).
Lemma Continuity_x2 x :
continuity_pt (fun x0 : R => (4 - 4 * x0) / (x0 ^ 2 - 2 * x0 + 2)) x.
Lemma Rint_x2 :
RInt (fun x => (4 - 4 * x) / (x ^ 2 - 2 * x + 2)) 0 1 = 2 * (ln 2 - ln 1).
Lemma Sqrt_x2 x : -sqrt 2 < x < sqrt 2 -> 0 < 2 - x ^ 2.
Lemma Sqrt2_m1 : 0 < sqrt 2 - 1.
Lemma X_between x : 0 <= x <= 1 -> - sqrt 2 < x < sqrt 2.
Lemma Derive_ln_x2 x :
-sqrt 2 < x < sqrt 2 ->
Derive (fun x => 2 * ln (2 - x ^ 2)) x = (4 * x) / (x ^ 2 - 2).
Lemma Continuity_inv_x2 x (Hx : 0 <= x <= 1) :
continuity_pt (fun x0 : R => 4 * x0 / (x0 ^ 2 - 2)) x.
Lemma Rint_inv_x2 :
RInt (fun x => (4 * x) / (x ^ 2 - 2)) 0 1 = 2 * (ln 1 - ln 2).
Lemma Derive_4atan x :
Derive (fun x => 4 * atan (x - 1)) x = 4 / (x ^ 2 - 2 * x + 2).
Lemma Continuity_4_inv_x2 x :
continuity_pt (fun x0 : R => 4 / (x0 ^ 2 - 2 * x0 + 2)) x.
Lemma Rint_4x2 : RInt (fun x => 4 / (x ^ 2 - 2 * x + 2)) 0 1 = PI.
Lemma Rint_x8_PI (f := fun (n : nat) x =>
sqrt 2 ^ (1 + n) * (x ^ n / (1 - x ^ 8))) :
RInt (fun x => (4 * f 0%nat x - 2 * f 3%nat x - f 4%nat x - f 5%nat x))
0 (/sqrt 2) = PI.
Lemma ex_RInt_x8 n (f := fun x => sqrt 2 ^ (1 + n) * (x ^ n / (1 - x ^ 8))) :
ex_RInt f 0 (/ sqrt 2).
Lemma Plouffe_PI (f := fun n i => / (16 ^ i * (8 * i + n)%:R)) :
PI = 4 * (Series (f 1%nat)) - 2 * (Series (f 4%nat))
- (Series (f 5%nat)) - (Series (f 6%nat)).
End Plouffe.
(******************************************************************************)
(* *)
(* PROOF OF THE PLOUFFE FORMULA *)
(* *)
(******************************************************************************)
Notation "a %:R " := (INR a) (at level 10).
Open Scope R_scope.
Ltac tlia := try lia.
Ltac tlra := try lra.
Ltac Rcotpatch :=
change (NormedModule.sort R_AbsRing R_NormedModule) with R; repeat
match goal with
|- context[scal ?X ?Y] => change (scal X Y) with (X * Y)
end.
Lemma pow_lt_1 x n: Rabs x < 1 -> Rabs (x ^ (S n)) < 1.
Section Hole.
Definition hole (k : nat) (a : nat -> R) (n : nat) :=
if Nat.eqb (n mod k) 0 then a (n / k)%nat else 0.
Lemma sum_f_R0_holes_small a g k n
(trm := fun i => hole k a i * g i) l :
(S l < k)%nat
-> sum_f_R0 (fun i => trm (S(k * n) + i)%nat) l = 0.
Lemma sumf_f_R0_add_r f m n : (n <> 0)%nat ->
sum_f_R0 f (m + n) =
sum_f_R0 f m + sum_f_R0 (fun i => f (S m + i)%nat) (pred n).
Lemma sum_f_R0_holes k a g n :
(k <> 0)%nat -> sum_f_R0 (fun i => hole k a i * g i) (k * n) =
sum_f_R0 (fun i => a i * g (k * i)%nat) n.
Lemma fill_holes : forall k a x, (k <> 0)%nat ->
ex_pseries a (x ^ k)->
PSeries (hole k a) x = Series (fun n => a n * x ^ (k * n)).
End Hole.
Section Radius.
Variable a : nat -> R.
Hypothesis aP : forall n, 0 <= a n <= 1.
Lemma PS_cv: forall x, 0 <= x -> x < 1 ->
Rbar.Rbar_lt (Rbar.Finite (Rabs x)) (CV_radius a) .
End Radius.
Section Plouffe.
Lemma PSerie_pow_add x k f a:
Series (fun i => (a i ) * x ^ (k + f i)) =
x ^ k * Series (fun i => (a i) * x ^ f i).
Lemma PSeries_hole x a d k : 0 <= x < 1 ->
Series (fun i : nat => (a * x ^ (d + (S k) * i))) =
PSeries (PS_incr_n (hole (S k) (fun _ => a)) d) x.
Lemma RInt_Series_hole n r: 0 <= r -> r < 1 -> (0 < n)%nat ->
RInt (fun x => Series (fun i => x ^ (n - 1 + 8 * i))) 0 r =
PSeries (PS_Int (PS_incr_n (hole 8 (fun _ => 1)) (n - 1))) r.
Lemma is_RInt_pow n b: is_RInt (fun x => (x ^ n)) 0 b (b ^ (S n)/ (S n)%:R).
Lemma RInt_pow n b : RInt (fun x => x ^ n) 0 b = b ^ (S n)/ (S n)%:R.
Lemma Series_16_Rint n (f := fun i => / (16 ^ i * (8 * i + n)%:R)) :
(0 < n)%nat ->
Series f =
sqrt 2 ^ n *
RInt (fun x => x ^ (n - 1) / (1 - x ^ 8)) 0 (/(sqrt 2)).
(* have hr2:= sqrt2_neq_0.*)
Lemma Pos_x2 x : 0 < x ^ 2 - 2 * x + 2.
Lemma Derive_x2 x :
Derive (fun x => -2 * ln (x ^ 2 - 2 * x + 2)) x =
(4 - 4 * x) / (x ^ 2 - 2 * x + 2).
Lemma Continuity_x2 x :
continuity_pt (fun x0 : R => (4 - 4 * x0) / (x0 ^ 2 - 2 * x0 + 2)) x.
Lemma Rint_x2 :
RInt (fun x => (4 - 4 * x) / (x ^ 2 - 2 * x + 2)) 0 1 = 2 * (ln 2 - ln 1).
Lemma Sqrt_x2 x : -sqrt 2 < x < sqrt 2 -> 0 < 2 - x ^ 2.
Lemma Sqrt2_m1 : 0 < sqrt 2 - 1.
Lemma X_between x : 0 <= x <= 1 -> - sqrt 2 < x < sqrt 2.
Lemma Derive_ln_x2 x :
-sqrt 2 < x < sqrt 2 ->
Derive (fun x => 2 * ln (2 - x ^ 2)) x = (4 * x) / (x ^ 2 - 2).
Lemma Continuity_inv_x2 x (Hx : 0 <= x <= 1) :
continuity_pt (fun x0 : R => 4 * x0 / (x0 ^ 2 - 2)) x.
Lemma Rint_inv_x2 :
RInt (fun x => (4 * x) / (x ^ 2 - 2)) 0 1 = 2 * (ln 1 - ln 2).
Lemma Derive_4atan x :
Derive (fun x => 4 * atan (x - 1)) x = 4 / (x ^ 2 - 2 * x + 2).
Lemma Continuity_4_inv_x2 x :
continuity_pt (fun x0 : R => 4 / (x0 ^ 2 - 2 * x0 + 2)) x.
Lemma Rint_4x2 : RInt (fun x => 4 / (x ^ 2 - 2 * x + 2)) 0 1 = PI.
Lemma Rint_x8_PI (f := fun (n : nat) x =>
sqrt 2 ^ (1 + n) * (x ^ n / (1 - x ^ 8))) :
RInt (fun x => (4 * f 0%nat x - 2 * f 3%nat x - f 4%nat x - f 5%nat x))
0 (/sqrt 2) = PI.
Lemma ex_RInt_x8 n (f := fun x => sqrt 2 ^ (1 + n) * (x ^ n / (1 - x ^ 8))) :
ex_RInt f 0 (/ sqrt 2).
Lemma Plouffe_PI (f := fun n i => / (16 ^ i * (8 * i + n)%:R)) :
PI = 4 * (Series (f 1%nat)) - 2 * (Series (f 4%nat))
- (Series (f 5%nat)) - (Series (f 6%nat)).
End Plouffe.