(* ========================================================================= *)
(* The real and complex gamma functions and Euler-Mascheroni constant. *)
(* ========================================================================= *)
needs "Multivariate/cauchy.ml";;
(* ------------------------------------------------------------------------- *)
(* Euler-Macheroni constant. *)
(* ------------------------------------------------------------------------- *)
(* ------------------------------------------------------------------------- *)
(* Simple-minded estimation of gamma using Euler-Maclaurin summation. *)
(* ------------------------------------------------------------------------- *)
let LOG2_APPROX_40 = prove
(`abs(log(&2) - &381061692393 / &549755813888) <= inv(&2 pow 40)`,
MP_TAC(SPECL [`41`; `Cx(--inv(&2))`]
TAYLOR_CLOG) THEN
SIMP_TAC[GSYM
CX_DIV; GSYM
CX_POW; GSYM
CX_NEG; GSYM
CX_ADD; GSYM
CX_MUL;
VSUM_CX;
COMPLEX_NORM_CX; GSYM
CX_LOG; GSYM
CX_SUB;
REAL_ARITH `&0 < &1 + --inv(&2)`] THEN
CONV_TAC NUM_REDUCE_CONV THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
REWRITE_TAC[REAL_ARITH `a * b / c:real = a / c * b`] THEN
CONV_TAC(ONCE_DEPTH_CONV EXPAND_SUM_CONV) THEN
CONV_TAC NUM_REDUCE_CONV THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
REWRITE_TAC[
real_div; REAL_MUL_LID] THEN
SIMP_TAC[
LOG_INV; REAL_ARITH `&0 < &2`] THEN REAL_ARITH_TAC);;
(* ------------------------------------------------------------------------- *)
(* Start with the log-gamma function. It's otherwise quite tedious to repeat *)
(* essentially the same argument when we want the logarithm of the gamma *)
(* function, since we can't just take the usual principal value of log. *)
(* ------------------------------------------------------------------------- *)
let lgamma = new_definition
`lgamma z = lim sequentially
(\n. z * clog(Cx(&n)) - clog z -
vsum(1..n) (\m. clog((Cx(&m) + z) / Cx(&m))))`;;
let LGAMMA,COMPLEX_DIFFERENTIABLE_AT_LGAMMA = (CONJ_PAIR o prove)
(`(!z. (!n. ~(z + Cx(&n) = Cx(&0)))
==> ((\n. z * clog(Cx(&n)) - clog z -
vsum(1..n) (\m. clog((Cx(&m) + z) / Cx(&m))))
--> lgamma(z)) sequentially) /\
(!z. (Im z = &0 ==> &0 < Re z) ==> lgamma complex_differentiable at z)`,
SUBGOAL_THEN `open {z | !n. ~(z + Cx(&n) = Cx(&0))}` ASSUME_TAC THENL
[REWRITE_TAC[SET_RULE `{z | !n. P n z} = UNIV DIFF {z | ?n. ~P n z}`] THEN
REWRITE_TAC[GSYM closed] THEN MATCH_MP_TAC DISCRETE_IMP_CLOSED THEN
EXISTS_TAC `&1` THEN REWRITE_TAC[REAL_LT_01; IMP_CONJ] THEN
REWRITE_TAC[IN_ELIM_THM; LEFT_IMP_EXISTS_THM] THEN
SIMP_TAC[COMPLEX_RING `x + y = Cx(&0) <=> x = --y`] THEN
REWRITE_TAC[COMPLEX_RING `--x - --y:complex = y - x`] THEN
REWRITE_TAC[COMPLEX_EQ_NEG2; CX_INJ; GSYM CX_SUB; COMPLEX_NORM_CX] THEN
SIMP_TAC[GSYM REAL_EQ_INTEGERS; INTEGER_CLOSED];
ALL_TAC] THEN
SUBGOAL_THEN
`!y. (!n. ~(y + Cx(&n) = Cx(&0)))
==> ?d l. &0 < d /\
!e. &0 < e
==> ?N. !n z. N <= n /\ z IN cball(y,d)
==> dist(z * clog(Cx(&n)) -
vsum(1..n)
(\m. clog((Cx(&m) + z) / Cx(&m))),
l z) < e`
MP_TAC THENL
[REPEAT STRIP_TAC THEN
FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [OPEN_CONTAINS_BALL]) THEN
REWRITE_TAC[OPEN_CONTAINS_CBALL] THEN
DISCH_THEN(MP_TAC o SPEC `y:complex`) THEN
ASM_REWRITE_TAC[IN_ELIM_THM] THEN
MATCH_MP_TAC MONO_EXISTS THEN X_GEN_TAC `d:real` THEN
STRIP_TAC THEN ASM_REWRITE_TAC[] THEN
REWRITE_TAC[UNIFORMLY_CONVERGENT_EQ_CAUCHY_ALT] THEN
X_GEN_TAC `e:real` THEN DISCH_TAC THEN
SUBGOAL_THEN
`summable (from 1)
(\m. Cx(&2 * ((norm(y:complex) + d) +
(norm(y) + d) pow 2)) / Cx(&m) pow 2)`
MP_TAC THENL
[REWRITE_TAC[complex_div; COMPLEX_MUL_ASSOC] THEN
MATCH_MP_TAC SUMMABLE_COMPLEX_LMUL THEN
MATCH_MP_TAC SUMMABLE_ZETA_INTEGER THEN REWRITE_TAC[LE_REFL];
ALL_TAC] THEN
REWRITE_TAC[summable; SERIES_CAUCHY; GE] THEN
DISCH_THEN(MP_TAC o SPEC `e:real`) THEN ASM_REWRITE_TAC[] THEN
DISCH_THEN(X_CHOOSE_THEN `M:num` (LABEL_TAC "M")) THEN
MP_TAC(SPEC `&2 * (norm(y:complex) + d) + &1` REAL_ARCH_SIMPLE) THEN
DISCH_THEN(X_CHOOSE_THEN `N:num` (LABEL_TAC "N")) THEN
EXISTS_TAC `MAX (MAX 1 2) (MAX M N)` THEN
MAP_EVERY X_GEN_TAC [`m:num`; `n:num`; `z:complex`] THEN
REWRITE_TAC[GE; ARITH_RULE `MAX a b <= c <=> a <= c /\ b <= c`] THEN
STRIP_TAC THEN FIRST_X_ASSUM(MP_TAC o SPECL [`m + 1`; `n:num`]) THEN
ANTS_TAC THENL [ASM_ARITH_TAC; ALL_TAC] THEN
REWRITE_TAC[FROM_INTER_NUMSEG_MAX; ARITH_RULE `MAX 1 (m + 1) = m + 1`] THEN
REWRITE_TAC[dist] THEN
SUBGOAL_THEN
`!n. 1 <= n
==> z * clog(Cx(&n)) - vsum(1..n) (\m. clog((Cx(&m) + z) / Cx(&m))) =
vsum(1..n) (\m. z * (clog(Cx(&(m + 1) - &1)) -
clog(Cx(&m - &1))) -
clog((Cx(&m) + z) / Cx(&m))) +
z * clog(Cx(&0))`
(fun th -> ASM_SIMP_TAC[th])
THENL
[REWRITE_TAC[VSUM_SUB_NUMSEG] THEN
ASM_SIMP_TAC[VSUM_COMPLEX_LMUL; FINITE_NUMSEG; VSUM_DIFFS_ALT] THEN
REWRITE_TAC[GSYM REAL_OF_NUM_ADD; REAL_ARITH `(x + &1) - &1 = x`;
REAL_SUB_REFL] THEN
REPEAT STRIP_TAC THEN CONV_TAC COMPLEX_RING;
ALL_TAC] THEN
MATCH_MP_TAC(REWRITE_RULE[IMP_CONJ] REAL_LET_TRANS) THEN
SUBGOAL_THEN `1 <= m + 1 /\ m <= n` MP_TAC THENL
[ASM_ARITH_TAC;
DISCH_THEN(fun th -> REWRITE_TAC[GSYM(MATCH_MP VSUM_COMBINE_R th)])] THEN
REWRITE_TAC[COMPLEX_RING `(x + a) - ((x + y) + a):complex = --y`] THEN
REWRITE_TAC[NORM_NEG] THEN MATCH_MP_TAC COMPLEX_NORM_VSUM_BOUND THEN
REWRITE_TAC[FINITE_NUMSEG; IN_NUMSEG] THEN X_GEN_TAC `k:num` THEN
STRIP_TAC THEN REWRITE_TAC[GSYM CX_POW; GSYM CX_DIV; REAL_CX; RE_CX] THEN
REWRITE_TAC[COMPLEX_RING `z * (a - b) - c:complex = --(z * (b - a) + c)`;
NORM_NEG; GSYM REAL_OF_NUM_ADD; REAL_ARITH `(x + &1) - &1 = x`] THEN
SUBGOAL_THEN `1 <= k /\ 1 < k /\ 2 <= k` STRIP_ASSUME_TAC THENL
[ASM_ARITH_TAC; ALL_TAC] THEN
SUBGOAL_THEN
`clog(Cx(&k - &1)) - clog(Cx(&k)) = clog(Cx(&1) - inv(Cx(&k)))`
SUBST1_TAC THENL
[ASM_SIMP_TAC[GSYM CX_LOG; REAL_OF_NUM_LE; REAL_SUB_LT; REAL_INV_LT_1;
REAL_ARITH `&2 <= x ==> &0 < x /\ &1 < x /\ &0 < x - &1`;
GSYM CX_SUB; GSYM CX_INV; GSYM LOG_DIV] THEN
AP_TERM_TAC THEN AP_TERM_TAC THEN UNDISCH_TAC `2 <= k` THEN
REWRITE_TAC[GSYM REAL_OF_NUM_LE] THEN CONV_TAC REAL_FIELD;
ALL_TAC] THEN
MP_TAC(ISPECL [`1`; `z / Cx(&k)`] TAYLOR_CLOG) THEN
MP_TAC(ISPECL [`1`; `--inv(Cx(&k))`] TAYLOR_CLOG) THEN
REWRITE_TAC[VSUM_SING_NUMSEG] THEN CONV_TAC NUM_REDUCE_CONV THEN
REWRITE_TAC[GSYM VECTOR_SUB; NORM_NEG] THEN
REWRITE_TAC[COMPLEX_NORM_INV; COMPLEX_NORM_DIV; COMPLEX_NORM_CX] THEN
REWRITE_TAC[COMPLEX_POW_NEG; ARITH; REAL_ABS_NUM; COMPLEX_POW_ONE] THEN
ASM_SIMP_TAC[REAL_INV_LT_1; REAL_OF_NUM_LT; COMPLEX_DIV_1] THEN
ASM_SIMP_TAC[REAL_LT_LDIV_EQ; REAL_OF_NUM_LT; LE_1; COMPLEX_POW_1] THEN
REWRITE_TAC[REAL_MUL_LID; COMPLEX_MUL_LID] THEN
DISCH_THEN(MP_TAC o SPEC `norm(z:complex)` o MATCH_MP
(REWRITE_RULE[IMP_CONJ_ALT] REAL_LE_LMUL)) THEN
REWRITE_TAC[NORM_POS_LE; GSYM COMPLEX_NORM_MUL] THEN
SUBGOAL_THEN `norm(z:complex) < &k` ASSUME_TAC THENL
[FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [IN_CBALL]) THEN
FIRST_X_ASSUM(K ALL_TAC o check (is_forall o concl)) THEN
REPEAT(POP_ASSUM MP_TAC) THEN
REWRITE_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_ADD] THEN
CONV_TAC NORM_ARITH;
ASM_REWRITE_TAC[COMPLEX_SUB_LDISTRIB]] THEN
ASM_SIMP_TAC[CX_INJ; REAL_OF_NUM_EQ; LE_1; COMPLEX_FIELD
`~(k = Cx(&0)) ==> (k + z) / k = Cx(&1) + z / k`] THEN
MATCH_MP_TAC(NORM_ARITH
`x' = --y' /\ d + e <= f
==> norm(x - x') <= d ==> norm(y - y') <= e
==> norm(x + y) <= f`) THEN
REWRITE_TAC[complex_div; COMPLEX_MUL_LNEG; COMPLEX_MUL_RNEG] THEN
REWRITE_TAC[REAL_POW_DIV; REAL_POW_INV; REAL_ARITH
`n * inv k / d + n pow 2 / k / e <= (&2 * (x + x pow 2)) / k <=>
(n * (&1 / d + n / e)) / k <= (x * (&2 + &2 * x)) / k`] THEN
ASM_SIMP_TAC[REAL_LE_DIV2_EQ; REAL_OF_NUM_LT; LE_1; REAL_POW_LT] THEN
MATCH_MP_TAC REAL_LE_MUL2 THEN REWRITE_TAC[NORM_POS_LE] THEN
SUBGOAL_THEN `norm(z:complex) <= norm(y:complex) + d` ASSUME_TAC THENL
[FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [IN_CBALL]) THEN
CONV_TAC NORM_ARITH;
ALL_TAC] THEN
REPEAT CONJ_TAC THENL
[FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [IN_CBALL]) THEN
CONV_TAC NORM_ARITH;
MATCH_MP_TAC REAL_LE_ADD THEN CONJ_TAC THEN
MATCH_MP_TAC REAL_LE_DIV THEN
REWRITE_TAC[REAL_SUB_LE; NORM_POS_LE; REAL_POS] THEN
REWRITE_TAC[REAL_ARITH `inv x = &1 / x`] THEN
ASM_SIMP_TAC[REAL_LE_LDIV_EQ; REAL_OF_NUM_LT; LE_1];
REWRITE_TAC[REAL_ARITH `&2 + &2 * x = &1 * &2 + x * &2`] THEN
ONCE_REWRITE_TAC[real_div] THEN MATCH_MP_TAC REAL_LE_ADD2 THEN
CONJ_TAC THEN MATCH_MP_TAC REAL_LE_MUL2 THEN
ASM_REWRITE_TAC[NORM_POS_LE; REAL_POS; REAL_LE_REFL; REAL_LE_INV_EQ] THEN
REWRITE_TAC[REAL_SUB_LE] THEN
REWRITE_TAC[REAL_ARITH `inv x = &1 / x`] THEN
ASM_SIMP_TAC[REAL_LE_LDIV_EQ; REAL_OF_NUM_LT; LE_1] THEN
ASM_REWRITE_TAC[REAL_MUL_LID; REAL_OF_NUM_LE] THEN
TRY(CONJ_TAC THENL
[RULE_ASSUM_TAC(REWRITE_RULE
[GSYM REAL_OF_NUM_LT; GSYM REAL_OF_NUM_ADD]) THEN
ASM_REAL_ARITH_TAC; ALL_TAC]) THEN
REWRITE_TAC[real_div; REAL_MUL_LID] THEN
GEN_REWRITE_TAC RAND_CONV [GSYM REAL_INV_INV] THEN
MATCH_MP_TAC REAL_LE_INV2 THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
REWRITE_TAC[REAL_ARITH `&1 / &2 <= &1 - x <=> x <= &1 / &2`] THEN
REWRITE_TAC[real_div; REAL_MUL_LID] THEN REWRITE_TAC[GSYM real_div] THEN
REWRITE_TAC[REAL_ARITH `inv x = &1 / x`] THEN
ASM_SIMP_TAC[REAL_LE_LDIV_EQ; REAL_OF_NUM_LT; LE_1]] THEN
FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [IN_CBALL]) THEN
FIRST_X_ASSUM(K ALL_TAC o check (is_forall o concl)) THEN
REPEAT(POP_ASSUM MP_TAC) THEN
REWRITE_TAC[GSYM REAL_OF_NUM_LE; GSYM REAL_OF_NUM_ADD] THEN
CONV_TAC NORM_ARITH;
GEN_REWRITE_TAC (LAND_CONV o REDEPTH_CONV) [RIGHT_IMP_EXISTS_THM] THEN
REWRITE_TAC[SKOLEM_THM; LEFT_IMP_EXISTS_THM] THEN MAP_EVERY X_GEN_TAC
[`dd:complex->real`; `ll:complex->complex->complex`] THEN
DISCH_THEN(LABEL_TAC "*")] THEN
SUBGOAL_THEN
`!z. (!n. ~(z + Cx(&n) = Cx(&0)))
==> ((\n. z * clog(Cx(&n)) -
vsum (1..n) (\m. clog((Cx(&m) + z) / Cx(&m)))) --> ll z z)
sequentially`
ASSUME_TAC THENL
[X_GEN_TAC `z:complex` THEN DISCH_TAC THEN
FIRST_X_ASSUM(MP_TAC o SPEC `z:complex`) THEN
ASM_REWRITE_TAC[LIM_SEQUENTIALLY; GSYM SKOLEM_THM] THEN
MESON_TAC[CENTRE_IN_CBALL; REAL_LT_IMP_LE];
ALL_TAC] THEN
MATCH_MP_TAC(TAUT `p /\ (p ==> q) ==> p /\ q`) THEN CONJ_TAC THENL
[REPEAT STRIP_TAC THEN
REWRITE_TAC[lgamma; lim] THEN CONV_TAC SELECT_CONV THEN
EXISTS_TAC `(ll:complex->complex->complex) z z - clog z` THEN
ONCE_REWRITE_TAC[COMPLEX_RING `w - z - v:complex = (w - v) - z`] THEN
MATCH_MP_TAC LIM_SUB THEN REWRITE_TAC[LIM_CONST] THEN ASM_MESON_TAC[];
DISCH_TAC] THEN
X_GEN_TAC `z:complex` THEN DISCH_TAC THEN
SUBGOAL_THEN `!n. ~(z + Cx(&n) = Cx(&0))` ASSUME_TAC THENL
[GEN_TAC THEN
REWRITE_TAC[COMPLEX_RING `z + x = Cx(&0) <=> z = --x`] THEN
DISCH_THEN SUBST_ALL_TAC THEN
FIRST_X_ASSUM(MP_TAC o check (is_imp o concl)) THEN
REWRITE_TAC[IM_NEG; RE_NEG; IM_CX; RE_CX] THEN REAL_ARITH_TAC;
ALL_TAC] THEN
MATCH_MP_TAC COMPLEX_DIFFERENTIABLE_TRANSFORM_AT THEN
EXISTS_TAC `\z. (ll:complex->complex->complex) z z - clog z` THEN
FIRST_X_ASSUM(MP_TAC o GEN_REWRITE_RULE I [OPEN_CONTAINS_BALL]) THEN
REWRITE_TAC[OPEN_CONTAINS_BALL] THEN
DISCH_THEN(MP_TAC o SPEC `z:complex`) THEN
ASM_REWRITE_TAC[IN_ELIM_THM; SUBSET; IN_BALL] THEN
MATCH_MP_TAC MONO_EXISTS THEN X_GEN_TAC `e:real` THEN STRIP_TAC THEN
ASM_REWRITE_TAC[] THEN CONJ_TAC THENL
[X_GEN_TAC `w:complex` THEN ONCE_REWRITE_TAC[DIST_SYM] THEN DISCH_TAC THEN
MATCH_MP_TAC(ISPEC `sequentially` LIM_UNIQUE) THEN
REWRITE_TAC[TRIVIAL_LIMIT_SEQUENTIALLY] THEN
EXISTS_TAC `\n. w * clog(Cx(&n)) - clog w -
vsum(1..n) (\m. clog((Cx(&m) + w) / Cx(&m)))` THEN
ASM_SIMP_TAC[] THEN
ONCE_REWRITE_TAC[COMPLEX_RING `w - z - v:complex = (w - v) - z`] THEN
MATCH_MP_TAC LIM_SUB THEN REWRITE_TAC[LIM_CONST] THEN ASM_MESON_TAC[];
ALL_TAC] THEN
MATCH_MP_TAC COMPLEX_DIFFERENTIABLE_SUB THEN
ASM_SIMP_TAC[COMPLEX_DIFFERENTIABLE_AT_CLOG] THEN
MATCH_MP_TAC COMPLEX_DIFFERENTIABLE_TRANSFORM_AT THEN
EXISTS_TAC `(ll:complex->complex->complex) z` THEN
EXISTS_TAC `min e (dd(z:complex))` THEN
ASM_SIMP_TAC[REAL_LT_MIN] THEN CONJ_TAC THENL
[X_GEN_TAC `w:complex` THEN ONCE_REWRITE_TAC[DIST_SYM] THEN
STRIP_TAC THEN
MATCH_MP_TAC(ISPEC `sequentially` LIM_UNIQUE) THEN
EXISTS_TAC `\n. w * clog(Cx(&n)) -
vsum(1..n) (\m. clog((Cx(&m) + w) / Cx(&m)))` THEN
REWRITE_TAC[TRIVIAL_LIMIT_SEQUENTIALLY; LIM_SEQUENTIALLY] THEN
CONJ_TAC THEN X_GEN_TAC `r:real` THEN STRIP_TAC THENL
[REMOVE_THEN "*" (MP_TAC o SPEC `z:complex`);
REMOVE_THEN "*" (MP_TAC o SPEC `w:complex`)] THEN
ASM_SIMP_TAC[] THEN
REWRITE_TAC[GSYM SKOLEM_THM; RIGHT_EXISTS_IMP_THM] THEN
DISCH_THEN(CONJUNCTS_THEN2 ASSUME_TAC (MP_TAC o SPEC `r:real`)) THEN
ASM_REWRITE_TAC[] THEN MATCH_MP_TAC MONO_EXISTS THEN
REPEAT STRIP_TAC THEN FIRST_X_ASSUM MATCH_MP_TAC THEN
ASM_SIMP_TAC[CENTRE_IN_CBALL; REAL_LT_IMP_LE] THEN
ASM_SIMP_TAC[IN_CBALL; REAL_LT_IMP_LE];
ALL_TAC] THEN
SUBGOAL_THEN `open {z | Im z = &0 ==> &0 < Re z}` MP_TAC THENL
[SUBGOAL_THEN
`{z | Im z = &0 ==> &0 < Re z} =
(:complex) DIFF ({z | Im z = &0} INTER {z | Re z <= &0})`
(fun th -> SIMP_TAC[th; CLOSED_HALFSPACE_RE_LE; GSYM closed;
CLOSED_HALFSPACE_IM_EQ; CLOSED_INTER]) THEN
REWRITE_TAC[IN_ELIM_THM; IN_UNIV; IN_DIFF; IN_INTER; EXTENSION] THEN
REAL_ARITH_TAC;
REWRITE_TAC[OPEN_CONTAINS_CBALL; IN_ELIM_THM; SUBSET] THEN
DISCH_THEN(MP_TAC o SPEC `z:complex`) THEN ASM_REWRITE_TAC[] THEN
DISCH_THEN(X_CHOOSE_THEN `r:real` STRIP_ASSUME_TAC)] THEN
SUBGOAL_THEN
`(ll:complex->complex->complex) z continuous_on cball(z,min r (dd z)) /\
ll z holomorphic_on ball(z,min r (dd z))`
MP_TAC THENL
[MATCH_MP_TAC(ISPEC `sequentially` HOLOMORPHIC_UNIFORM_LIMIT) THEN
EXISTS_TAC `\n z. z * clog(Cx(&n)) -
vsum(1..n) (\m. clog((Cx(&m) + z) / Cx(&m)))` THEN
ASM_REWRITE_TAC[TRIVIAL_LIMIT_SEQUENTIALLY; EVENTUALLY_SEQUENTIALLY] THEN
REWRITE_TAC[CBALL_MIN_INTER; IN_INTER] THEN
CONJ_TAC THENL [ALL_TAC; SIMP_TAC[GSYM dist] THEN ASM_MESON_TAC[]] THEN
EXISTS_TAC `1` THEN X_GEN_TAC `n:num` THEN DISCH_TAC THEN
MATCH_MP_TAC(MESON[HOLOMORPHIC_ON_IMP_CONTINUOUS_ON;
HOLOMORPHIC_ON_SUBSET]
`t SUBSET s /\ f holomorphic_on s
==> f continuous_on s /\ f holomorphic_on t`) THEN
REWRITE_TAC[BALL_SUBSET_CBALL; GSYM CBALL_MIN_INTER] THEN
MATCH_MP_TAC HOLOMORPHIC_ON_SUB THEN
SIMP_TAC[HOLOMORPHIC_ON_MUL; HOLOMORPHIC_ON_ID;
HOLOMORPHIC_ON_CONST] THEN
MATCH_MP_TAC HOLOMORPHIC_ON_VSUM THEN
REWRITE_TAC[FINITE_NUMSEG; IN_NUMSEG] THEN
X_GEN_TAC `m:num` THEN STRIP_TAC THEN
GEN_REWRITE_TAC LAND_CONV [GSYM o_DEF] THEN
MATCH_MP_TAC HOLOMORPHIC_ON_COMPOSE THEN
SIMP_TAC[HOLOMORPHIC_ON_ADD; HOLOMORPHIC_ON_ID;
HOLOMORPHIC_ON_CONST; complex_div; HOLOMORPHIC_ON_MUL] THEN
MATCH_MP_TAC HOLOMORPHIC_ON_CLOG THEN
REWRITE_TAC[FORALL_IN_IMAGE; GSYM complex_div; IMP_CONJ] THEN
ASM_SIMP_TAC[RE_DIV_CX; IM_DIV_CX; REAL_DIV_EQ_0; RE_ADD; IM_ADD] THEN
X_GEN_TAC `w:complex` THEN DISCH_TAC THEN
ASM_SIMP_TAC[REAL_LT_IMP_NZ; REAL_OF_NUM_LT; LE_1] THEN
REWRITE_TAC[IM_CX; RE_CX; REAL_ADD_LID] THEN DISCH_TAC THEN
MATCH_MP_TAC REAL_LT_DIV THEN
ASM_SIMP_TAC[REAL_LT_IMP_NZ; REAL_OF_NUM_LT; LE_1] THEN
MATCH_MP_TAC(REAL_ARITH `&0 < x ==> &0 < &m + x`) THEN
RULE_ASSUM_TAC(REWRITE_RULE[CBALL_MIN_INTER; IN_INTER]) THEN
ASM_MESON_TAC[];
SIMP_TAC[HOLOMORPHIC_ON_OPEN; OPEN_BALL; complex_differentiable] THEN
DISCH_THEN(MATCH_MP_TAC o CONJUNCT2) THEN
ASM_SIMP_TAC[CENTRE_IN_BALL; REAL_LT_MIN]]);;
let LGAMMA_ALT = prove
(`!z. (!n. ~(z + Cx(&n) = Cx(&0)))
==> ((\n. (z * clog(Cx(&n)) + clog(Cx(&(
FACT n)))) -
vsum(0..n) (\m. clog(Cx(&m) + z)))
--> lgamma(z)) sequentially`,
let LGAMMA_ALT2 = prove
(`!z. (!n. ~(z + Cx(&n) = Cx(&0)))
==> ((\n. vsum(1..n) (\m. z * clog(Cx(&1) + Cx(&1) / Cx(&m)) -
clog(Cx(&1) + z / Cx(&m))) -
clog(z))
--> lgamma(z)) sequentially`,
(* ------------------------------------------------------------------------- *)
(* The complex gamma function (defined using the Gauss/Euler product). *)
(* Note that this totalizes it with the value zero at the poles. *)
(* ------------------------------------------------------------------------- *)
let [CGAMMA;CGAMMA_EQ_0;CGAMMA_LGAMMA] = (CONJUNCTS o prove)
(`(!z. ((\n. (Cx(&n) cpow z * Cx(&(FACT n))) / cproduct(0..n) (\m. z + Cx(&m)))
--> cgamma(z)) sequentially) /\
(!z. cgamma(z) = Cx(&0) <=> ?n. z + Cx(&n) = Cx(&0)) /\
(!z. cgamma(z) = if ?n. z + Cx(&n) = Cx(&0) then Cx(&0)
else cexp(lgamma z))`,
REWRITE_TAC[AND_FORALL_THM] THEN X_GEN_TAC `y:complex` THEN
ASM_CASES_TAC `?n. y + Cx(&n) = Cx(&0)` THENL
[ASM_REWRITE_TAC[GSYM NOT_EXISTS_THM] THEN
FIRST_X_ASSUM(X_CHOOSE_TAC `N:num`) THEN
REWRITE_TAC[cgamma; lim] THEN MATCH_MP_TAC(MESON[LIM_UNIQUE]
`~trivial_limit net /\ (f --> a) net
==> (f --> @a. (f --> a) net) net /\ (@a. (f --> a) net) = a`) THEN
REWRITE_TAC[TRIVIAL_LIMIT_SEQUENTIALLY] THEN
MATCH_MP_TAC LIM_EVENTUALLY THEN REWRITE_TAC[EVENTUALLY_SEQUENTIALLY] THEN
EXISTS_TAC `N:num` THEN X_GEN_TAC `n:num` THEN REWRITE_TAC[GE] THEN
DISCH_TAC THEN
SIMP_TAC[COMPLEX_ENTIRE; COMPLEX_DIV_EQ_0; CPRODUCT_EQ_0; FINITE_NUMSEG;
IN_NUMSEG; LE_0] THEN
ASM_MESON_TAC[LE_REFL];
ASM_REWRITE_TAC[] THEN RULE_ASSUM_TAC(REWRITE_RULE[NOT_EXISTS_THM])] THEN
SUBGOAL_THEN
`((\n. (Cx(&n) cpow y * Cx(&(FACT n))) / cproduct(0..n) (\m. y + Cx(&m)))
--> cexp(lgamma y)) sequentially`
ASSUME_TAC THENL
[MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC `\n. cexp(y * clog(Cx(&n)) - clog y -
vsum(1..n) (\m. clog((Cx(&m) + y) / Cx(&m))))` THEN
CONJ_TAC THENL
[REWRITE_TAC[EVENTUALLY_SEQUENTIALLY] THEN EXISTS_TAC `1` THEN
X_GEN_TAC `n:num` THEN DISCH_TAC THEN
ASM_SIMP_TAC[CEXP_SUB; cpow; CX_INJ; REAL_OF_NUM_EQ; LE_1] THEN
REWRITE_TAC[complex_div; GSYM COMPLEX_MUL_ASSOC] THEN AP_TERM_TAC THEN
SIMP_TAC[GSYM NPRODUCT_FACT; REAL_OF_NUM_NPRODUCT; FINITE_NUMSEG] THEN
SIMP_TAC[CX_PRODUCT; FINITE_NUMSEG; GSYM CPRODUCT_INV] THEN
SIMP_TAC[CPRODUCT_CLAUSES_LEFT; LE_0] THEN
GEN_REWRITE_TAC RAND_CONV [COMPLEX_RING
`a * b * c:complex = b * a * c`] THEN
REWRITE_TAC[COMPLEX_ADD_RID] THEN BINOP_TAC THENL
[ASM_MESON_TAC[COMPLEX_ADD_RID; CEXP_CLOG]; ALL_TAC] THEN
SIMP_TAC[ADD_CLAUSES; GSYM CPRODUCT_MUL; FINITE_NUMSEG] THEN
REWRITE_TAC[GSYM CEXP_NEG; GSYM VSUM_NEG] THEN
SIMP_TAC[CEXP_VSUM; FINITE_NUMSEG] THEN MATCH_MP_TAC CPRODUCT_EQ THEN
REWRITE_TAC[IN_NUMSEG] THEN X_GEN_TAC `m:num` THEN STRIP_TAC THEN
REWRITE_TAC[CEXP_NEG] THEN
GEN_REWRITE_TAC (RAND_CONV o LAND_CONV) [GSYM COMPLEX_INV_INV] THEN
REWRITE_TAC[GSYM COMPLEX_INV_MUL] THEN AP_TERM_TAC THEN
GEN_REWRITE_TAC RAND_CONV [COMPLEX_MUL_SYM] THEN
GEN_REWRITE_TAC (RAND_CONV o ONCE_DEPTH_CONV) [COMPLEX_ADD_SYM] THEN
MATCH_MP_TAC CEXP_CLOG THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN
ASM_REWRITE_TAC[COMPLEX_ENTIRE; COMPLEX_INV_EQ_0] THEN
REWRITE_TAC[CX_INJ; REAL_OF_NUM_EQ] THEN ASM_ARITH_TAC;
MATCH_MP_TAC(ISPEC `cexp` LIM_CONTINUOUS_FUNCTION) THEN
ASM_SIMP_TAC[LGAMMA; CONTINUOUS_AT_CEXP]];
MATCH_MP_TAC(TAUT `p /\ (p ==> q) ==> p /\ q`) THEN CONJ_TAC THENL
[REWRITE_TAC[cgamma; lim] THEN CONV_TAC SELECT_CONV THEN ASM_MESON_TAC[];
DISCH_TAC] THEN
MATCH_MP_TAC(TAUT `q /\ (q ==> p) ==> p /\ q`) THEN CONJ_TAC THENL
[MATCH_MP_TAC(ISPEC `sequentially` LIM_UNIQUE) THEN
ASM_MESON_TAC[TRIVIAL_LIMIT_SEQUENTIALLY];
SIMP_TAC[CEXP_NZ]]]);;
let CGAMMA_LEGENDRE_ALT = prove
(`!z. cgamma(z) * cgamma(z + Cx(&1) / Cx(&2)) =
Cx(&2) cpow (Cx(&1) - Cx(&2) * z) *
cgamma(Cx(&1) / Cx(&2)) * cgamma(Cx(&2) * z)`,
let CGAMMA_HALF = prove
(`cgamma(Cx(&1) / Cx(&2)) = Cx(sqrt pi)`,
MP_TAC(SPEC `Cx(&1) / Cx(&2)`
CGAMMA_REFLECTION) THEN
REWRITE_TAC[COMPLEX_RING `Cx(&1) - Cx(&1) / Cx(&2) = Cx(&1) / Cx(&2)`] THEN
REWRITE_TAC[GSYM
CX_DIV; GSYM
CX_MUL; GSYM
CX_SIN] THEN
REWRITE_TAC[REAL_ARITH `x * &1 / &2 = x / &2`;
SIN_PI2;
REAL_DIV_1] THEN
SUBGOAL_THEN `Cx pi = Cx(sqrt pi) pow 2` SUBST1_TAC THENL
[REWRITE_TAC[GSYM
CX_POW] THEN AP_TERM_TAC THEN CONV_TAC SYM_CONV THEN
REWRITE_TAC[
SQRT_POW2;
PI_POS_LE];
REWRITE_TAC[COMPLEX_RING
`a * a:complex = b pow 2 <=> a = b \/ a = --b`] THEN
STRIP_TAC THEN ASM_REWRITE_TAC[] THEN
MP_TAC(SPEC `Cx(&1 / &2)`
RE_POS_CGAMMA_REAL) THEN
ASM_REWRITE_TAC[
REAL_CX;
RE_CX;
RE_NEG] THEN
CONV_TAC REAL_RAT_REDUCE_CONV THEN MATCH_MP_TAC(TAUT `~p ==> p ==> q`) THEN
REWRITE_TAC[REAL_ARITH `~(&0 <= --x) <=> &0 < x`] THEN
MATCH_MP_TAC
SQRT_POS_LT THEN REWRITE_TAC[
PI_POS]]);;
let CGAMMA_LEGENDRE = prove
(`!z. cgamma(z) * cgamma(z + Cx(&1) / Cx(&2)) =
Cx(&2) cpow (Cx(&1) - Cx(&2) * z) * Cx(sqrt pi) * cgamma(Cx(&2) * z)`,
(* ------------------------------------------------------------------------- *)
(* Thw Weierstrass product definition. *)
(* ------------------------------------------------------------------------- *)
(* ------------------------------------------------------------------------- *)
(* Sometimes the reciprocal function is convenient, since it's entire. *)
(* ------------------------------------------------------------------------- *)
let RECIP_CGAMMA = prove
(`!z. ((\n. cproduct(0..n) (\m. z + Cx(&m)) /
(Cx(&n) cpow z * Cx(&(
FACT n)))) --> inv(cgamma z))
sequentially`,
(* ------------------------------------------------------------------------- *)
(* The real gamma function. *)
(* ------------------------------------------------------------------------- *)
let GAMMA = prove
(`!x. ((\n. (&n rpow x * &(
FACT n)) / product(0..n) (\m. x + &m))
---> gamma(x)) sequentially`,
let GAMMA_LEGENDRE = prove
(`!x. gamma(x) * gamma(x + &1 / &2) =
&2 rpow (&1 - &2 * x) * sqrt pi * gamma(&2 * x)`,
(* ------------------------------------------------------------------------- *)
(* Characterization of the real gamma function using log-convexity. *)
(* ------------------------------------------------------------------------- *)
let GAMMA_REAL_LOG_CONVEX_UNIQUE = prove
(`!f:real->real.
f(&1) = &1 /\ (!x. &0 < x ==> f(x + &1) = x * f(x)) /\
(!x. &0 < x ==> &0 < f x) /\ f
real_log_convex_on {x | &0 < x}
==> !x. &0 < x ==> f x = gamma x`,
GEN_TAC THEN STRIP_TAC THEN
SUBGOAL_THEN `!x. &0 < x /\ x <= &1 ==> f x = gamma x` ASSUME_TAC THENL
[ALL_TAC;
SUBGOAL_THEN
`!y. &0 < y /\ y <= &1 ==> !n. f(&n + y) = gamma(&n + y)`
ASSUME_TAC THENL
[GEN_TAC THEN STRIP_TAC THEN
INDUCT_TAC THEN ASM_SIMP_TAC[REAL_ADD_LID] THEN
REWRITE_TAC[GSYM
REAL_OF_NUM_SUC] THEN
REWRITE_TAC[REAL_ARITH `(n + &1) + x = (n + x) + &1`] THEN
ASM_SIMP_TAC[
GAMMA_RECURRENCE;
REAL_LET_ADD;
REAL_POS] THEN
ASM_REAL_ARITH_TAC;
X_GEN_TAC `x:real` THEN DISCH_TAC THEN
MP_TAC(SPEC `x:real`
FLOOR_POS) THEN ASM_SIMP_TAC[
REAL_LT_IMP_LE] THEN
DISCH_THEN(X_CHOOSE_TAC `n:num`) THEN
MP_TAC(SPEC `x:real`
FLOOR_FRAC) THEN ASM_REWRITE_TAC[] THEN
ASM_CASES_TAC `frac x = &0` THEN ASM_REWRITE_TAC[
REAL_LT_LE] THENL
[ASM_CASES_TAC `n = 0` THEN
ASM_SIMP_TAC[
REAL_ADD_RID;
REAL_LT_IMP_NZ] THEN
SUBGOAL_THEN `&(n - 1) + &1 = &n`
(fun th -> ASM_MESON_TAC[th;
REAL_LE_REFL;
REAL_LT_01]) THEN
ASM_SIMP_TAC[GSYM
REAL_OF_NUM_SUB;
LE_1] THEN REAL_ARITH_TAC;
STRIP_TAC THEN FIRST_X_ASSUM(MP_TAC o SPEC `frac x`) THEN
ANTS_TAC THENL [ASM_REAL_ARITH_TAC; ASM_MESON_TAC[]]]]] THEN
X_GEN_TAC `x:real` THEN STRIP_TAC THEN
ASM_CASES_TAC `x = &1` THEN ASM_REWRITE_TAC[
GAMMA_1] THEN
SUBGOAL_THEN `&0 < x` ASSUME_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
MATCH_MP_TAC(REAL_FIELD `&0 < g /\ f / g = &1 ==> f = g`) THEN
ASM_SIMP_TAC[
GAMMA_POS_LT] THEN
MATCH_MP_TAC(ISPEC `sequentially`
REALLIM_UNIQUE) THEN
EXISTS_TAC
`\n. f(x) / ((&n rpow x * &(
FACT n)) / product (0..n) (\m. x + &m))` THEN
REWRITE_TAC[
TRIVIAL_LIMIT_SEQUENTIALLY] THEN
ASM_SIMP_TAC[
REALLIM_DIV;
REALLIM_CONST;
GAMMA_POS_LT;
GAMMA;
REAL_LT_IMP_NZ] THEN
ONCE_REWRITE_TAC[
REALLIM_NULL] THEN
MATCH_MP_TAC
REALLIM_NULL_COMPARISON THEN
EXISTS_TAC `\n. x * inv(&n)` THEN
SIMP_TAC[
REALLIM_NULL_LMUL;
REALLIM_1_OVER_N] THEN
REWRITE_TAC[
EVENTUALLY_SEQUENTIALLY] THEN EXISTS_TAC `2` THEN
SUBGOAL_THEN
`!n. &2 <= &n
==> log(f(&n)) - log(f(&n - &1))
<= (log(f(&n + x)) - log(f(&n))) / x /\
(log(f(&n + x)) - log(f(&n))) / x
<= log(f(&n + &1)) - log(f(&n))`
MP_TAC THENL
[MP_TAC(SPECL [`f:real->real`; `{x | &0 < x}`]
REAL_LOG_CONVEX_ON) THEN
ASM_REWRITE_TAC[
IN_ELIM_THM] THEN ANTS_TAC THENL
[REWRITE_TAC[
is_realinterval;
IN_ELIM_THM] THEN REAL_ARITH_TAC;
DISCH_TAC] THEN
MAP_EVERY (MP_TAC o SPECL [`log o f:real->real`; `{x | &0 < x}`])
[
REAL_CONVEX_ON_LEFT_SECANT;
REAL_CONVEX_ON_RIGHT_SECANT] THEN
ASM_REWRITE_TAC[] THEN
DISCH_THEN(LABEL_TAC "L") THEN DISCH_THEN(LABEL_TAC "R") THEN
REPEAT STRIP_TAC THENL
[USE_THEN "L" (MP_TAC o SPECL [`&n - &1`; `&n + x`; `&n`]) THEN
USE_THEN "R" (MP_TAC o SPECL [`&n - &1`; `&n + x`; `&n`]) THEN
REWRITE_TAC[
IN_ELIM_THM;
IN_REAL_SEGMENT;
o_THM] THEN
ANTS_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
REWRITE_TAC[REAL_ARITH `abs(x - (x - &1)) = &1`;
REAL_DIV_1] THEN
DISCH_TAC THEN ANTS_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
ASM_SIMP_TAC[REAL_ARITH `&0 < x ==> abs((n + x) - n) = x`] THEN
ASM_REAL_ARITH_TAC;
ASM_CASES_TAC `x = &1` THEN
ASM_REWRITE_TAC[
REAL_LE_REFL;
REAL_DIV_1] THEN
USE_THEN "R" (MP_TAC o SPECL [`&n`; `&n + &1`; `&n + x`]) THEN
REWRITE_TAC[
IN_ELIM_THM;
IN_REAL_SEGMENT;
o_THM] THEN
ANTS_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
REWRITE_TAC[REAL_ARITH `abs((n + &1) - n) = &1`;
REAL_DIV_1] THEN
ASM_SIMP_TAC[REAL_ARITH `&0 < x ==> abs((n + x) - n) = x`] THEN
ASM_REAL_ARITH_TAC];
ALL_TAC] THEN
ASM_SIMP_TAC[GSYM
LOG_DIV;
REAL_LT_MUL;
REAL_ARITH `&2 <= x ==> &0 < x /\ &0 < x + &1 /\ &0 < x - &1`;
REAL_FIELD `&0 < x ==> (n * x) / x = n`] THEN
SUBGOAL_THEN
`!n. &2 <= n ==> f(n - &1) = f(n) / (n - &1)` (fun th -> SIMP_TAC[th])
THENL
[REPEAT STRIP_TAC THEN
SUBGOAL_THEN `(f:real->real) n = f((n - &1) + &1)` SUBST1_TAC THENL
[AP_TERM_TAC THEN REAL_ARITH_TAC;
ASM_SIMP_TAC[REAL_ARITH `&2 <= n ==> &0 < n - &1`] THEN
UNDISCH_TAC `&2 <= n` THEN CONV_TAC REAL_FIELD];
ASM_SIMP_TAC[REAL_ARITH `&2 <= x ==> &0 < x`; REAL_FIELD
`&0 < x /\ &2 <= y ==> x / (x / (y - &1)) = y - &1`] THEN
ASM_SIMP_TAC[
REAL_LE_LDIV_EQ;
REAL_LE_RDIV_EQ]] THEN
SIMP_TAC[REAL_OF_NUM_LE;
REAL_OF_NUM_SUB;
ARITH_RULE `2 <= n ==> 1 <= n`] THEN
REWRITE_TAC[
REAL_LE_SUB_LADD;
REAL_LE_SUB_RADD] THEN
SUBGOAL_THEN
`!n. 0 < n ==> f(&n) = &(
FACT(n - 1))`
(fun th -> SIMP_TAC[ARITH_RULE `2 <= n ==> 0 < n /\ 0 < n - 1`; th])
THENL
[INDUCT_TAC THEN REWRITE_TAC[ARITH;
LT_0] THEN
ASM_CASES_TAC `n = 0` THEN
ASM_REWRITE_TAC[REAL_ADD_LID;
FACT; ARITH] THEN
ASM_SIMP_TAC[ARITH_RULE `~(n = 0) ==> SUC n - 1 = SUC(n - 1)`] THEN
ASM_SIMP_TAC[
REAL_OF_NUM_LT; GSYM
REAL_OF_NUM_SUC;
LE_1;
FACT] THEN
ASM_SIMP_TAC[ARITH_RULE `~(n = 0) ==> SUC(n - 1) = n`] THEN
REWRITE_TAC[REAL_OF_NUM_MUL;
MULT_SYM];
ALL_TAC] THEN
GEN_REWRITE_TAC (LAND_CONV o BINDER_CONV o RAND_CONV o BINOP_CONV)
[GSYM
REAL_EXP_MONO_LE] THEN
REWRITE_TAC[GSYM REAL_OF_NUM_LE] THEN
ASM_SIMP_TAC[
REAL_EXP_ADD;
EXP_LOG;
REAL_OF_NUM_LT;
LE_1;
FACT_NZ;
ARITH_RULE `&2 <= &n /\ &0 < x ==> &0 < &n + x`] THEN
REWRITE_TAC[REAL_OF_NUM_LE] THEN
SUBGOAL_THEN
`!n. 0 < n ==> f(&n + x) = product(0..n-1) (\k. x + &k) * f(x)`
(fun th -> SIMP_TAC[ARITH_RULE `2 <= n ==> 0 < n`; th])
THENL
[INDUCT_TAC THEN REWRITE_TAC[ARITH] THEN
DISCH_TAC THEN REWRITE_TAC[GSYM
REAL_OF_NUM_SUC] THEN
ASM_CASES_TAC `n = 0` THENL
[ASM_SIMP_TAC[REAL_ARITH `(&0 + &1) + x = x + &1`; ARITH] THEN
REWRITE_TAC[
PRODUCT_SING_NUMSEG;
REAL_ADD_RID];
ASM_SIMP_TAC[REAL_ARITH `(n + &1) + x = (n + x) + &1`;
REAL_LT_ADD;
REAL_OF_NUM_LT;
LE_1] THEN
ASM_SIMP_TAC[ARITH_RULE `~(n = 0) ==> SUC n - 1 = SUC(n - 1)`] THEN
REWRITE_TAC[
PRODUCT_CLAUSES_NUMSEG;
LE_0] THEN
ASM_SIMP_TAC[ARITH_RULE `~(n = 0) ==> SUC(n - 1) = n`] THEN
REWRITE_TAC[
REAL_MUL_AC;
REAL_ADD_AC]];
DISCH_TAC] THEN
X_GEN_TAC `n:num` THEN DISCH_TAC THEN
MATCH_MP_TAC(REAL_ARITH
`&1 <= x /\ x <= &1 + e ==> abs(x - &1) <= e`) THEN
ASM_SIMP_TAC[REAL_OF_NUM_LE; REAL_FIELD
`&2 <= n ==> &1 + x * inv n = (x + n) / n`] THEN
REWRITE_TAC[
real_div;
REAL_INV_MUL;
REAL_INV_INV] THEN
ONCE_REWRITE_TAC[REAL_ARITH `f * (r * n) * p:real = (p * f) * r * n`] THEN
REWRITE_TAC[GSYM
REAL_INV_MUL] THEN REWRITE_TAC[GSYM
real_div] THEN
SUBGOAL_THEN `&0 < &n rpow x * &(
FACT n)` ASSUME_TAC THENL
[ASM_SIMP_TAC[
REAL_LT_MUL;
REAL_OF_NUM_LT;
LE_1;
FACT_NZ;
RPOW_POS_LT; ARITH_RULE `2 <= x ==> 0 < x`];
ALL_TAC] THEN
CONJ_TAC THENL
[ASM_SIMP_TAC[
REAL_LE_RDIV_EQ; REAL_MUL_LID] THEN
FIRST_X_ASSUM(MP_TAC o SPEC `n + 1`) THEN
ANTS_TAC THENL [ASM_ARITH_TAC; DISCH_THEN(MP_TAC o CONJUNCT1)] THEN
REWRITE_TAC[
ADD_SUB] THEN
ASM_SIMP_TAC[rpow; REAL_OF_NUM_LE; REAL_ARITH `&2 <= x ==> &0 < x`] THEN
REWRITE_TAC[
REAL_MUL_AC];
ASM_SIMP_TAC[
PRODUCT_CLAUSES_RIGHT;
LE_0;
ARITH_RULE `2 <= n ==> 0 < n`] THEN
REWRITE_TAC[
real_div; GSYM REAL_MUL_ASSOC] THEN
GEN_REWRITE_TAC LAND_CONV
[REAL_ARITH `a * b * c * d:real = b * (a * c) * d`] THEN
MATCH_MP_TAC
REAL_LE_LMUL THEN
CONJ_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
ASM_SIMP_TAC[GSYM
real_div;
REAL_LE_LDIV_EQ] THEN
FIRST_X_ASSUM(MP_TAC o SPEC `n:num`) THEN ASM_REWRITE_TAC[] THEN
DISCH_THEN(MP_TAC o CONJUNCT2) THEN
MATCH_MP_TAC(REWRITE_RULE[
IMP_CONJ_ALT]
REAL_LE_TRANS) THEN
SUBGOAL_THEN `
FACT n =
FACT(SUC(n - 1))` SUBST1_TAC THENL
[AP_TERM_TAC THEN ASM_ARITH_TAC; REWRITE_TAC[
FACT]] THEN
ASM_SIMP_TAC[ARITH_RULE `2 <= n ==> SUC(n - 1) = n`] THEN
MATCH_MP_TAC
REAL_EQ_IMP_LE THEN
SUBGOAL_THEN `&0 < &n` MP_TAC THENL
[REWRITE_TAC[
REAL_OF_NUM_LT] THEN ASM_ARITH_TAC;
SIMP_TAC[rpow; GSYM REAL_OF_NUM_MUL;
REAL_MUL_AC] THEN
CONV_TAC REAL_FIELD]]);;
(* ------------------------------------------------------------------------- *)
(* The integral definition, the current usual one and Euler's original one. *)
(* ------------------------------------------------------------------------- *)
let EULER_INTEGRAL = prove
(`!z. &0 < Re z
==> integral {t | &0 <= drop t}
(\t. Cx(drop t) cpow (z - Cx(&1)) / cexp(Cx(drop t))) =
cgamma(z)`,
(* ------------------------------------------------------------------------- *)
(* Stirling's approximation. *)
(* ------------------------------------------------------------------------- *)
let LGAMMA_STIRLING_INTEGRALS_EXIST,LGAMMA_STIRLING = (CONJ_PAIR o prove)
(`(!z n. 1 <= n /\ (Im z = &0 ==> &0 < Re z)
==> (\t. Cx(bernoulli n (frac (drop t))) / (z + Cx(drop t)) pow n)
integrable_on {t | &0 <= drop t}) /\
(!z p.
(Im z = &0 ==> &0 < Re z)
==> lgamma(z) =
((z - Cx(&1) / Cx(&2)) * clog(z) - z + Cx(log(&2 * pi) / &2)) +
vsum(1..p) (\k. Cx(bernoulli (2 * k) (&0) /
(&4 * &k pow 2 - &2 * &k)) / z pow (2 * k - 1)) -
integral {t | &0 <= drop t}
(\t. Cx(bernoulli (2 * p + 1) (frac(drop t))) /
(z + Cx(drop t)) pow (2 * p + 1)) /
Cx(&(2 * p + 1)))`,
let lemma1 = prove
(`!p n z. (Im z = &0 ==> &0 < Re z)
==> (\x. Cx(bernoulli (2 * p + 1) (frac(drop x))) *
Cx(&(
FACT(2 * p))) / (z + Cx(drop x)) pow (2 * p + 1))
integrable_on interval [lift(&0),lift(&n)] /\
vsum(0..n) (\m. clog(Cx(&m) + z)) =
(z + Cx(&n) + Cx(&1) / Cx(&2)) * clog(z + Cx(&n)) -
(z - Cx(&1) / Cx(&2)) * clog(z) - Cx(&n) +
vsum(1..p)
(\k. Cx(bernoulli (2 * k) (&0) / &(
FACT(2 * k))) *
(Cx(&(
FACT(2 * k - 2))) / (z + Cx(&n)) pow (2 * k - 1) -
Cx(&(
FACT(2 * k - 2))) / z pow (2 * k - 1))) +
integral (interval[lift(&0),lift(&n)])
(\x. Cx(bernoulli (2 * p + 1) (frac(drop x))) *
Cx(&(
FACT(2 * p))) / (z + Cx(drop x)) pow (2 * p + 1)) /
Cx(&(
FACT(2 * p + 1)))`,
REPEAT GEN_TAC THEN STRIP_TAC THEN
REWRITE_TAC[COMPLEX_RING `Cx(&n) + z = z + Cx(&n)`] THEN MP_TAC(ISPECL
[`\n w. if n = 0 then (z + w) * clog(z + w) - (z + w)
else if n = 1 then clog(z + w)
else Cx(--(&1) pow n * &(
FACT(n - 2))) /
(z + w) pow (n - 1)`;
`0`; `n:num`; `p:num`]
COMPLEX_EULER_MACLAURIN_ANTIDERIVATIVE) THEN
ASM_REWRITE_TAC[
ADD_EQ_0;
MULT_EQ_0;
MULT_EQ_1] THEN
CONV_TAC NUM_REDUCE_CONV THEN ANTS_TAC THENL
[REWRITE_TAC[
IN_REAL_INTERVAL;
LE_0] THEN
MAP_EVERY X_GEN_TAC [`k:num`; `x:real`] THEN REPEAT STRIP_TAC THEN
SUBGOAL_THEN `~(z + Cx x = Cx(&0))` ASSUME_TAC THENL
[REWRITE_TAC[
COMPLEX_EQ;
IM_ADD;
RE_ADD;
IM_CX;
RE_CX] THEN
ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[
o_THM; ARITH_RULE `k + 1 = 1 <=> k = 0`] THEN
ASM_CASES_TAC `k = 0` THEN ASM_REWRITE_TAC[] THENL
[COMPLEX_DIFF_TAC THEN CONJ_TAC THENL
[REWRITE_TAC[
RE_ADD;
IM_ADD;
RE_CX;
IM_CX] THEN ASM_REAL_ARITH_TAC;
UNDISCH_TAC `~(z + Cx x = Cx(&0))` THEN CONV_TAC COMPLEX_FIELD];
ALL_TAC] THEN
ASM_CASES_TAC `k = 1` THEN ASM_REWRITE_TAC[] THEN COMPLEX_DIFF_TAC THEN
CONV_TAC NUM_REDUCE_CONV THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
REWRITE_TAC[
complex_div;
COMPLEX_ADD_LID;
COMPLEX_MUL_LID;
COMPLEX_POW_1;
IM_ADD;
RE_ADD;
IM_CX;
RE_CX] THEN
(CONJ_TAC ORELSE ASM_REAL_ARITH_TAC) THEN
ASM_REWRITE_TAC[
COMPLEX_POW_EQ_0] THEN
REWRITE_TAC[
REAL_POW_ADD;
REAL_POW_1;
COMPLEX_MUL_LZERO] THEN
REWRITE_TAC[
REAL_MUL_RNEG;
REAL_MUL_RID;
REAL_NEG_NEG;
REAL_MUL_LNEG] THEN
REWRITE_TAC[
COMPLEX_MUL_RID;
CX_NEG;
COMPLEX_POW_POW] THEN
REWRITE_TAC[GSYM
COMPLEX_MUL_ASSOC;
COMPLEX_SUB_LZERO;
COMPLEX_NEG_NEG;
COMPLEX_MUL_LNEG] THEN
REWRITE_TAC[GSYM
complex_div] THEN AP_TERM_TAC THEN
FIRST_X_ASSUM(K ALL_TAC o check (is_imp o concl)) THEN
ASM_SIMP_TAC[GSYM
complex_div;
COMPLEX_DIV_POW2] THEN COND_CASES_TAC THENL
[MATCH_MP_TAC(TAUT `F ==> p`) THEN ASM_ARITH_TAC; ALL_TAC] THEN
ASM_SIMP_TAC[
ADD_SUB; ARITH_RULE
`~(k = 0) /\ ~(k = 1) ==> (k - 1) * 2 - (k - 1 - 1) = k`] THEN
REWRITE_TAC[
COMPLEX_MUL_ASSOC;
complex_div] THEN
AP_THM_TAC THEN AP_TERM_TAC THEN REWRITE_TAC[GSYM
CX_MUL] THEN
AP_TERM_TAC THEN REWRITE_TAC[GSYM REAL_MUL_ASSOC] THEN AP_TERM_TAC THEN
ASM_SIMP_TAC[ARITH_RULE `~(k = 0) ==> (k + 1) - 2 = k - 1`] THEN
ASM_SIMP_TAC[ARITH_RULE
`~(k = 0) /\ ~(k = 1) ==> k - 1 = SUC(k - 2)`] THEN
REWRITE_TAC[
FACT; REAL_OF_NUM_MUL] THEN REWRITE_TAC[
MULT_SYM];
REWRITE_TAC[ARITH_RULE `~(2 * p + 2 = 1)`] THEN
DISCH_THEN(CONJUNCTS_THEN2 MP_TAC SUBST1_TAC)] THEN
SIMP_TAC[
LE_1] THEN
REWRITE_TAC[ARITH_RULE `(2 * p + 2) - 1 = 2 * p + 1`;
ADD_SUB] THEN
REWRITE_TAC[
REAL_POW_NEG;
EVEN_ADD;
EVEN_MULT;
ARITH_EVEN] THEN
REWRITE_TAC[
REAL_POW_ONE; REAL_MUL_LID] THEN
DISCH_TAC THEN ASM_REWRITE_TAC[
COMPLEX_ADD_RID] THEN
CONV_TAC COMPLEX_RING) in
let lemma2 = prove
(`!z n. 2 <= n /\ (Im z = &0 ==> &0 < Re z)
==> (\t. inv(z + Cx(drop t)) pow n)
absolutely_integrable_on {t | &0 <= drop t}`,
REPEAT STRIP_TAC THEN REWRITE_TAC[COMPLEX_POW_INV] THEN
MATCH_MP_TAC
MEASURABLE_BOUNDED_BY_INTEGRABLE_IMP_ABSOLUTELY_INTEGRABLE THEN
EXISTS_TAC `\t. lift(inv(max (abs(Im z)) (Re z + drop t)) pow n)` THEN
REWRITE_TAC[] THEN REPEAT CONJ_TAC THENL
[MATCH_MP_TAC CONTINUOUS_IMP_MEASURABLE_ON_LEBESGUE_MEASURABLE_SUBSET THEN
CONJ_TAC THENL
[MATCH_MP_TAC CONTINUOUS_ON_COMPLEX_INV THEN
SIMP_TAC[CONTINUOUS_ON_COMPLEX_POW; CONTINUOUS_ON_ADD;
CONTINUOUS_ON_CONST; CONTINUOUS_ON_ID; CONTINUOUS_ON_CX_DROP] THEN
REWRITE_TAC[IN_ELIM_THM; GSYM FORALL_DROP] THEN
REWRITE_TAC[COMPLEX_POW_EQ_0] THEN
REWRITE_TAC[COMPLEX_EQ; RE_ADD; IM_ADD; RE_CX; IM_CX] THEN
ASM_REAL_ARITH_TAC;
MATCH_MP_TAC LEBESGUE_MEASURABLE_CONVEX THEN
MATCH_MP_TAC IS_INTERVAL_CONVEX THEN
REWRITE_TAC[IS_INTERVAL_1; IN_ELIM_THM] THEN REAL_ARITH_TAC];
ALL_TAC;
REWRITE_TAC[FORALL_IN_GSPEC; GSYM FORALL_DROP; LIFT_DROP] THEN
X_GEN_TAC `x:real` THEN DISCH_TAC THEN
REWRITE_TAC[GSYM COMPLEX_POW_INV; COMPLEX_NORM_POW] THEN
MATCH_MP_TAC REAL_POW_LE2 THEN
REWRITE_TAC[NORM_POS_LE; COMPLEX_NORM_INV] THEN
MATCH_MP_TAC REAL_LE_INV2 THEN
MP_TAC(SPEC `z + Cx x` COMPLEX_NORM_GE_RE_IM) THEN
REWRITE_TAC[RE_ADD; IM_ADD; RE_CX; IM_CX] THEN ASM_REAL_ARITH_TAC] THEN
REWRITE_TAC[REAL_ARITH `max m n = if n <= m then m else n`] THEN
REWRITE_TAC[COND_RAND; COND_RATOR] THEN
MATCH_MP_TAC INTEGRABLE_CASES THEN REWRITE_TAC[IN_ELIM_THM] THEN
CONJ_TAC THENL
[SUBGOAL_THEN
`{t | &0 <= drop t /\ Re z + drop t <= abs (Im z)} =
interval[vec 0,lift(abs(Im z) - Re z)]`
(fun th -> REWRITE_TAC[th; INTEGRABLE_CONST]) THEN
SIMP_TAC[EXTENSION; IN_ELIM_THM; IN_INTERVAL_1; DROP_VEC; LIFT_DROP] THEN
REAL_ARITH_TAC;
ALL_TAC] THEN
MATCH_MP_TAC(REWRITE_RULE[IMP_IMP] INTEGRABLE_SPIKE_SET) THEN
EXISTS_TAC
`{t | abs(Im z) - Re z <= drop t} INTER {t | &0 <= drop t}` THEN
CONJ_TAC THENL
[MATCH_MP_TAC NEGLIGIBLE_SUBSET THEN
EXISTS_TAC `{lift(abs(Im z) - Re z)}` THEN
REWRITE_TAC[NEGLIGIBLE_SING] THEN
REWRITE_TAC[SUBSET; IN_DIFF; IN_UNION;
IN_ELIM_THM; IN_SING; IN_INTER] THEN
REWRITE_TAC[GSYM DROP_EQ; LIFT_DROP] THEN ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[GSYM INTEGRABLE_RESTRICT_INTER] THEN
REWRITE_TAC[integrable_on] THEN
ONCE_REWRITE_TAC[HAS_INTEGRAL_LIM_AT_POSINFINITY] THEN
REWRITE_TAC[INTEGRABLE_RESTRICT_INTER; INTEGRAL_RESTRICT_INTER] THEN
SUBGOAL_THEN
`!a. {t | abs(Im z) - Re z <= drop t} INTER interval[vec 0,a] =
interval[lift(max (&0) (abs (Im z) - Re z)),a]`
(fun th -> REWRITE_TAC[th])
THENL
[REWRITE_TAC[EXTENSION; IN_INTER; IN_ELIM_THM; IN_INTERVAL_1] THEN
REWRITE_TAC[LIFT_DROP; DROP_VEC] THEN REAL_ARITH_TAC;
ALL_TAC] THEN
SUBGOAL_THEN
`!b. &0 <= drop b /\ abs (Im z) - Re z <= drop b
==> ((\x. lift(inv (Re z + drop x) pow n)) has_integral
lift
(inv((&1 - &n) * (Re z + drop b) pow (n - 1)) -
inv((&1 - &n) *
(Re z + max(&0) (abs(Im z) - Re z)) pow (n - 1))))
(interval[lift(max (&0) (abs (Im z) - Re z)),b])`
ASSUME_TAC THENL
[REPEAT STRIP_TAC THEN
MP_TAC(ISPECL
[`\x. lift(inv((&1 - &n) * (Re z + drop x) pow (n - 1)))`;
`\x. lift(inv(Re z + drop x) pow n)`;
`lift(max (&0) (abs (Im z) - Re z))`;
`b:real^1`] FUNDAMENTAL_THEOREM_OF_CALCULUS) THEN
ASM_REWRITE_TAC[IN_INTERVAL_1; LIFT_DROP; REAL_MAX_LE] THEN
ANTS_TAC THENL
[REWRITE_TAC[FORALL_LIFT; LIFT_DROP; INTERVAL_REAL_INTERVAL] THEN
REWRITE_TAC[REWRITE_RULE[o_DEF]
(GSYM HAS_REAL_VECTOR_DERIVATIVE_WITHIN)] THEN
REWRITE_TAC[REAL_INV_MUL; REAL_INV_POW] THEN REPEAT STRIP_TAC THEN
REAL_DIFF_TAC THEN
ASM_SIMP_TAC[GSYM REAL_OF_NUM_SUB; ARITH_RULE `2 <= n ==> 1 <= n`] THEN
CONJ_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
ASM_SIMP_TAC[real_div; REAL_OF_NUM_LE; REAL_FIELD
`&2 <= x ==> inv(&1 - x) * (x - &1) * y * --(&0 + &1) * p =
y * p`] THEN
REWRITE_TAC[GSYM REAL_INV_MUL; GSYM REAL_INV_POW] THEN
AP_TERM_TAC THEN REWRITE_TAC[GSYM REAL_POW_ADD] THEN
AP_TERM_TAC THEN ASM_ARITH_TAC;
REWRITE_TAC[LIFT_DROP; GSYM LIFT_SUB]];
ALL_TAC] THEN
EXISTS_TAC `lift(inv
((&n - &1) *
(Re z + max (&0) (abs(Im z) - Re z)) pow (n - 1)))` THEN
CONJ_TAC THENL
[X_GEN_TAC `b:real^1` THEN
ASM_CASES_TAC `&0 <= drop b /\ abs (Im z) - Re z <= drop b` THENL
[ASM_MESON_TAC[integrable_on]; MATCH_MP_TAC INTEGRABLE_ON_NULL] THEN
REWRITE_TAC[CONTENT_EQ_0_1; LIFT_DROP] THEN ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC
`\b. lift(inv((&1 - &n) * (Re z + b) pow (n - 1)) -
inv((&1 - &n) *
(Re z + max (&0) (abs(Im z) - Re z)) pow (n - 1)))` THEN
CONJ_TAC THENL
[REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN
EXISTS_TAC `max(&0) (abs(Im z) - Re z)` THEN
REWRITE_TAC[FORALL_DROP; LIFT_DROP] THEN
REPEAT STRIP_TAC THEN CONV_TAC SYM_CONV THEN
MATCH_MP_TAC INTEGRAL_UNIQUE THEN FIRST_X_ASSUM MATCH_MP_TAC THEN
ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[LIFT_SUB; VECTOR_SUB] THEN
GEN_REWRITE_TAC LAND_CONV [GSYM VECTOR_ADD_LID] THEN
MATCH_MP_TAC LIM_ADD THEN
REWRITE_TAC[GSYM LIFT_NEG; GSYM REAL_INV_NEG; GSYM REAL_MUL_LNEG] THEN
REWRITE_TAC[REAL_NEG_SUB; LIM_CONST] THEN
REWRITE_TAC[REAL_INV_MUL; LIFT_CMUL] THEN
MATCH_MP_TAC LIM_NULL_CMUL THEN
REWRITE_TAC[GSYM LIFT_NUM; GSYM LIM_CX_LIFT] THEN
REWRITE_TAC[CX_INV; CX_POW; GSYM COMPLEX_POW_INV; CX_ADD] THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN
MATCH_MP_TAC LIM_INV_X_POW_OFFSET THEN ASM_ARITH_TAC) in
let lemma3 = prove
(`!z n. 2 <= n /\ (Im z = &0 ==> &0 < Re z)
==> (\t. Cx(bernoulli n (frac (drop t))) / (z + Cx(drop t)) pow n)
integrable_on {t | &0 <= drop t}`,
REPEAT STRIP_TAC THEN
REWRITE_TAC[complex_div; GSYM COMPLEX_CMUL] THEN
MATCH_MP_TAC ABSOLUTELY_INTEGRABLE_IMP_INTEGRABLE THEN
MATCH_MP_TAC ABSOLUTELY_INTEGRABLE_ON_MUL_BERNOULLI_FRAC THEN
ASM_SIMP_TAC[GSYM COMPLEX_POW_INV; lemma2]) in
let lemma4 = prove
(`!z p. (Im z = &0 ==> &0 < Re z) /\ 1 <= p
==> ((\t. Cx(bernoulli 1 (frac(drop t))) / (z + Cx(drop t)))
has_integral
(integral {t | &0 <= drop t}
(\t. Cx(bernoulli (2 * p + 1) (frac(drop t))) /
(z + Cx(drop t)) pow (2 * p + 1)) /
Cx(&(2 * p + 1)) -
vsum(1..p) (\k. Cx(bernoulli (2 * k) (&0) /
(&4 * &k pow 2 - &2 * &k)) /
z pow (2 * k - 1))))
{t | &0 <= drop t}`,
REPEAT STRIP_TAC THEN
MATCH_MP_TAC HAS_INTEGRAL_LIM_SEQUENTIALLY THEN CONJ_TAC THENL
[REWRITE_TAC[o_DEF; LIFT_DROP; complex_div; COMPLEX_VEC_0] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_LMUL_BOUNDED THEN
EXISTS_TAC `&1 / &2` THEN CONJ_TAC THENL
[MATCH_MP_TAC ALWAYS_EVENTUALLY THEN X_GEN_TAC `x:real` THEN
REWRITE_TAC[BERNOULLI_CONV `bernoulli 1 x`] THEN DISJ1_TAC THEN
REWRITE_TAC[COMPLEX_NORM_CX] THEN
MP_TAC(SPEC `x:real` FLOOR_FRAC) THEN REAL_ARITH_TAC;
MATCH_MP_TAC(REWRITE_RULE[o_DEF] LIM_INFINITY_POSINFINITY_CX) THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN REWRITE_TAC[LIM_INV_Z_OFFSET]];
ALL_TAC] THEN
MP_TAC(GEN `n:num` (CONJ
(SPECL [`0`; `n:num`; `z:complex`] lemma1)
(SPECL [`p:num`; `n:num`; `z:complex`] lemma1))) THEN
ASM_REWRITE_TAC[FORALL_AND_THM] THEN ONCE_REWRITE_TAC[IMP_CONJ] THEN
REWRITE_TAC[VSUM_CLAUSES_NUMSEG] THEN CONV_TAC NUM_REDUCE_CONV THEN
REWRITE_TAC[SIMPLE_COMPLEX_ARITH `a * Cx(&1) / b = a / b`] THEN
REWRITE_TAC[LIFT_NUM; COMPLEX_POW_1; COMPLEX_DIV_1] THEN
DISCH_THEN(fun th -> ONCE_REWRITE_TAC[th]) THEN
REWRITE_TAC[VECTOR_ARITH `a + vec 0 + x:real^N = a + y <=> x = y`] THEN
DISCH_THEN(fun th -> ONCE_REWRITE_TAC[th]) THEN
GEN_REWRITE_TAC LAND_CONV [COMPLEX_RING `a - b:complex = --b + a`] THEN
MATCH_MP_TAC LIM_ADD THEN CONJ_TAC THENL
[REWRITE_TAC[GSYM VSUM_NEG] THEN MATCH_MP_TAC LIM_VSUM THEN
REWRITE_TAC[FINITE_NUMSEG; IN_NUMSEG] THEN
X_GEN_TAC `k:num` THEN STRIP_TAC THEN
REWRITE_TAC[CX_DIV; complex_div; GSYM COMPLEX_MUL_RNEG] THEN
REWRITE_TAC[GSYM COMPLEX_MUL_ASSOC] THEN
MATCH_MP_TAC LIM_COMPLEX_LMUL THEN
REWRITE_TAC[SIMPLE_COMPLEX_ARITH
`inv x * (y * w - y * z):complex = y / x * (w - z)`] THEN
SUBGOAL_THEN
`Cx(&(FACT(2 * k - 2))) / Cx(&(FACT(2 * k))) =
inv(Cx(&4 * &k pow 2 - &2 * &k))`
(fun th -> ONCE_REWRITE_TAC[th])
THENL
[MATCH_MP_TAC(COMPLEX_FIELD
`~(y = Cx(&0)) /\ ~(z = Cx(&0)) /\ x * z = y
==> x / y = inv(z)`) THEN
REWRITE_TAC[CX_INJ; REAL_OF_NUM_EQ; FACT_NZ] THEN
REWRITE_TAC[REAL_ENTIRE; REAL_ARITH
`&4 * x pow 2 - &2 * x = (&2 * x) * (&2 * x - &1)`] THEN
ASM_SIMP_TAC[ARITH_RULE `1 <= k ==> 1 <= 2 * k`; REAL_OF_NUM_SUB;
REAL_OF_NUM_MUL; GSYM CX_MUL; CX_INJ; REAL_OF_NUM_EQ] THEN
CONJ_TAC THENL [ASM_ARITH_TAC; ALL_TAC] THEN
UNDISCH_TAC `1 <= k` THEN SPEC_TAC(`k:num`,`k:num`) THEN
INDUCT_TAC THEN REWRITE_TAC[ADD_CLAUSES; MULT_CLAUSES] THEN
CONV_TAC NUM_REDUCE_CONV THEN DISCH_THEN (K ALL_TAC) THEN
REWRITE_TAC[ARITH_RULE `(2 + k) - 2 = k /\ (2 + k) - 1 = k + 1`] THEN
REWRITE_TAC[ARITH_RULE `2 + k = SUC(SUC k)`; FACT] THEN ARITH_TAC;
MATCH_MP_TAC LIM_COMPLEX_LMUL THEN
ONCE_REWRITE_TAC[COMPLEX_RING `--z = Cx(&0) - z`] THEN
MATCH_MP_TAC LIM_SUB THEN REWRITE_TAC[LIM_CONST] THEN
REWRITE_TAC[GSYM COMPLEX_POW_INV] THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN
MATCH_MP_TAC LIM_INV_N_POW_OFFSET THEN ASM_ARITH_TAC];
MP_TAC(SPECL [`z:complex`; `2 * p + 1`] lemma3) THEN
ASM_REWRITE_TAC[] THEN ANTS_TAC THENL [ASM_ARITH_TAC; DISCH_TAC] THEN
FIRST_ASSUM(MP_TAC o SPEC `Cx(&(FACT(2 * p)))` o
MATCH_MP INTEGRABLE_COMPLEX_LMUL) THEN
DISCH_THEN(MP_TAC o MATCH_MP INTEGRABLE_INTEGRAL) THEN
ASM_SIMP_TAC[INTEGRAL_COMPLEX_LMUL] THEN
REWRITE_TAC[HAS_INTEGRAL_LIM_AT_POSINFINITY] THEN
DISCH_THEN(MP_TAC o
MATCH_MP LIM_POSINFINITY_SEQUENTIALLY o CONJUNCT2) THEN
DISCH_THEN(MP_TAC o SPEC `inv(Cx(&(FACT(2 * p + 1))))` o
MATCH_MP LIM_COMPLEX_RMUL) THEN
MATCH_MP_TAC EQ_IMP THEN AP_THM_TAC THEN BINOP_TAC THENL
[REWRITE_TAC[LIFT_NUM; FUN_EQ_THM] THEN
X_GEN_TAC `n:num` THEN REWRITE_TAC[complex_div] THEN
AP_THM_TAC THEN AP_TERM_TAC THEN MATCH_MP_TAC INTEGRAL_EQ THEN
REWRITE_TAC[complex_div; COMPLEX_MUL_AC];
REWRITE_TAC[complex_div] THEN MATCH_MP_TAC(COMPLEX_RING
`x * y:complex = z ==> (x * i) * y = i * z`) THEN
REWRITE_TAC[GSYM ADD1; FACT; GSYM REAL_OF_NUM_MUL; CX_MUL] THEN
MATCH_MP_TAC(COMPLEX_FIELD
`~(a = Cx(&0)) /\ ~(b = Cx(&0)) ==> a * inv(b * a) = inv b`) THEN
REWRITE_TAC[CX_INJ; REAL_OF_NUM_EQ; FACT_NZ; NOT_SUC]]]) in
let lemma5 = prove
(`!z n. 1 <= n /\ (Im z = &0 ==> &0 < Re z)
==> (\t. Cx(bernoulli n (frac (drop t))) / (z + Cx(drop t)) pow n)
integrable_on {t | &0 <= drop t}`,
REPEAT STRIP_TAC THEN FIRST_X_ASSUM(DISJ_CASES_TAC o MATCH_MP (ARITH_RULE
`1 <= n ==> n = 1 \/ 2 <= n`)) THEN
ASM_SIMP_TAC[lemma3; COMPLEX_POW_1] THEN
REWRITE_TAC[integrable_on] THEN
ASM_MESON_TAC[LE_REFL; lemma4]) in
let lemma6 = prove
(`!z. (Im z = &0 ==> &0 < Re z)
==> lgamma(z) =
((z - Cx(&1) / Cx(&2)) * clog(z) - z + Cx(&1)) +
(integral {t | &0 <= drop t}
(\t. Cx(bernoulli 1 (frac(drop t))) /
(Cx(&1) + Cx(drop t))) -
integral {t | &0 <= drop t}
(\t. Cx(bernoulli 1 (frac(drop t))) /
(z + Cx(drop t))))`,
REPEAT STRIP_TAC THEN
SUBGOAL_THEN `!n. ~(z + Cx(&n) = Cx(&0))` ASSUME_TAC THENL
[REWRITE_TAC[COMPLEX_EQ; IM_ADD; RE_ADD; IM_CX; RE_CX] THEN
ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
SUBGOAL_THEN `~(z = Cx(&0))` ASSUME_TAC THENL
[ASM_MESON_TAC[COMPLEX_ADD_RID]; ALL_TAC] THEN
FIRST_ASSUM(MP_TAC o MATCH_MP LGAMMA_ALT) THEN
MATCH_MP_TAC(REWRITE_RULE[IMP_CONJ_ALT]
(REWRITE_RULE[TRIVIAL_LIMIT_SEQUENTIALLY]
(ISPEC `sequentially` LIM_UNIQUE))) THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN EXISTS_TAC
`\n. z * clog (Cx(&n)) - clog (Cx(&n + &1)) +
((Cx(&1) + Cx(&n) + Cx(&1) / Cx(&2)) * clog (Cx(&1) + Cx(&n)) -
(Cx(&1) - Cx(&1) / Cx(&2)) * clog (Cx(&1)) -
Cx(&n) +
integral (interval [lift (&0),lift (&n)])
(\x. Cx(bernoulli 1 (frac (drop x))) *
Cx(&1) / (Cx(&1) + Cx(drop x)) pow 1) /
Cx(&1)) -
((z + Cx(&n) + Cx(&1) / Cx(&2)) * clog (z + Cx(&n)) -
(z - Cx(&1) / Cx(&2)) * clog z -
Cx(&n) +
integral (interval [lift (&0),lift (&n)])
(\x. Cx(bernoulli 1 (frac (drop x))) *
Cx(&1) / (z + Cx(drop x)) pow 1) /
Cx(&1))` THEN
REWRITE_TAC[] THEN CONJ_TAC THENL
[REWRITE_TAC[EVENTUALLY_SEQUENTIALLY] THEN
EXISTS_TAC `1` THEN X_GEN_TAC `n:num` THEN DISCH_TAC THEN
SUBGOAL_THEN
`clog(Cx(&(FACT n))) =
vsum(0..n) (\m. clog(Cx(&m) + Cx(&1))) - clog(Cx(&n + &1))`
SUBST1_TAC THENL
[REWRITE_TAC[REAL_OF_NUM_ADD; GSYM CX_ADD] THEN
REWRITE_TAC[GSYM(SPEC `1` VSUM_OFFSET); ADD_CLAUSES] THEN
SIMP_TAC[GSYM NPRODUCT_FACT; REAL_OF_NUM_NPRODUCT; FINITE_NUMSEG;
GSYM CX_LOG; LOG_PRODUCT; PRODUCT_POS_LT; IN_NUMSEG;
REAL_OF_NUM_LT; LE_1; GSYM VSUM_CX] THEN
REWRITE_TAC[GSYM ADD1; VSUM_CLAUSES_NUMSEG] THEN
REWRITE_TAC[ARITH_RULE `1 <= SUC n`; REAL_OF_NUM_SUC] THEN
SIMP_TAC[CX_LOG; REAL_OF_NUM_LT; LT_0] THEN
REWRITE_TAC[COMPLEX_RING `(a + b) - b:complex = a`];
ONCE_REWRITE_TAC[COMPLEX_RING
`(a + b - c) - d:complex = (a - c) + (b - d)`]] THEN
MP_TAC(SPECL [`0`; `n:num`] lemma1) THEN
CONV_TAC NUM_REDUCE_CONV THEN
ASM_REWRITE_TAC[VSUM_CLAUSES_NUMSEG; ARITH; VECTOR_ADD_LID] THEN
ASM_SIMP_TAC[RE_CX; REAL_LT_01];
ALL_TAC] THEN
REWRITE_TAC[COMPLEX_POW_1; COMPLEX_DIV_1] THEN
ONCE_REWRITE_TAC[COMPLEX_RING
`a + (b - c - n + x) - (d - e - n + y):complex =
(a + (b - c) - (d - e)) + (x - y)`] THEN
MATCH_MP_TAC LIM_ADD THEN CONJ_TAC THENL
[REWRITE_TAC[CLOG_1; COMPLEX_MUL_RZERO; COMPLEX_SUB_RZERO] THEN
ONCE_REWRITE_TAC[LIM_NULL_COMPLEX] THEN REWRITE_TAC[] THEN
ONCE_REWRITE_TAC[REAL_ADD_SYM] THEN REWRITE_TAC[CX_ADD] THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN EXISTS_TAC
`\n. z * (clog(Cx(&n)) - clog(z + Cx(&n))) +
(Cx(&n) + Cx(&1) / Cx(&2)) *
(clog(Cx(&1) + Cx(&n)) - clog(z + Cx(&n))) -
(Cx(&1) - z)` THEN
CONJ_TAC THENL
[MATCH_MP_TAC ALWAYS_EVENTUALLY THEN CONV_TAC COMPLEX_RING;
ALL_TAC] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_ADD THEN CONJ_TAC THENL
[MATCH_MP_TAC LIM_NULL_COMPLEX_LMUL THEN
GEN_REWRITE_TAC (RATOR_CONV o LAND_CONV o BINDER_CONV o
LAND_CONV o RAND_CONV) [GSYM COMPLEX_ADD_LID] THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN REWRITE_TAC[LIM_SUB_CLOG];
ONCE_REWRITE_TAC[COMPLEX_RING
`(a + h) * x - y:complex = h * x + a * x - y`] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_ADD THEN REWRITE_TAC[] THEN
ONCE_REWRITE_TAC[COMPLEX_ADD_SYM] THEN
SIMP_TAC[LIM_NULL_COMPLEX_LMUL; LIM_SUB_CLOG] THEN
REWRITE_TAC[GSYM LIM_NULL_COMPLEX; LIM_N_MUL_SUB_CLOG]];
REWRITE_TAC[complex_div; COMPLEX_MUL_LID] THEN
REWRITE_TAC[LIFT_NUM] THEN MATCH_MP_TAC LIM_SUB THEN
CONJ_TAC THEN REWRITE_TAC[GSYM LIFT_NUM] THEN
MATCH_MP_TAC LIM_POSINFINITY_SEQUENTIALLY THEN
REWRITE_TAC[LIFT_NUM] THEN
MATCH_MP_TAC(MATCH_MP (TAUT `(p <=> q /\ r) ==> (p ==> r)`)
(SPEC_ALL HAS_INTEGRAL_LIM_AT_POSINFINITY)) THEN
REWRITE_TAC[GSYM HAS_INTEGRAL_INTEGRAL; GSYM complex_div] THEN
ONCE_REWRITE_TAC[COMPLEX_RING `z + x:complex = (z + x) pow 1`] THEN
MATCH_MP_TAC lemma5 THEN
ASM_REWRITE_TAC[LE_REFL; RE_CX; REAL_LT_01]]) in
let lemma7 = prove
(`((\y. integral {t | &0 <= drop t}
(\t. Cx (bernoulli 1 (frac (drop t))) / (ii * Cx y + Cx(drop t))))
--> Cx(&0)) at_posinfinity`,
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC
`\y. integral {t | &0 <= drop t}
(\t. Cx(bernoulli 3 (frac (drop t))) /
(ii * Cx y + Cx(drop t)) pow 3) / Cx(&3) -
Cx(bernoulli 2 (&0) / &2) / (ii * Cx y)` THEN
CONJ_TAC THENL
[REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN
EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
SUBGOAL_THEN `&0 < y` ASSUME_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
MP_TAC(ISPECL [`ii * Cx y`; `1`] lemma4) THEN
ASM_SIMP_TAC[IM_MUL_II; RE_CX; REAL_LT_IMP_NZ; LE_REFL] THEN
REWRITE_TAC[VSUM_SING_NUMSEG] THEN
CONV_TAC NUM_REDUCE_CONV THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
REWRITE_TAC[COMPLEX_POW_1] THEN
DISCH_THEN(SUBST1_TAC o MATCH_MP INTEGRAL_UNIQUE) THEN REFL_TAC;
ALL_TAC] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_SUB THEN CONJ_TAC THENL
[ALL_TAC;
REWRITE_TAC[complex_div] THEN MATCH_MP_TAC LIM_NULL_COMPLEX_LMUL THEN
REWRITE_TAC[COMPLEX_INV_MUL] THEN MATCH_MP_TAC LIM_NULL_COMPLEX_LMUL THEN
REWRITE_TAC[LIM_INV_X]] THEN
ONCE_REWRITE_TAC[complex_div] THEN MATCH_MP_TAC LIM_NULL_COMPLEX_RMUL THEN
MATCH_MP_TAC LIM_NULL_COMPARISON_COMPLEX THEN
EXISTS_TAC
`\y. Cx(&1 / &2 / y) *
integral {t | &0 <= drop t}
(\t. inv(Cx(&1) + Cx(drop t)) pow 2)` THEN
CONJ_TAC THENL
[ALL_TAC;
MATCH_MP_TAC LIM_NULL_COMPLEX_RMUL THEN
REWRITE_TAC[real_div; CX_MUL] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_LMUL THEN
REWRITE_TAC[LIM_INV_X; CX_INV]] THEN
REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN
EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
SIMP_TAC[GSYM INTEGRAL_COMPLEX_LMUL; lemma2;
ABSOLUTELY_INTEGRABLE_IMP_INTEGRABLE; ARITH;
RE_CX; REAL_OF_NUM_LT] THEN
MATCH_MP_TAC(REAL_ARITH
`abs(Re z) <= norm z /\ x <= Re z ==> x <= norm z`) THEN
REWRITE_TAC[COMPLEX_NORM_GE_RE_IM] THEN REWRITE_TAC[RE_DEF] THEN
MATCH_MP_TAC INTEGRAL_NORM_BOUND_INTEGRAL_COMPONENT THEN
REWRITE_TAC[DIMINDEX_2; ARITH] THEN REPEAT CONJ_TAC THENL
[MATCH_MP_TAC lemma5 THEN
REWRITE_TAC[ARITH; IM_MUL_II; RE_CX] THEN ASM_REAL_ARITH_TAC;
MATCH_MP_TAC INTEGRABLE_COMPLEX_LMUL THEN
MATCH_MP_TAC ABSOLUTELY_INTEGRABLE_IMP_INTEGRABLE THEN
MATCH_MP_TAC lemma2 THEN
REWRITE_TAC[RE_CX; REAL_LT_01] THEN ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[FORALL_LIFT; IN_ELIM_THM; LIFT_DROP; GSYM RE_DEF] THEN
X_GEN_TAC `x:real` THEN DISCH_TAC THEN
REWRITE_TAC[GSYM CX_ADD; GSYM CX_INV; GSYM CX_MUL; GSYM CX_POW] THEN
REWRITE_TAC[COMPLEX_NORM_DIV; COMPLEX_NORM_CX; RE_CX] THEN
REWRITE_TAC[REAL_ARITH `&1 / &2 / y * x = (&1 / &4) * (&2 / y * x)`] THEN
REWRITE_TAC[real_div] THEN MATCH_MP_TAC REAL_LE_MUL2 THEN
REWRITE_TAC[REAL_ABS_POS; REAL_LE_INV_EQ; NORM_POS_LE] THEN CONJ_TAC THENL
[MP_TAC(SPECL [`3`; `frac x`] BERNOULLI_BOUND) THEN
CONV_TAC NUM_REDUCE_CONV THEN
REWRITE_TAC[BERNOULLI_CONV `bernoulli 2 (&0)`] THEN
ANTS_TAC THENL [ALL_TAC; REAL_ARITH_TAC] THEN
SIMP_TAC[IN_REAL_INTERVAL; REAL_LT_IMP_LE; FLOOR_FRAC];
ALL_TAC] THEN
REWRITE_TAC[REAL_POW_INV; GSYM REAL_INV_MUL; GSYM REAL_MUL_ASSOC] THEN
GEN_REWRITE_TAC (RAND_CONV o LAND_CONV) [GSYM REAL_INV_INV] THEN
REWRITE_TAC[GSYM REAL_INV_MUL] THEN MATCH_MP_TAC REAL_LE_INV2 THEN
CONJ_TAC THENL
[MATCH_MP_TAC REAL_LT_MUL THEN CONV_TAC REAL_RAT_REDUCE_CONV THEN
MATCH_MP_TAC REAL_LT_MUL THEN
CONJ_TAC THENL [ASM_REAL_ARITH_TAC; MATCH_MP_TAC REAL_POW_LT] THEN
ASM_REAL_ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[COMPLEX_NORM_POW; REAL_ARITH
`x pow 3 = (x:real) * x pow 2`] THEN
ONCE_REWRITE_TAC[REAL_ARITH `inv(&2) * y * x = y * x / &2`] THEN
MATCH_MP_TAC REAL_LE_MUL2 THEN REPEAT CONJ_TAC THENL
[ASM_REAL_ARITH_TAC;
W(MP_TAC o PART_MATCH (rand o rand) COMPLEX_NORM_GE_RE_IM o
rand o snd) THEN
REWRITE_TAC[IM_ADD; IM_MUL_II; RE_CX; IM_CX] THEN REAL_ARITH_TAC;
REWRITE_TAC[REAL_ARITH `&0 <= x / &2 <=> &0 <= x`] THEN
MATCH_MP_TAC REAL_POW_LE THEN ASM_REAL_ARITH_TAC;
REWRITE_TAC[COMPLEX_SQNORM] THEN
REWRITE_TAC[RE_ADD; IM_ADD; RE_MUL_II; IM_MUL_II; IM_CX; RE_CX] THEN
REWRITE_TAC[REAL_ADD_LID; REAL_ADD_RID; REAL_NEG_0] THEN
MATCH_MP_TAC(REAL_ARITH
`&1 pow 2 <= y pow 2 /\ &0 <= (x - &1) pow 2
==> (&1 + x) pow 2 / &2 <= x pow 2 + y pow 2`) THEN
REWRITE_TAC[REAL_LE_POW_2] THEN MATCH_MP_TAC REAL_POW_LE2 THEN
ASM_REAL_ARITH_TAC]) in
let lemma8 = prove
(`integral {t | &0 <= drop t}
(\t. Cx (bernoulli 1 (frac (drop t))) / (Cx(&1) + Cx(drop t))) =
Cx(log(&2 * pi) / &2 - &1)`,
MATCH_MP_TAC(MESON[REAL]
`real z /\ Re z = y ==> z = Cx y`) THEN
CONJ_TAC THENL
[MATCH_MP_TAC REAL_COMPLEX_INTEGRAL THEN CONJ_TAC THENL
[ONCE_REWRITE_TAC[COMPLEX_RING `z + x:complex = (z + x) pow 1`] THEN
MATCH_MP_TAC lemma5 THEN
ASM_REWRITE_TAC[LE_REFL; RE_CX; REAL_LT_01];
REWRITE_TAC[GSYM CX_ADD; GSYM CX_DIV; REAL_CX]];
GEN_REWRITE_TAC I [GSYM CX_INJ]] THEN
MATCH_MP_TAC(ISPEC `at_posinfinity` LIM_UNIQUE) THEN
EXISTS_TAC
`\y:real. Cx(Re(integral {t | &0 <= drop t}
(\t. Cx(bernoulli 1 (frac (drop t))) / (Cx(&1) + Cx(drop t)))))` THEN
REWRITE_TAC[TRIVIAL_LIMIT_AT_POSINFINITY; LIM_CONST] THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC
`\y. Cx(Re(lgamma(ii * Cx y) -
((ii * Cx y - Cx(&1) / Cx(&2)) * clog(ii * Cx y) -
ii * Cx y + Cx(&1)) +
integral {t | &0 <= drop t}
(\t. Cx(bernoulli 1 (frac(drop t))) /
(ii * Cx y + Cx(drop t)))))` THEN
REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN CONJ_TAC THENL
[EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
AP_TERM_TAC THEN AP_TERM_TAC THEN
MP_TAC(SPEC `ii * Cx y` lemma6) THEN
ASM_SIMP_TAC[IM_MUL_II; RE_CX; REAL_LT_IMP_NZ] THEN
ANTS_TAC THENL [ASM_REAL_ARITH_TAC; DISCH_THEN SUBST1_TAC] THEN
REWRITE_TAC[COMPLEX_RING `(s + i - j) - s + j:complex = i`];
ALL_TAC] THEN
REWRITE_TAC[RE_ADD; RE_SUB; CX_ADD; CX_SUB] THEN
GEN_REWRITE_TAC LAND_CONV [GSYM COMPLEX_ADD_RID] THEN
MATCH_MP_TAC LIM_ADD THEN CONJ_TAC THENL
[ALL_TAC;
MATCH_MP_TAC LIM_NULL_COMPARISON_COMPLEX THEN
EXISTS_TAC `\y. integral {t | &0 <= drop t}
(\t. Cx(bernoulli 1 (frac(drop t))) / (ii * Cx y + Cx(drop t)))` THEN
REWRITE_TAC[lemma7] THEN MATCH_MP_TAC ALWAYS_EVENTUALLY THEN
REWRITE_TAC[COMPLEX_NORM_CX; COMPLEX_NORM_GE_RE_IM]] THEN
REWRITE_TAC[RE_MUL_II; IM_CX; REAL_NEG_0; COMPLEX_SUB_RZERO; RE_CX] THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC
`\y. Cx(log(norm(cgamma(ii * Cx y)))) -
(Cx(--(pi * y + log y) / &2) + Cx(&1))` THEN
REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN CONJ_TAC THENL
[EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
SUBGOAL_THEN `&0 < y` ASSUME_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
BINOP_TAC THENL
[MP_TAC(SPEC `ii * Cx y` CGAMMA_LGAMMA) THEN COND_CASES_TAC THENL
[FIRST_X_ASSUM(CHOOSE_THEN (MP_TAC o AP_TERM `Im`)) THEN
REWRITE_TAC[IM_ADD; IM_MUL_II; IM_CX; RE_CX] THEN ASM_REAL_ARITH_TAC;
DISCH_THEN SUBST1_TAC THEN REWRITE_TAC[NORM_CEXP; LOG_EXP]];
AP_THM_TAC THEN AP_TERM_TAC THEN
AP_TERM_TAC THEN CONV_TAC SYM_CONV THEN
REWRITE_TAC[COMPLEX_SUB_RDISTRIB; RE_SUB; GSYM CX_DIV] THEN
REWRITE_TAC[GSYM COMPLEX_MUL_ASSOC; RE_MUL_II; IM_MUL_II; RE_MUL_CX;
IM_MUL_CX] THEN
ASM_SIMP_TAC[RE_CX; IM_CX; REAL_LT_IMP_LE; CX_INJ; REAL_LT_IMP_NZ;
GSYM CX_LOG; CLOG_MUL_II] THEN
REWRITE_TAC[RE_ADD; IM_ADD; RE_CX; IM_CX; RE_MUL_II; IM_MUL_II] THEN
REAL_ARITH_TAC];
ALL_TAC] THEN
MATCH_MP_TAC LIM_TRANSFORM_EVENTUALLY THEN
EXISTS_TAC
`\y. Cx((log(&2 * pi) - log(&1 - inv(exp(pi * y) pow 2))) / &2 - &1)` THEN
REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN CONJ_TAC THENL
[EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
SUBGOAL_THEN `&0 < y` ASSUME_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
MP_TAC(SPEC `ii * Cx y + Cx(&1)` CGAMMA_REFLECTION) THEN
REWRITE_TAC[CGAMMA_RECURRENCE; COMPLEX_ENTIRE; II_NZ; CX_INJ] THEN
ASM_SIMP_TAC[REAL_LT_IMP_NZ] THEN
REWRITE_TAC[COMPLEX_RING `Cx(&1) - (z + Cx(&1)) = --z`] THEN
MP_TAC(SPEC `cgamma(ii * Cx y)` COMPLEX_NORM_POW_2) THEN
REWRITE_TAC[CNJ_MUL; CNJ_II; CNJ_CX; CNJ_CGAMMA] THEN
REWRITE_TAC[COMPLEX_MUL_LNEG; GSYM COMPLEX_MUL_ASSOC] THEN
DISCH_THEN(SUBST1_TAC o SYM) THEN
REWRITE_TAC[COMPLEX_RING `w * (z + Cx(&1)) = w * z + w`] THEN
REWRITE_TAC[CSIN_ADD; GSYM CX_SIN; GSYM CX_COS; SIN_PI; COS_PI] THEN
REWRITE_TAC[COMPLEX_MUL_RZERO; CX_NEG; COMPLEX_MUL_RNEG] THEN
REWRITE_TAC[COMPLEX_ADD_RID; COMPLEX_MUL_RID] THEN
REWRITE_TAC[csin; COMPLEX_RING `--ii * x * ii * y = x * y`] THEN
REWRITE_TAC[COMPLEX_RING `ii * x * ii * y = --(x * y)`] THEN
REWRITE_TAC[complex_div; COMPLEX_INV_INV; COMPLEX_INV_MUL;
COMPLEX_INV_NEG; COMPLEX_MUL_RNEG] THEN
ASM_SIMP_TAC[CX_INJ; REAL_LT_IMP_NZ; COMPLEX_FIELD
`~(y = Cx(&0))
==> (ii * y * z = --(p * i * Cx(&2) * ii) <=>
z = --(Cx(&2) * p) * inv y * i)`] THEN
REWRITE_TAC[GSYM COMPLEX_INV_MUL; GSYM CX_MUL;
GSYM CX_EXP; CEXP_NEG] THEN
REWRITE_TAC[COMPLEX_MUL_LNEG; GSYM COMPLEX_MUL_RNEG;
COMPLEX_NEG_INV] THEN
REWRITE_TAC[COMPLEX_NEG_SUB] THEN
REWRITE_TAC[GSYM CX_INV; GSYM CX_SUB; GSYM CX_MUL; GSYM CX_INV] THEN
REWRITE_TAC[CX_INJ; GSYM CX_POW] THEN
DISCH_THEN(MP_TAC o AP_TERM `sqrt`) THEN
REWRITE_TAC[POW_2_SQRT_ABS; REAL_ABS_NORM] THEN
DISCH_THEN SUBST1_TAC THEN REWRITE_TAC[REAL_INV_MUL] THEN
SUBGOAL_THEN `&0 < exp(pi * y) - inv(exp(pi * y))` ASSUME_TAC THENL
[MATCH_MP_TAC(REAL_ARITH `x < &1 /\ &1 < y ==> &0 < y - x`) THEN
CONJ_TAC THENL [MATCH_MP_TAC REAL_INV_LT_1; ALL_TAC] THEN
MATCH_MP_TAC REAL_EXP_LT_1 THEN MATCH_MP_TAC REAL_LT_MUL THEN
ASM_REWRITE_TAC[PI_POS];
ASM_SIMP_TAC[LOG_MUL; REAL_LT_INV_EQ; PI_POS; REAL_LT_MUL;
REAL_ARITH `&0 < &2 * x <=> &0 < x`; LOG_SQRT; LOG_INV]] THEN
REWRITE_TAC[GSYM CX_ADD; GSYM CX_SUB] THEN
REWRITE_TAC[REAL_ARITH
`(l2 + --l + x) / &2 - (--(py + l) / &2 + &1) =
(l2 + py + x) / &2 - &1`] THEN
SIMP_TAC[REAL_EXP_NZ; REAL_FIELD
`~(e = &0) ==> e - inv e = e * (&1 - inv(e pow 2))`] THEN
SUBGOAL_THEN `inv (exp (pi * y) pow 2) < &1` ASSUME_TAC THENL
[MATCH_MP_TAC REAL_INV_LT_1 THEN
MATCH_MP_TAC REAL_POW_LT_1 THEN
REWRITE_TAC[ARITH] THEN
MATCH_MP_TAC REAL_EXP_LT_1 THEN
MATCH_MP_TAC REAL_LT_MUL THEN
ASM_SIMP_TAC[REAL_LT_MUL; PI_POS];
ASM_SIMP_TAC[LOG_MUL; REAL_EXP_POS_LT; REAL_SUB_LT; LOG_EXP]] THEN
REWRITE_TAC[CX_INJ] THEN REAL_ARITH_TAC;
ALL_TAC] THEN
REWRITE_TAC[CX_SUB; CX_DIV] THEN MATCH_MP_TAC LIM_SUB THEN
REWRITE_TAC[LIM_CONST] THEN
MATCH_MP_TAC LIM_COMPLEX_DIV THEN REWRITE_TAC[LIM_CONST] THEN
CONJ_TAC THENL [ALL_TAC; CONV_TAC COMPLEX_RING] THEN
GEN_REWRITE_TAC LAND_CONV [GSYM COMPLEX_SUB_RZERO] THEN
MATCH_MP_TAC LIM_SUB THEN REWRITE_TAC[LIM_CONST] THEN
MP_TAC(ISPECL
[`clog`; `at_posinfinity`;
`\y. Cx(&1 - inv(exp(pi * y) pow 2))`;
`Cx(&1)`] LIM_CONTINUOUS_FUNCTION) THEN
REWRITE_TAC[CLOG_1] THEN ANTS_TAC THENL
[SIMP_TAC[CONTINUOUS_AT_CLOG; RE_CX; REAL_LT_01; CX_SUB] THEN
GEN_REWRITE_TAC LAND_CONV [GSYM COMPLEX_SUB_RZERO] THEN
MATCH_MP_TAC LIM_SUB THEN REWRITE_TAC[LIM_CONST] THEN
REWRITE_TAC[CX_INV; CX_EXP; CX_POW; GSYM COMPLEX_POW_INV] THEN
MATCH_MP_TAC LIM_NULL_COMPLEX_POW THEN CONV_TAC NUM_REDUCE_CONV THEN
REWRITE_TAC[LIM_AT_POSINFINITY] THEN
X_GEN_TAC `e:real` THEN DISCH_TAC THEN
EXISTS_TAC `&1 + inv(e)` THEN REWRITE_TAC[dist; real_ge] THEN
X_GEN_TAC `x:real` THEN DISCH_TAC THEN
REWRITE_TAC[COMPLEX_SUB_RZERO; COMPLEX_NORM_INV] THEN
MATCH_MP_TAC REAL_LT_LINV THEN ASM_REWRITE_TAC[NORM_CEXP; RE_CX] THEN
FIRST_ASSUM(MATCH_MP_TAC o MATCH_MP (REAL_ARITH
`&1 + e <= x
==> &1 + pi * x <= z /\ &0 <= x * (pi - &1) ==> e < z`)) THEN
REWRITE_TAC[REAL_EXP_LE_X] THEN MATCH_MP_TAC REAL_LE_MUL THEN
MP_TAC(SPEC `e:real` REAL_LT_INV_EQ) THEN
MP_TAC PI_APPROX_32 THEN ASM_REAL_ARITH_TAC;
MATCH_MP_TAC(REWRITE_RULE[IMP_CONJ] LIM_TRANSFORM_EVENTUALLY) THEN
REWRITE_TAC[EVENTUALLY_AT_POSINFINITY; real_ge] THEN
EXISTS_TAC `&1` THEN X_GEN_TAC `y:real` THEN DISCH_TAC THEN
SUBGOAL_THEN `&0 < y` ASSUME_TAC THENL [ASM_REAL_ARITH_TAC; ALL_TAC] THEN
MATCH_MP_TAC(GSYM CX_LOG) THEN REWRITE_TAC[REAL_SUB_LT] THEN
MATCH_MP_TAC REAL_INV_LT_1 THEN MATCH_MP_TAC REAL_POW_LT_1 THEN
REWRITE_TAC[ARITH] THEN MATCH_MP_TAC REAL_EXP_LT_1 THEN
MATCH_MP_TAC REAL_LT_MUL THEN ASM_SIMP_TAC[REAL_LT_MUL; PI_POS]]) in
CONJ_TAC THENL [MATCH_ACCEPT_TAC lemma5; ALL_TAC] THEN
REPEAT STRIP_TAC THEN ASM_SIMP_TAC[lemma6] THEN
REWRITE_TAC[lemma8; CX_SUB] THEN
REWRITE_TAC[COMPLEX_RING
`(x + Cx(&1)) + (a - Cx(&1)) - b = (x + a) - b`] THEN
REWRITE_TAC[complex_sub; GSYM COMPLEX_ADD_ASSOC] THEN AP_TERM_TAC THEN
AP_TERM_TAC THEN AP_TERM_TAC THEN ASM_CASES_TAC `p = 0` THENL
[ASM_REWRITE_TAC[VSUM_CLAUSES_NUMSEG] THEN
CONV_TAC NUM_REDUCE_CONV THEN
REWRITE_TAC[COMPLEX_POW_1; COMPLEX_DIV_1; VECTOR_ADD_LID];
MP_TAC(SPECL [`z:complex`; `p:num`] lemma4) THEN
ASM_SIMP_TAC[LE_1] THEN
DISCH_THEN(SUBST1_TAC o MATCH_MP INTEGRAL_UNIQUE) THEN
REWRITE_TAC[complex_sub; COMPLEX_NEG_ADD; VECTOR_NEG_NEG] THEN
REWRITE_TAC[COMPLEX_ADD_SYM]]);;
let LOG_GAMMA_STIRLING = prove
(`!x p. &0 < x
==> log(gamma x) =
((x - &1 / &2) * log(x) - x + log(&2 * pi) / &2) +
sum(1..p) (\k. bernoulli (2 * k) (&0) /
(&4 * &k pow 2 - &2 * &k) / x pow (2 * k - 1)) -
real_integral {t | &0 <= t}
(\t. bernoulli (2 * p + 1) (frac t) /
(x + t) pow (2 * p + 1)) /
&(2 * p + 1)`,
(* ------------------------------------------------------------------------- *)
(* Some other mathematical facts that don't directly involve the gamma *)
(* function can now be proved relatively easily: Euler's product for sin, *)
(* Wallis's product for pi, and the Gaussian integral. *)
(* ------------------------------------------------------------------------- *)
let CSIN_PRODUCT = prove
(`!z. ((\n. z * cproduct(1..n) (\m. Cx(&1) - (z / Cx(pi * &m)) pow 2))
--> csin(z)) sequentially`,
let SIN_PRODUCT = prove
(`!x. ((\n. x * product(1..n) (\m. &1 - (x / (pi * &m)) pow 2)) ---> sin(x))
sequentially`,
let WALLIS_ALT = prove
(`((\n. product(1..n) (\k. (&2 * &k) / (&2 * &k - &1) *
(&2 * &k) / (&2 * &k + &1))) ---> pi / &2)
sequentially`,