1(* ------------------------------------------------------------------------- *)
2(* Mobius Function and Inversion.                                            *)
3(* ------------------------------------------------------------------------- *)
4
5(*===========================================================================*)
6
7(* add all dependent libraries for script *)
8open HolKernel boolLib bossLib Parse;
9
10(* declare new theory at start *)
11val _ = new_theory "Mobius";
12
13(* ------------------------------------------------------------------------- *)
14
15
16
17(* val _ = load "jcLib"; *)
18open jcLib;
19
20(* Get dependent theories in lib *)
21(*
22(* val _ = load "helperNumTheory"; *)
23(* val _ = load "helperSetTheory"; *)
24(* val _ = load "helperListTheory"; *)
25*)
26(* val _ = load "GaussTheory"; *)
27
28(* open dependent theories *)
29open pred_setTheory listTheory;
30(* (* val _ = load "dividesTheory"; -- in helperNumTheory *) *)
31(* (* val _ = load "gcdTheory"; -- in helperNumTheory *) *)
32open arithmeticTheory dividesTheory gcdTheory;
33
34open helperNumTheory helperSetTheory helperListTheory;
35
36open GaussTheory;
37open EulerTheory;
38
39
40(* ------------------------------------------------------------------------- *)
41(* Mobius Function and Inversion Documentation                               *)
42(* ------------------------------------------------------------------------- *)
43(* Overloading:
44   sq_free s          = {n | n IN s /\ square_free n}
45   non_sq_free s      = {n | n IN s /\ ~(square_free n)}
46   even_sq_free s     = {n | n IN (sq_free s) /\ EVEN (CARD (prime_factors n))}
47   odd_sq_free s      = {n | n IN (sq_free s) /\ ODD (CARD (prime_factors n))}
48   less_divisors n    = {x | x IN (divisors n) /\ x <> n}
49   proper_divisors n  = {x | x IN (divisors n) /\ x <> 1 /\ x <> n}
50*)
51(* Definitions and Theorems (# are exported):
52
53   Helper Theorems:
54
55   Square-free Number and Square-free Sets:
56   square_free_def     |- !n. square_free n <=> !p. prime p /\ p divides n ==> ~(p * p divides n)
57   square_free_1       |- square_free 1
58   square_free_prime   |- !n. prime n ==> square_free n
59
60   sq_free_element     |- !s n. n IN sq_free s <=> n IN s /\ square_free n
61   sq_free_subset      |- !s. sq_free s SUBSET s
62   sq_free_finite      |- !s. FINITE s ==> FINITE (sq_free s)
63   non_sq_free_element |- !s n. n IN non_sq_free s <=> n IN s /\ ~square_free n
64   non_sq_free_subset  |- !s. non_sq_free s SUBSET s
65   non_sq_free_finite  |- !s. FINITE s ==> FINITE (non_sq_free s)
66   sq_free_split       |- !s. (s = sq_free s UNION non_sq_free s) /\
67                              (sq_free s INTER non_sq_free s = {})
68   sq_free_union       |- !s. s = sq_free s UNION non_sq_free s
69   sq_free_inter       |- !s. sq_free s INTER non_sq_free s = {}
70   sq_free_disjoint    |- !s. DISJOINT (sq_free s) (non_sq_free s)
71
72   Prime Divisors of a Number and Partitions of Square-free Set:
73   prime_factors_def      |- !n. prime_factors n = {p | prime p /\ p IN divisors n}
74   prime_factors_element  |- !n p. p IN prime_factors n <=> prime p /\ p <= n /\ p divides n
75   prime_factors_subset   |- !n. prime_factors n SUBSET divisors n
76   prime_factors_finite   |- !n. FINITE (prime_factors n)
77
78   even_sq_free_element    |- !s n. n IN even_sq_free s <=> n IN s /\ square_free n /\ EVEN (CARD (prime_factors n))
79   even_sq_free_subset     |- !s. even_sq_free s SUBSET s
80   even_sq_free_finite     |- !s. FINITE s ==> FINITE (even_sq_free s)
81   odd_sq_free_element     |- !s n. n IN odd_sq_free s <=> n IN s /\ square_free n /\ ODD (CARD (prime_factors n))
82   odd_sq_free_subset      |- !s. odd_sq_free s SUBSET s
83   odd_sq_free_finite      |- !s. FINITE s ==> FINITE (odd_sq_free s)
84   sq_free_split_even_odd  |- !s. (sq_free s = even_sq_free s UNION odd_sq_free s) /\
85                                  (even_sq_free s INTER odd_sq_free s = {})
86   sq_free_union_even_odd  |- !s. sq_free s = even_sq_free s UNION odd_sq_free s
87   sq_free_inter_even_odd  |- !s. even_sq_free s INTER odd_sq_free s = {}
88   sq_free_disjoint_even_odd  |- !s. DISJOINT (even_sq_free s) (odd_sq_free s)
89
90   Less Divisors of a number:
91   less_divisors_element       |- !n x. x IN less_divisors n <=> x < n /\ x divides n
92   less_divisors_0             |- less_divisors 0 = {}
93   less_divisors_1             |- less_divisors 1 = {}
94   less_divisors_subset_divisors   |- !n. less_divisors n SUBSET divisors n
95   less_divisors_finite        |- !n. FINITE (less_divisors n)
96   less_divisors_prime         |- !n. prime n ==> (less_divisors n = {1})
97   less_divisors_has_one       |- !n. 1 < n ==> 1 IN less_divisors n
98   less_divisors_nonzero       |- !n x. x IN less_divisors n ==> 0 < x
99   less_divisors_has_cofactor  |- !n. 0 < n ==>
100                                  !d. 1 < d /\ d IN less_divisors n ==> n DIV d IN less_divisors n
101
102   Proper Divisors of a number:
103   proper_divisors_element     |- !n x. x IN proper_divisors n <=> 1 < x /\ x < n /\ x divides n
104   proper_divisors_0           |- proper_divisors 0 = {}
105   proper_divisors_1           |- proper_divisors 1 = {}
106   proper_divisors_subset      |- !n. proper_divisors n SUBSET less_divisors n
107   proper_divisors_finite      |- !n. FINITE (proper_divisors n)
108   proper_divisors_not_one     |- !n. 1 NOTIN proper_divisors n
109   proper_divisors_by_less_divisors   |- !n. proper_divisors n = less_divisors n DELETE 1
110   proper_divisors_prime       |- !n. prime n ==> (proper_divisors n = {})
111   proper_divisors_has_cofactor   |- !n d. d IN proper_divisors n ==> n DIV d IN proper_divisors n
112   proper_divisors_min_gt_1    |- !n. proper_divisors n <> {} ==> 1 < MIN_SET (proper_divisors n)
113   proper_divisors_max_min     |- !n. proper_divisors n <> {} ==>
114                                      (MAX_SET (proper_divisors n) = n DIV MIN_SET (proper_divisors n)) /\
115                                      (MIN_SET (proper_divisors n) = n DIV MAX_SET (proper_divisors n))
116
117   Useful Properties of Less Divisors:
118   less_divisors_min             |- !n. 1 < n ==> (MIN_SET (less_divisors n) = 1)
119   less_divisors_max             |- !n. MAX_SET (less_divisors n) <= n DIV 2
120   less_divisors_subset_natural  |- !n. less_divisors n SUBSET natural (n DIV 2)
121
122   Properties of Summation equals Perfect Power:
123   perfect_power_special_inequality  |- !p. 1 < p ==> !n. p * (p ** n - 1) < (p - 1) * (2 * p ** n)
124   perfect_power_half_inequality_1   |- !p n. 1 < p /\ 0 < n ==> 2 * p ** (n DIV 2) <= p ** n
125   perfect_power_half_inequality_2   |- !p n. 1 < p /\ 0 < n ==>
126                                        (p ** (n DIV 2) - 2) * p ** (n DIV 2) <= p ** n - 2 * p ** (n DIV 2)
127   sigma_eq_perfect_power_bounds_1   |- !p. 1 < p ==>
128                          !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
129                              (!n. 0 < n ==> n * f n <= p ** n) /\
130                               !n. 0 < n ==> p ** n - 2 * p ** (n DIV 2) < n * f n
131   sigma_eq_perfect_power_bounds_2   |- !p. 1 < p ==>
132                          !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
133                              (!n. 0 < n ==> n * f n <= p ** n) /\
134                               !n. 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) < n * f n
135
136*)
137
138(* ------------------------------------------------------------------------- *)
139(* Helper Theorems                                                           *)
140(* ------------------------------------------------------------------------- *)
141
142(* ------------------------------------------------------------------------- *)
143(* Mobius Function and Inversion                                             *)
144(* ------------------------------------------------------------------------- *)
145
146
147(* ------------------------------------------------------------------------- *)
148(* Square-free Number and Square-free Sets                                   *)
149(* ------------------------------------------------------------------------- *)
150
151(* Define square-free number *)
152val square_free_def = Define`
153    square_free n = !p. prime p /\ p divides n ==> ~(p * p divides n)
154`;
155
156(* Theorem: square_free 1 *)
157(* Proof:
158       square_free 1
159   <=> !p. prime p /\ p divides 1 ==> ~(p * p divides 1)    by square_free_def
160   <=> prime 1 ==> ~(1 * 1 divides 1)                       by DIVIDES_ONE
161   <=> F ==> ~(1 * 1 divides 1)                             by NOT_PRIME_1
162   <=> T                                                    by false assumption
163*)
164val square_free_1 = store_thm(
165  "square_free_1",
166  ``square_free 1``,
167  rw[square_free_def]);
168
169(* Theorem: prime n ==> square_free n *)
170(* Proof:
171       square_free n
172   <=> !p. prime p /\ p divides n ==> ~(p * p divides n)   by square_free_def
173   By contradiction, suppose (p * p divides n).
174   Since p divides n ==> (p = n) \/ (p = 1)                by prime_def
175     and p * p divides  ==> (p * p = n) \/ (p * p = 1)     by prime_def
176     but p <> 1                                            by prime_def
177      so p * p <> 1              by MULT_EQ_1
178    Thus p * p = n = p,
179      or p = 0 \/ p = 1          by SQ_EQ_SELF
180     But p <> 0                  by NOT_PRIME_0
181     and p <> 1                  by NOT_PRIME_1
182    Thus there is a contradiction.
183*)
184val square_free_prime = store_thm(
185  "square_free_prime",
186  ``!n. prime n ==> square_free n``,
187  rw_tac std_ss[square_free_def] >>
188  spose_not_then strip_assume_tac >>
189  `p * p = p` by metis_tac[prime_def, MULT_EQ_1] >>
190  metis_tac[SQ_EQ_SELF, NOT_PRIME_0, NOT_PRIME_1]);
191
192(* Overload square-free filter of a set *)
193val _ = overload_on("sq_free", ``\s. {n | n IN s /\ square_free n}``);
194
195(* Overload non-square-free filter of a set *)
196val _ = overload_on("non_sq_free", ``\s. {n | n IN s /\ ~(square_free n)}``);
197
198(* Theorem: n IN sq_free s <=> n IN s /\ square_free n *)
199(* Proof: by notation. *)
200val sq_free_element = store_thm(
201  "sq_free_element",
202  ``!s n. n IN sq_free s <=> n IN s /\ square_free n``,
203  rw[]);
204
205(* Theorem: sq_free s SUBSET s *)
206(* Proof: by SUBSET_DEF *)
207val sq_free_subset = store_thm(
208  "sq_free_subset",
209  ``!s. sq_free s SUBSET s``,
210  rw[SUBSET_DEF]);
211
212(* Theorem: FINITE s ==> FINITE (sq_free s) *)
213(* Proof: by sq_free_subset, SUBSET_FINITE *)
214val sq_free_finite = store_thm(
215  "sq_free_finite",
216  ``!s. FINITE s ==> FINITE (sq_free s)``,
217  metis_tac[sq_free_subset, SUBSET_FINITE]);
218
219(* Theorem: n IN non_sq_free s <=> n IN s /\ ~(square_free n) *)
220(* Proof: by notation. *)
221val non_sq_free_element = store_thm(
222  "non_sq_free_element",
223  ``!s n. n IN non_sq_free s <=> n IN s /\ ~(square_free n)``,
224  rw[]);
225
226(* Theorem: non_sq_free s SUBSET s *)
227(* Proof: by SUBSET_DEF *)
228val non_sq_free_subset = store_thm(
229  "non_sq_free_subset",
230  ``!s. non_sq_free s SUBSET s``,
231  rw[SUBSET_DEF]);
232
233(* Theorem: FINITE s ==> FINITE (non_sq_free s) *)
234(* Proof: by non_sq_free_subset, SUBSET_FINITE *)
235val non_sq_free_finite = store_thm(
236  "non_sq_free_finite",
237  ``!s. FINITE s ==> FINITE (non_sq_free s)``,
238  metis_tac[non_sq_free_subset, SUBSET_FINITE]);
239
240(* Theorem: (s = (sq_free s) UNION (non_sq_free s)) /\ ((sq_free s) INTER (non_sq_free s) = {}) *)
241(* Proof:
242   This is to show:
243   (1) s = (sq_free s) UNION (non_sq_free s)
244       True by EXTENSION, IN_UNION.
245   (2) (sq_free s) INTER (non_sq_free s) = {}
246       True by EXTENSION, IN_INTER
247*)
248val sq_free_split = store_thm(
249  "sq_free_split",
250  ``!s. (s = (sq_free s) UNION (non_sq_free s)) /\ ((sq_free s) INTER (non_sq_free s) = {})``,
251  (rw[EXTENSION] >> metis_tac[]));
252
253(* Theorem: s = (sq_free s) UNION (non_sq_free s) *)
254(* Proof: extract from sq_free_split. *)
255val sq_free_union = save_thm("sq_free_union", sq_free_split |> SPEC_ALL |> CONJUNCT1 |> GEN_ALL);
256(* val sq_free_union = |- !s. s = sq_free s UNION non_sq_free s: thm *)
257
258(* Theorem: (sq_free s) INTER (non_sq_free s) = {} *)
259(* Proof: extract from sq_free_split. *)
260val sq_free_inter = save_thm("sq_free_inter", sq_free_split |> SPEC_ALL |> CONJUNCT2 |> GEN_ALL);
261(* val sq_free_inter = |- !s. sq_free s INTER non_sq_free s = {}: thm *)
262
263(* Theorem: DISJOINT (sq_free s) (non_sq_free s) *)
264(* Proof: by DISJOINT_DEF, sq_free_inter. *)
265val sq_free_disjoint = store_thm(
266  "sq_free_disjoint",
267  ``!s. DISJOINT (sq_free s) (non_sq_free s)``,
268  rw_tac std_ss[DISJOINT_DEF, sq_free_inter]);
269
270(* ------------------------------------------------------------------------- *)
271(* Prime Divisors of a Number and Partitions of Square-free Set              *)
272(* ------------------------------------------------------------------------- *)
273
274(* Define the prime divisors of a number *)
275val prime_factors_def = Define`
276    prime_factors n = {p | prime p /\ p IN (divisors n)}
277`;
278(* prime_divisors is used in triangle.hol *)
279
280(* Theorem: p IN prime_factors n <=> prime p /\ p <= n /\ p divides n *)
281(* Proof: by prime_factors_def, divisors_def *)
282val prime_factors_element = store_thm(
283  "prime_factors_element",
284  ``!n p. p IN prime_factors n <=> prime p /\ p <= n /\ p divides n``,
285  rw[prime_factors_def, divisors_def]);
286
287(* Theorem: (prime_factors n) SUBSET (divisors n) *)
288(* Proof:
289       p IN (prime_factors n)
290   ==> p IN (divisors n)              by prime_factors_def
291   Hence (prime_factors n) SUBSET (divisors n)   by SUBSET_DEF
292*)
293val prime_factors_subset = store_thm(
294  "prime_factors_subset",
295  ``!n. (prime_factors n) SUBSET (divisors n)``,
296  rw[prime_factors_def, SUBSET_DEF]);
297
298(* Theorem: FINITE (prime_factors n) *)
299(* Proof:
300   Since (prime_factors n) SUBSET (divisors n)    by prime_factors_subset
301     and FINITE (divisors n)           by divisors_finite
302    Thus FINITE (prime_factors n)                 by SUBSET_FINITE
303*)
304val prime_factors_finite = store_thm(
305  "prime_factors_finite",
306  ``!n. FINITE (prime_factors n)``,
307  metis_tac[prime_factors_subset, divisors_finite, SUBSET_FINITE]);
308
309(* Overload even square-free filter of a set *)
310val _ = overload_on("even_sq_free", ``\s. {n | n IN (sq_free s) /\ EVEN (CARD (prime_factors n))}``);
311
312(* Overload odd square-free filter of a set *)
313val _ = overload_on("odd_sq_free", ``\s. {n | n IN (sq_free s) /\ ODD (CARD (prime_factors n))}``);
314
315(* Theorem: n IN even_sq_free s <=> n IN s /\ square_free n /\ EVEN (CARD (prime_factors n)) *)
316(* Proof: by notation. *)
317val even_sq_free_element = store_thm(
318  "even_sq_free_element",
319  ``!s n. n IN even_sq_free s <=> n IN s /\ square_free n /\ EVEN (CARD (prime_factors n))``,
320  (rw[] >> metis_tac[]));
321
322(* Theorem: even_sq_free s SUBSET s *)
323(* Proof: by SUBSET_DEF *)
324val even_sq_free_subset = store_thm(
325  "even_sq_free_subset",
326  ``!s. even_sq_free s SUBSET s``,
327  rw[SUBSET_DEF]);
328
329(* Theorem: FINITE s ==> FINITE (even_sq_free s) *)
330(* Proof: by even_sq_free_subset, SUBSET_FINITE *)
331val even_sq_free_finite = store_thm(
332  "even_sq_free_finite",
333  ``!s. FINITE s ==> FINITE (even_sq_free s)``,
334  metis_tac[even_sq_free_subset, SUBSET_FINITE]);
335
336(* Theorem: n IN odd_sq_free s <=> n IN s /\ square_free n /\ ODD (CARD (prime_factors n)) *)
337(* Proof: by notation. *)
338val odd_sq_free_element = store_thm(
339  "odd_sq_free_element",
340  ``!s n. n IN odd_sq_free s <=> n IN s /\ square_free n /\ ODD (CARD (prime_factors n))``,
341  (rw[] >> metis_tac[]));
342
343(* Theorem: odd_sq_free s SUBSET s *)
344(* Proof: by SUBSET_DEF *)
345val odd_sq_free_subset = store_thm(
346  "odd_sq_free_subset",
347  ``!s. odd_sq_free s SUBSET s``,
348  rw[SUBSET_DEF]);
349
350(* Theorem: FINITE s ==> FINITE (odd_sq_free s) *)
351(* Proof: by odd_sq_free_subset, SUBSET_FINITE *)
352val odd_sq_free_finite = store_thm(
353  "odd_sq_free_finite",
354  ``!s. FINITE s ==> FINITE (odd_sq_free s)``,
355  metis_tac[odd_sq_free_subset, SUBSET_FINITE]);
356
357(* Theorem: (sq_free s = (even_sq_free s) UNION (odd_sq_free s)) /\
358            ((even_sq_free s) INTER (odd_sq_free s) = {}) *)
359(* Proof:
360   This is to show:
361   (1) sq_free s = even_sq_free s UNION odd_sq_free s
362       True by EXTENSION, IN_UNION, EVEN_ODD.
363   (2) even_sq_free s INTER odd_sq_free s = {}
364       True by EXTENSION, IN_INTER, EVEN_ODD.
365*)
366val sq_free_split_even_odd = store_thm(
367  "sq_free_split_even_odd",
368  ``!s. (sq_free s = (even_sq_free s) UNION (odd_sq_free s)) /\
369       ((even_sq_free s) INTER (odd_sq_free s) = {})``,
370  (rw[EXTENSION] >> metis_tac[EVEN_ODD]));
371
372(* Theorem: sq_free s = (even_sq_free s) UNION (odd_sq_free s) *)
373(* Proof: extract from sq_free_split_even_odd. *)
374val sq_free_union_even_odd =
375    save_thm("sq_free_union_even_odd", sq_free_split_even_odd |> SPEC_ALL |> CONJUNCT1 |> GEN_ALL);
376(* val sq_free_union_even_odd =
377   |- !s. sq_free s = even_sq_free s UNION odd_sq_free s: thm *)
378
379(* Theorem: (even_sq_free s) INTER (odd_sq_free s) = {} *)
380(* Proof: extract from sq_free_split_even_odd. *)
381val sq_free_inter_even_odd =
382    save_thm("sq_free_inter_even_odd", sq_free_split_even_odd |> SPEC_ALL |> CONJUNCT2 |> GEN_ALL);
383(* val sq_free_inter_even_odd =
384   |- !s. even_sq_free s INTER odd_sq_free s = {}: thm *)
385
386(* Theorem: DISJOINT (even_sq_free s) (odd_sq_free s) *)
387(* Proof: by DISJOINT_DEF, sq_free_inter_even_odd. *)
388val sq_free_disjoint_even_odd = store_thm(
389  "sq_free_disjoint_even_odd",
390  ``!s. DISJOINT (even_sq_free s) (odd_sq_free s)``,
391  rw_tac std_ss[DISJOINT_DEF, sq_free_inter_even_odd]);
392
393(* ------------------------------------------------------------------------- *)
394(* Less Divisors of a number.                                                *)
395(* ------------------------------------------------------------------------- *)
396
397(* Overload the set of divisors less than n *)
398val _ = overload_on("less_divisors", ``\n. {x | x IN (divisors n) /\ x <> n}``);
399
400(* Theorem: x IN (less_divisors n) <=> (x < n /\ x divides n) *)
401(* Proof: by divisors_element. *)
402val less_divisors_element = store_thm(
403  "less_divisors_element",
404  ``!n x. x IN (less_divisors n) <=> (x < n /\ x divides n)``,
405  rw[divisors_element, EQ_IMP_THM]);
406
407(* Theorem: less_divisors 0 = {} *)
408(* Proof: by divisors_element. *)
409val less_divisors_0 = store_thm(
410  "less_divisors_0",
411  ``less_divisors 0 = {}``,
412  rw[divisors_element]);
413
414(* Theorem: less_divisors 1 = {} *)
415(* Proof: by divisors_element. *)
416val less_divisors_1 = store_thm(
417  "less_divisors_1",
418  ``less_divisors 1 = {}``,
419  rw[divisors_element]);
420
421(* Theorem: (less_divisors n) SUBSET (divisors n) *)
422(* Proof: by SUBSET_DEF *)
423val less_divisors_subset_divisors = store_thm(
424  "less_divisors_subset_divisors",
425  ``!n. (less_divisors n) SUBSET (divisors n)``,
426  rw[SUBSET_DEF]);
427
428(* Theorem: FINITE (less_divisors n) *)
429(* Proof:
430   Since (less_divisors n) SUBSET (divisors n)   by less_divisors_subset_divisors
431     and FINITE (divisors n)                     by divisors_finite
432      so FINITE (proper_divisors n)              by SUBSET_FINITE
433*)
434val less_divisors_finite = store_thm(
435  "less_divisors_finite",
436  ``!n. FINITE (less_divisors n)``,
437  metis_tac[divisors_finite, less_divisors_subset_divisors, SUBSET_FINITE]);
438
439(* Theorem: prime n ==> (less_divisors n = {1}) *)
440(* Proof:
441   Since prime n
442     ==> !b. b divides n ==> (b = n) \/ (b = 1)   by prime_def
443   But (less_divisors n) excludes n               by less_divisors_element
444   and 1 < n                                      by ONE_LT_PRIME
445   Hence less_divisors n = {1}
446*)
447val less_divisors_prime = store_thm(
448  "less_divisors_prime",
449  ``!n. prime n ==> (less_divisors n = {1})``,
450  rpt strip_tac >>
451  `!b. b divides n ==> (b = n) \/ (b = 1)` by metis_tac[prime_def] >>
452  rw[less_divisors_element, EXTENSION, EQ_IMP_THM] >| [
453    `x <> n` by decide_tac >>
454    metis_tac[],
455    rw[ONE_LT_PRIME]
456  ]);
457
458(* Theorem: 1 < n ==> 1 IN (less_divisors n) *)
459(* Proof:
460       1 IN (less_divisors n)
461   <=> 1 < n /\ 1 divides n     by less_divisors_element
462   <=> T                        by ONE_DIVIDES_ALL
463*)
464val less_divisors_has_one = store_thm(
465  "less_divisors_has_one",
466  ``!n. 1 < n ==> 1 IN (less_divisors n)``,
467  rw[less_divisors_element]);
468
469(* Theorem: x IN (less_divisors n) ==> 0 < x *)
470(* Proof:
471   Since x IN (less_divisors n),
472         x < n /\ x divides n                 by less_divisors_element
473   By contradiction, if x = 0, then n = 0     by ZERO_DIVIDES
474   This contradicts x < n.
475*)
476val less_divisors_nonzero = store_thm(
477  "less_divisors_nonzero",
478  ``!n x. x IN (less_divisors n) ==> 0 < x``,
479  metis_tac[less_divisors_element, ZERO_DIVIDES, NOT_ZERO_LT_ZERO]);
480
481(* Theorem: 0 < n ==> !d. 1 < d /\ d IN (less_divisors n) ==> (n DIV d) IN (less_divisors n) *)
482(* Proof:
483      d IN (less_divisors n)
484  ==> d IN (divisors n)                   by notation
485  ==> (n DIV d) IN (divisors n)           by divisors_has_cofactor
486   Also n DIV d < n                       by DIV_LESS, 0 < n /\ 1 < d
487   thus n DIV d <> n                      by LESS_NOT_EQ
488  Hence (n DIV d) IN (less_divisors n)    by notation
489*)
490val less_divisors_has_cofactor = store_thm(
491  "less_divisors_has_cofactor",
492  ``!n. 0 < n ==> !d. 1 < d /\ d IN (less_divisors n) ==> (n DIV d) IN (less_divisors n)``,
493  rw[divisors_has_cofactor, DIV_LESS, LESS_NOT_EQ]);
494
495(* ------------------------------------------------------------------------- *)
496(* Proper Divisors of a number.                                              *)
497(* ------------------------------------------------------------------------- *)
498
499(* Overload the set of proper divisors of n *)
500val _ = overload_on("proper_divisors", ``\n. {x | x IN (divisors n) /\ x <> 1 /\ x <> n}``);
501
502(* Theorem: x IN (proper_divisors n) <=> (1 < x /\ x < n /\ x divides n) *)
503(* Proof:
504   Since x IN (divisors n)
505     ==> x <= n /\ x divides n       by divisors_element
506   Since x <= n but x <> n, x < n.
507   If x = 0, x divides n ==> n = 0   by ZERO_DIVIDES
508   But x <> n, so x <> 0.
509   With x <> 0 /\ x <> 1 ==> 1 < x.
510*)
511val proper_divisors_element = store_thm(
512  "proper_divisors_element",
513  ``!n x. x IN (proper_divisors n) <=> (1 < x /\ x < n /\ x divides n)``,
514  rw[divisors_element, EQ_IMP_THM] >>
515  `x <> 0` by metis_tac[ZERO_DIVIDES] >>
516  decide_tac);
517
518(* Theorem: proper_divisors 0 = {} *)
519(* Proof: by proper_divisors_element. *)
520val proper_divisors_0 = store_thm(
521  "proper_divisors_0",
522  ``proper_divisors 0 = {}``,
523  rw[proper_divisors_element, EXTENSION]);
524
525(* Theorem: proper_divisors 1 = {} *)
526(* Proof: by proper_divisors_element. *)
527val proper_divisors_1 = store_thm(
528  "proper_divisors_1",
529  ``proper_divisors 1 = {}``,
530  rw[proper_divisors_element, EXTENSION]);
531
532(* Theorem: (proper_divisors n) SUBSET (less_divisors n) *)
533(* Proof: by SUBSET_DEF *)
534val proper_divisors_subset = store_thm(
535  "proper_divisors_subset",
536  ``!n. (proper_divisors n) SUBSET (less_divisors n)``,
537  rw[SUBSET_DEF]);
538
539(* Theorem: FINITE (proper_divisors n) *)
540(* Proof:
541   Since (proper_divisors n) SUBSET (less_divisors n) by proper_divisors_subset
542     and FINITE (less_divisors n)                     by less_divisors_finite
543      so FINITE (proper_divisors n)                   by SUBSET_FINITE
544*)
545val proper_divisors_finite = store_thm(
546  "proper_divisors_finite",
547  ``!n. FINITE (proper_divisors n)``,
548  metis_tac[less_divisors_finite, proper_divisors_subset, SUBSET_FINITE]);
549
550(* Theorem: 1 NOTIN (proper_divisors n) *)
551(* Proof: proper_divisors_element *)
552val proper_divisors_not_one = store_thm(
553  "proper_divisors_not_one",
554  ``!n. 1 NOTIN (proper_divisors n)``,
555  rw[proper_divisors_element]);
556
557(* Theorem: proper_divisors n = (less_divisors n) DELETE 1 *)
558(* Proof:
559      proper_divisors n
560    = {x | x IN (divisors n) /\ x <> 1 /\ x <> n}   by notation
561    = {x | x IN (divisors n) /\ x <> n} DELETE 1    by IN_DELETE
562    = (less_divisors n) DELETE 1
563*)
564val proper_divisors_by_less_divisors = store_thm(
565  "proper_divisors_by_less_divisors",
566  ``!n. proper_divisors n = (less_divisors n) DELETE 1``,
567  rw[divisors_element, EXTENSION, EQ_IMP_THM]);
568
569(* Theorem: prime n ==> (proper_divisors n = {}) *)
570(* Proof:
571      proper_divisors n
572    = (less_divisors n) DELETE 1  by proper_divisors_by_less_divisors
573    = {1} DELETE 1                by less_divisors_prime, prime n
574    = {}                          by SING_DELETE
575*)
576val proper_divisors_prime = store_thm(
577  "proper_divisors_prime",
578  ``!n. prime n ==> (proper_divisors n = {})``,
579  rw[proper_divisors_by_less_divisors, less_divisors_prime]);
580
581(* Theorem: d IN (proper_divisors n) ==> (n DIV d) IN (proper_divisors n) *)
582(* Proof:
583   Let e = n DIV d.
584   Since d IN (proper_divisors n)
585     ==> 1 < d /\ d < n                   by proper_divisors_element
586     and d IN (less_divisors n)           by proper_divisors_subset
587      so e IN (less_divisors n)           by less_divisors_has_cofactor
588     and 0 < e                            by less_divisors_nonzero
589   Since d divides n                      by less_divisors_element
590      so n = e * d                        by DIV_MULT_EQ, 0 < d
591    thus e <> 1 since n <> d              by MULT_LEFT_1
592    With 0 < e /\ e <> 1
593     ==> e IN (proper_divisors n)         by proper_divisors_by_less_divisors, IN_DELETE
594*)
595val proper_divisors_has_cofactor = store_thm(
596  "proper_divisors_has_cofactor",
597  ``!n d. d IN (proper_divisors n) ==> (n DIV d) IN (proper_divisors n)``,
598  rpt strip_tac >>
599  qabbrev_tac `e = n DIV d` >>
600  `1 < d /\ d < n` by metis_tac[proper_divisors_element] >>
601  `d IN (less_divisors n)` by metis_tac[proper_divisors_subset, SUBSET_DEF] >>
602  `e IN (less_divisors n)` by rw[less_divisors_has_cofactor, Abbr`e`] >>
603  `0 < e` by metis_tac[less_divisors_nonzero] >>
604  `0 < d /\ n <> d` by decide_tac >>
605  `e <> 1` by metis_tac[less_divisors_element, DIV_MULT_EQ, MULT_LEFT_1] >>
606  metis_tac[proper_divisors_by_less_divisors, IN_DELETE]);
607
608(* Theorem: (proper_divisors n) <> {} ==> 1 < MIN_SET (proper_divisors n) *)
609(* Proof:
610   Let s = proper_divisors n.
611   Since !x. x IN s ==> 1 < x        by proper_divisors_element
612     But MIN_SET s IN s              by MIN_SET_IN_SET
613   Hence 1 < MIN_SET s               by above
614*)
615val proper_divisors_min_gt_1 = store_thm(
616  "proper_divisors_min_gt_1",
617  ``!n. (proper_divisors n) <> {} ==> 1 < MIN_SET (proper_divisors n)``,
618  metis_tac[MIN_SET_IN_SET, proper_divisors_element]);
619
620(* Theorem: (proper_divisors n) <> {} ==>
621            (MAX_SET (proper_divisors n) = n DIV (MIN_SET (proper_divisors n))) /\
622            (MIN_SET (proper_divisors n) = n DIV (MAX_SET (proper_divisors n)))     *)
623(* Proof:
624   Let s = proper_divisors n, b = MIN_SET s.
625   By MAX_SET_ELIM, this is to show:
626   (1) FINITE s, true                     by proper_divisors_finite
627   (2) s <> {} /\ x IN s /\ !y. y IN s ==> y <= x ==> x = n DIV b /\ b = n DIV x
628       Note s <> {} ==> n <> 0            by proper_divisors_0
629        Let m = n DIV b.
630       Note n DIV x IN s                  by proper_divisors_has_cofactor, 0 < n, 1 < b.
631       Also b IN s /\ b <= x              by MIN_SET_IN_SET, s <> {}
632       thus 1 < b                         by proper_divisors_min_gt_1
633         so m IN s                        by proper_divisors_has_cofactor, 0 < n, 1 < x.
634         or 1 < m                         by proper_divisors_nonzero
635        and m <= x                        by implication, x = MAX_SET s.
636       Thus n DIV x <= n DIV m            by DIV_LE_MONOTONE_REVERSE [1], 0 < x, 0 < m.
637        But n DIV m
638          = n DIV (n DIV b) = b           by divide_by_cofactor, b divides n.
639         so n DIV x <= b                  by [1]
640      Since b <= n DIV x                  by MIN_SET_PROPERTY, b = MIN_SET s, n DIV x IN s.
641         so n DIV x = b                   by LESS_EQUAL_ANTISYM (gives second subgoal)
642      Hence m = n DIV b
643              = n DIV (n DIV x) = x       by divide_by_cofactor, x divides n (gives first subgoal)
644*)
645val proper_divisors_max_min = store_thm(
646  "proper_divisors_max_min",
647  ``!n. (proper_divisors n) <> {} ==>
648       (MAX_SET (proper_divisors n) = n DIV (MIN_SET (proper_divisors n))) /\
649       (MIN_SET (proper_divisors n) = n DIV (MAX_SET (proper_divisors n)))``,
650  ntac 2 strip_tac >>
651  qabbrev_tac `s = proper_divisors n` >>
652  qabbrev_tac `b = MIN_SET s` >>
653  DEEP_INTRO_TAC MAX_SET_ELIM >>
654  strip_tac >-
655  rw[proper_divisors_finite, Abbr`s`] >>
656  ntac 3 strip_tac >>
657  `n <> 0` by metis_tac[proper_divisors_0] >>
658  `b IN s /\ b <= x` by rw[MIN_SET_IN_SET, Abbr`b`] >>
659  `1 < b` by rw[proper_divisors_min_gt_1, Abbr`s`, Abbr`b`] >>
660  `0 < n /\ 1 < x` by decide_tac >>
661  qabbrev_tac `m = n DIV b` >>
662  `m IN s /\ (n DIV x) IN s` by rw[proper_divisors_has_cofactor, Abbr`m`, Abbr`s`] >>
663  `1 < m` by metis_tac[proper_divisors_element] >>
664  `0 < x /\ 0 < m` by decide_tac >>
665  `n DIV x <= n DIV m` by rw[DIV_LE_MONOTONE_REVERSE] >>
666  `b divides n /\ x divides n` by metis_tac[proper_divisors_element] >>
667  `n DIV m = b` by rw[divide_by_cofactor, Abbr`m`] >>
668  `b <= n DIV x` by rw[MIN_SET_PROPERTY, Abbr`b`] >>
669  `b = n DIV x` by rw[LESS_EQUAL_ANTISYM] >>
670  `m = x` by rw[divide_by_cofactor, Abbr`m`] >>
671  decide_tac);
672
673(* This is a milestone theorem. *)
674
675(* ------------------------------------------------------------------------- *)
676(* Useful Properties of Less Divisors                                        *)
677(* ------------------------------------------------------------------------- *)
678
679(* Theorem: 1 < n ==> (MIN_SET (less_divisors n) = 1) *)
680(* Proof:
681   Let s = less_divisors n.
682   Since 1 < n ==> 1 IN s         by less_divisors_has_one
683      so s <> {}                  by MEMBER_NOT_EMPTY
684     and !y. y IN s ==> 0 < y     by less_divisors_nonzero
685      or !y. y IN s ==> 1 <= y    by LESS_EQ
686   Hence 1 = MIN_SET s            by MIN_SET_TEST
687*)
688val less_divisors_min = store_thm(
689  "less_divisors_min",
690  ``!n. 1 < n ==> (MIN_SET (less_divisors n) = 1)``,
691  metis_tac[less_divisors_has_one, MEMBER_NOT_EMPTY,
692             MIN_SET_TEST, less_divisors_nonzero, LESS_EQ, ONE]);
693
694(* Theorem: MAX_SET (less_divisors n) <= n DIV 2 *)
695(* Proof:
696   Let s = less_divisors n, m = MAX_SET s.
697   If s = {},
698      Then m = MAX_SET {} = 0          by MAX_SET_EMPTY
699       and 0 <= n DIV 2 is trivial.
700   If s <> {},
701      Then n <> 0 /\ n <> 1            by less_divisors_0, less_divisors_1
702   Note 1 IN s                         by less_divisors_has_one
703   Consider t = s DELETE 1.
704   Then t = proper_divisors n          by proper_divisors_by_less_divisors
705   If t = {},
706      Then s = {1}                     by DELETE_EQ_SING
707       and m = 1                       by SING_DEF, IN_SING (same as MAX_SET_SING)
708     Since 2 <= n                      by 1 < n
709      thus n DIV n <= n DIV 2          by DIV_LE_MONOTONE_REVERSE
710        or n DIV n = 1 = m <= n DIV 2  by DIVMOD_ID, 0 < n
711   If t <> {},
712      Let b = MIN_SET t
713      Then MAX_SET t = n DIV b         by proper_divisors_max_min, t <> {}
714     Since MIN_SET s = 1               by less_divisors_min, 1 < n
715       and FINITE s                    by less_divisors_finite
716       and s <> {1}                    by DELETE_EQ_SING
717      thus m = MAX_SET t               by MAX_SET_DELETE, s <> {1}
718
719       Now 1 < b                       by proper_divisors_min_gt_1
720        so 2 <= b                      by LESS_EQ, 1 < b
721     Hence n DIV b <= n DIV 2          by DIV_LE_MONOTONE_REVERSE
722       or        m <= n DIV 2          by m = MAX_SET t = n DIV b
723*)
724
725Theorem less_divisors_max:
726  !n. MAX_SET (less_divisors n) <= n DIV 2
727Proof
728  rpt strip_tac >>
729  qabbrev_tac `s = less_divisors n` >>
730  qabbrev_tac `m = MAX_SET s` >>
731  Cases_on `s = {}` >- rw[MAX_SET_EMPTY, Abbr`m`] >>
732  `n <> 0 /\ n <> 1` by metis_tac[less_divisors_0, less_divisors_1] >>
733  `1 < n` by decide_tac >>
734  `1 IN s` by rw[less_divisors_has_one, Abbr`s`] >>
735  qabbrev_tac `t = proper_divisors n` >>
736  `t = s DELETE 1`  by rw[proper_divisors_by_less_divisors, Abbr`t`, Abbr`s`] >>
737  Cases_on `t = {}` >| [
738    `s = {1}` by rfs[] >>
739    `m = 1` by rw[MAX_SET_SING, Abbr`m`] >>
740    `(2 <= n) /\ (0 < 2) /\ (0 < n) /\ (n DIV n = 1)` by rw[] >>
741    metis_tac[DIV_LE_MONOTONE_REVERSE],
742    qabbrev_tac `b = MIN_SET t` >>
743    `MAX_SET t = n DIV b` by metis_tac[proper_divisors_max_min] >>
744    `MIN_SET s = 1` by rw[less_divisors_min, Abbr`s`] >>
745    `FINITE s` by rw[less_divisors_finite, Abbr`s`] >>
746    `s <> {1}` by metis_tac[DELETE_EQ_SING] >>
747    `m = MAX_SET t` by metis_tac[MAX_SET_DELETE] >>
748    `1 < b` by rw[proper_divisors_min_gt_1, Abbr`b`, Abbr`t`] >>
749    `2 <= b /\ (0 < b) /\ (0 < 2)` by decide_tac >>
750    `n DIV b <= n DIV 2` by rw[DIV_LE_MONOTONE_REVERSE] >>
751    decide_tac
752  ]
753QED
754
755(* Theorem: (less_divisors n) SUBSET (natural (n DIV 2)) *)
756(* Proof:
757   Let s = less_divisors n
758   If n = 0 or n - 1,
759   Then s = {}                        by less_divisors_0, less_divisors_1
760    and {} SUBSET t, for any t.       by EMPTY_SUBSET
761   If n <> 0 and n <> 1, 1 < n.
762   Note FINITE s                      by less_divisors_finite
763    and x IN s ==> x <= MAX_SET s     by MAX_SET_PROPERTY, FINITE s
764    But MAX_SET s <= n DIV 2          by less_divisors_max
765   Thus x IN s ==> x <= n DIV 2       by LESS_EQ_TRANS
766   Note s <> {}                       by MEMBER_NOT_EMPTY, x IN s
767    and x IN s ==> MIN_SET s <= x     by MIN_SET_PROPERTY, s <> {}
768  Since 1 = MIN_SET s, 1 <= x         by less_divisors_min, 1 < n
769   Thus 0 < x <= n DIV 2              by LESS_EQ
770     or x IN (natural (n DIV 2))      by natural_element
771*)
772val less_divisors_subset_natural = store_thm(
773  "less_divisors_subset_natural",
774  ``!n. (less_divisors n) SUBSET (natural (n DIV 2))``,
775  rpt strip_tac >>
776  qabbrev_tac `s = less_divisors n` >>
777  qabbrev_tac `m = n DIV 2` >>
778  Cases_on `(n = 0) \/ (n = 1)` >-
779  metis_tac[less_divisors_0, less_divisors_1, EMPTY_SUBSET] >>
780  `1 < n` by decide_tac >>
781  rw_tac std_ss[SUBSET_DEF] >>
782  `s <> {}` by metis_tac[MEMBER_NOT_EMPTY] >>
783  `FINITE s` by rw[less_divisors_finite, Abbr`s`] >>
784  `x <= MAX_SET s` by rw[MAX_SET_PROPERTY] >>
785  `MIN_SET s <= x` by rw[MIN_SET_PROPERTY] >>
786  `MAX_SET s <= m` by rw[less_divisors_max, Abbr`s`, Abbr`m`] >>
787  `MIN_SET s = 1` by rw[less_divisors_min, Abbr`s`] >>
788  `0 < x /\ x <= m` by decide_tac >>
789  rw[natural_element]);
790
791(* ------------------------------------------------------------------------- *)
792(* Properties of Summation equals Perfect Power                              *)
793(* ------------------------------------------------------------------------- *)
794
795(* Idea for the theorem below (for m = n DIV 2 when applied in bounds):
796      p * (p ** m - 1) / (p - 1)
797   <  p * p ** m / (p - 1)        discard subtraction
798   <= p * p ** m / (p / 2)        replace by smaller denominator
799    = 2 * p ** m                  double division and cancel p
800   or p * (p ** m - 1) < (p - 1) * 2 * p ** m
801*)
802
803(* Theorem: 1 < p ==> !n. p * (p ** n - 1) < (p - 1) * (2 * p ** n) *)
804(* Proof:
805   Let q = p ** n
806   Then 1 <= q                       by ONE_LE_EXP, 0 < p
807     so p <= p * q                   by LE_MULT_LCANCEL, p <> 0
808   Also 1 < p ==> 2 <= p             by LESS_EQ
809     so 2 * q <= p * q               by LE_MULT_RCANCEL, q <> 0
810   Thus   LHS
811        = p * (q - 1)
812        = p * q - p                  by LEFT_SUB_DISTRIB
813    And   RHS
814        = (p - 1) * (2 * q)
815        = p * (2 * q) - 2 * q        by RIGHT_SUB_DISTRIB
816        = 2 * (p * q) - 2 * q        by MULT_ASSOC, MULT_COMM
817        = (p * q + p * q) - 2 * q    by TIMES2
818        = (p * q - p + p + p * q) - 2 * q  by SUB_ADD, p <= p * q
819        = LHS + p + p * q - 2 * q    by above
820        = LHS + p + (p * q - 2 * q)  by LESS_EQ_ADD_SUB, 2 * q <= p * q
821        = LHS + p + (p - 2) * q      by RIGHT_SUB_DISTRIB
822
823    Since 0 < p                      by 1 < p
824      and 0 <= (p - 2) * q           by 2 <= p
825    Hence LHS < RHS                  by discarding positive terms
826*)
827val perfect_power_special_inequality = store_thm(
828  "perfect_power_special_inequality",
829  ``!p. 1 < p ==> !n. p * (p ** n - 1) < (p - 1) * (2 * p ** n)``,
830  rpt strip_tac >>
831  qabbrev_tac `q = p ** n` >>
832  `p <> 0 /\ 2 <= p` by decide_tac >>
833  `1 <= q` by rw[ONE_LE_EXP, Abbr`q`] >>
834  `p <= p * q` by rw[] >>
835  `2 * q <= p * q` by rw[] >>
836  qabbrev_tac `l = p * (q - 1)` >>
837  qabbrev_tac `r = (p - 1) * (2 * q)` >>
838  `l = p * q - p` by rw[Abbr`l`] >>
839  `r = p * (2 * q) - 2 * q` by rw[Abbr`r`] >>
840  `_ = 2 * (p * q) - 2 * q` by rw[] >>
841  `_ = (p * q + p * q) - 2 * q` by rw[] >>
842  `_ = (p * q - p + p + p * q) - 2 * q` by rw[] >>
843  `_ = l + p + p * q - 2 * q` by rw[] >>
844  `_ = l + p + (p * q - 2 * q)` by rw[] >>
845  `_ = l + p + (p - 2) * q` by rw[] >>
846  decide_tac);
847
848(* Theorem: 1 < p /\ 1 < n ==>
849            p ** (n DIV 2) * p ** (n DIV 2) <= p ** n /\
850            2 * p ** (n DIV 2) <= p ** (n DIV 2) * p ** (n DIV 2) *)
851(* Proof:
852   Let m = n DIV 2, q = p ** m.
853   The goal becomes: q * q <= p ** n /\ 2 * q <= q * q.
854      Note 1 < p ==> 0 < p.
855   First goal: q * q <= p ** n
856      Then 0 < q                    by EXP_POS, 0 < p
857       and 2 * m <= n               by DIV_MULT_LE, 0 < 2.
858      thus p ** (2 * m) <= p ** n   by EXP_BASE_LE_MONO, 1 < p.
859     Since p ** (2 * m)
860         = p ** (m + m)             by TIMES2
861         = q * q                    by EXP_ADD
862      Thus q * q <= p ** n          by above
863
864   Second goal: 2 * q <= q * q
865     Since 1 < n, so 2 <= n         by LESS_EQ
866        so 2 DIV 2 <= n DIV 2       by DIV_LE_MONOTONE, 0 < 2.
867        or 1 <= m, i.e. 0 < m       by DIVMOD_ID, 0 < 2.
868      Thus 1 < q                    by ONE_LT_EXP, 1 < p, 0 < m.
869        so 2 <= q                   by LESS_EQ
870       and 2 * q <= q * q           by MULT_RIGHT_CANCEL, q <> 0.
871     Hence 2 * q <= p ** n          by LESS_EQ_TRANS
872*)
873val perfect_power_half_inequality_lemma = prove(
874  ``!p n. 1 < p /\ 1 < n ==>
875         p ** (n DIV 2) * p ** (n DIV 2) <= p ** n /\
876         2 * p ** (n DIV 2) <= p ** (n DIV 2) * p ** (n DIV 2)``,
877  ntac 3 strip_tac >>
878  qabbrev_tac `m = n DIV 2` >>
879  qabbrev_tac `q = p ** m` >>
880  strip_tac >| [
881    `0 < p /\ 0 < 2` by decide_tac >>
882    `0 < q /\ q <> 0` by rw[EXP_POS, Abbr`q`] >>
883    `2 * m <= n` by metis_tac[DIV_MULT_LE, MULT_COMM] >>
884    `p ** (2 * m) <= p ** n` by rw[EXP_BASE_LE_MONO] >>
885    `p ** (2 * m) = p ** (m + m)` by rw[] >>
886    `_ = q * q` by rw[EXP_ADD, Abbr`q`] >>
887    decide_tac,
888    `2 <= n /\ 0 < 2` by decide_tac >>
889    `1 <= m` by metis_tac[DIV_LE_MONOTONE, DIVMOD_ID] >>
890    `0 < m` by decide_tac >>
891    `1 < q` by rw[ONE_LT_EXP, Abbr`q`] >>
892    rw[]
893  ]);
894
895(* Theorem: 1 < p /\ 0 < n ==> 2 * p ** (n DIV 2) <= p ** n *)
896(* Proof:
897   Let m = n DIV 2, q = p ** m.
898   The goal becomes: 2 * q <= p ** n
899   If n = 1,
900      Then m = 0                    by ONE_DIV, 0 < 2.
901       and q = 1                    by EXP
902       and p ** n = p               by EXP_1
903     Since 1 < p ==> 2 <= p         by LESS_EQ
904     Hence 2 * q <= p = p ** n      by MULT_RIGHT_1
905   If n <> 1, 1 < n.
906      Then q * q <= p ** n /\
907           2 * q <= q * q           by perfect_power_half_inequality_lemma
908     Hence 2 * q <= p ** n          by LESS_EQ_TRANS
909*)
910val perfect_power_half_inequality_1 = store_thm(
911  "perfect_power_half_inequality_1",
912  ``!p n. 1 < p /\ 0 < n ==> 2 * p ** (n DIV 2) <= p ** n``,
913  rpt strip_tac >>
914  qabbrev_tac `m = n DIV 2` >>
915  qabbrev_tac `q = p ** m` >>
916  Cases_on `n = 1` >| [
917    `m = 0` by rw[Abbr`m`] >>
918    `(q = 1) /\ (p ** n = p)` by rw[Abbr`q`] >>
919    `2 <= p` by decide_tac >>
920    rw[],
921    `1 < n` by decide_tac >>
922    `q * q <= p ** n /\ 2 * q <= q * q` by rw[perfect_power_half_inequality_lemma, Abbr`q`, Abbr`m`] >>
923    decide_tac
924  ]);
925
926(* Theorem: 1 < p /\ 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) <= p ** n - 2 * p ** (n DIV 2) *)
927(* Proof:
928   Let m = n DIV 2, q = p ** m.
929   The goal becomes: (q - 2) * q <= p ** n - 2 * q
930   If n = 1,
931      Then m = 0                    by ONE_DIV, 0 < 2.
932       and q = 1                    by EXP
933       and p ** n = p               by EXP_1
934     Since 1 < p ==> 2 <= p         by LESS_EQ
935        or 0 <= p - 2               by SUB_LEFT_LESS_EQ
936     Hence (q - 2) * q = 0 <= p - 2
937   If n <> 1, 1 < n.
938      Then q * q <= p ** n /\ 2 * q <= q * q   by perfect_power_half_inequality_lemma
939      Thus q * q - 2 * q <= p ** n - 2 * q     by LE_SUB_RCANCEL, 2 * q <= q * q
940        or   (q - 2) * q <= p ** n - 2 * q     by RIGHT_SUB_DISTRIB
941*)
942val perfect_power_half_inequality_2 = store_thm(
943  "perfect_power_half_inequality_2",
944  ``!p n. 1 < p /\ 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) <= p ** n - 2 * p ** (n DIV 2)``,
945  rpt strip_tac >>
946  qabbrev_tac `m = n DIV 2` >>
947  qabbrev_tac `q = p ** m` >>
948  Cases_on `n = 1` >| [
949    `m = 0` by rw[Abbr`m`] >>
950    `(q = 1) /\ (p ** n = p)` by rw[Abbr`q`] >>
951    `0 <= p - 2 /\ (1 - 2 = 0)` by decide_tac >>
952    rw[],
953    `1 < n` by decide_tac >>
954    `q * q <= p ** n /\ 2 * q <= q * q` by rw[perfect_power_half_inequality_lemma, Abbr`q`, Abbr`m`] >>
955    decide_tac
956  ]);
957
958(* Already in pred_setTheory:
959SUM_IMAGE_SUBSET_LE;
960!f s t. FINITE s /\ t SUBSET s ==> SIGMA f t <= SIGMA f s: thm
961SUM_IMAGE_MONO_LESS_EQ;
962|- !s. FINITE s ==> (!x. x IN s ==> f x <= g x) ==> SIGMA f s <= SIGMA g s: thm
963*)
964
965(* Theorem: 1 < p ==> !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
966            (!n. 0 < n ==> n * (f n) <= p ** n) /\
967            (!n. 0 < n ==> p ** n - 2 * p ** (n DIV 2) < n * (f n)) *)
968(* Proof:
969   Step 1: prove a specific lemma for sum decomposition
970   Claim: !n. 0 < n ==> (divisors n DIFF {n}) SUBSET (natural (n DIV 2)) /\
971          (p ** n = SIGMA (\d. d * f d) (divisors n)) ==>
972          (p ** n = n * f n + SIGMA (\d. d * f d) (divisors n DIFF {n}))
973   Proof: Let s = divisors n, a = {n}, b = s DIFF a, m = n DIV 2.
974          Then b = less_divisors n        by EXTENSION,IN_DIFF
975           and b SUBSET (natural m)       by less_divisors_subset_natural
976          This gives the first part.
977          For the second part:
978          Note a SUBSET s                 by divisors_has_last, SUBSET_DEF
979           and b SUBSET s                 by DIFF_SUBSET
980          Thus s = b UNION a              by UNION_DIFF, a SUBSET s
981           and DISJOINT b a               by DISJOINT_DEF, EXTENSION
982           Now FINITE s                   by divisors_finite
983            so FINITE a /\ FINITE b       by SUBSET_FINITE, by a SUBSEt s /\ b SUBSET s
984
985               p ** n
986             = SIGMA (\d. d * f d) s              by implication
987             = SIGMA (\d. d * f d) (b UNION a)    by above, s = b UNION a
988             = SIGMA (\d. d * f d) b + SIGMA (\d. d * f d) a   by SUM_IMAGE_DISJOINT, FINITE a /\ FINITE b
989             = SIGMA (\d. d * f d) b + n * f n    by SUM_IMAGE_SING
990             = n * f n + SIGMA (\d. d * f d) b    by ADD_COMM
991          This gives the second part.
992
993   Step 2: Upper bound, to show: !n. 0 < n ==> n * f n <= p ** n
994           Let b = divisors n DIFF {n}
995           Since n * f n + SIGMA (\d. d * f d) b = p ** n    by lemma
996           Hence n * f n <= p ** n                           by 0 <= SIGMA (\d. d * f d) b
997
998   Step 3: Lower bound, to show: !n. 0 < n ==> p ** n - p ** (n DIV 2) <= n * f n
999           Let s = divisors n, a = {n}, b = s DIFF a, m = n DIV 2.
1000            Note b SUBSET (natural m) /\
1001                 (p ** n = n * f n + SIGMA (\d. d * f d) b)  by lemma
1002           Since FINITE (natural m)                          by natural_finite
1003            thus SIGMA (\d. d * f d) b
1004              <= SIGMA (\d. d * f d) (natural m)             by SUM_IMAGE_SUBSET_LE [1]
1005            Also !d. d IN (natural m) ==> 0 < d              by natural_element
1006             and !d. 0 < d ==> d * f d <= p ** d             by upper bound (Step 2)
1007            thus !d. d IN (natural m) ==> d * f d <= p ** d  by implication
1008           Hence SIGMA (\d. d * f d) (natural m)
1009              <= SIGMA (\d. p ** d) (natural m)              by SUM_IMAGE_MONO_LESS_EQ [2]
1010             Now 1 < p ==> 0 < p /\ (p - 1) <> 0             by arithmetic
1011
1012             (p - 1) * SIGMA (\d. d * f d) b
1013          <= (p - 1) * SIGMA (\d. d * f d) (natural m)       by LE_MULT_LCANCEL, [1]
1014          <= (p - 1) * SIGMA (\d. p ** d) (natural m)        by LE_MULT_LCANCEL, [2]
1015           = p * (p ** m - 1)                                by sigma_geometric_natural_eqn
1016           < (p - 1) * (2 * p ** m)                          by perfect_power_special_inequality
1017
1018             (p - 1) * SIGMA (\d. d * f d) b < (p - 1) * (2 * p ** m)   by LESS_EQ_LESS_TRANS
1019             or        SIGMA (\d. d * f d) b < 2 * p ** m               by LT_MULT_LCANCEL, (p - 1) <> 0
1020
1021            But 2 * p ** m <= p ** n                         by perfect_power_half_inequality_1, 1 < p, 0 < n
1022           Thus p ** n = p ** n - 2 * p ** m + 2 * p ** m    by SUB_ADD, 2 * p ** m <= p ** n
1023       Combinig with lemma,
1024           p ** n - 2 * p ** m + 2 * p ** m < n * f n + 2 * p ** m
1025             or         p ** n - 2 * p ** m < n * f n        by LESS_MONO_ADD_EQ, no condition
1026*)
1027Theorem sigma_eq_perfect_power_bounds_1:
1028  !p.
1029    1 < p ==>
1030    !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
1031        (!n. 0 < n ==> n * (f n) <= p ** n) /\
1032        (!n. 0 < n ==> p ** n - 2 * p ** (n DIV 2) < n * (f n))
1033Proof
1034  ntac 4 strip_tac >>
1035  ������n. 0 < n ==>
1036       (divisors n DIFF {n}) SUBSET (natural (n DIV 2)) /\
1037       (p ** n = SIGMA (\d. d * f d) (divisors n) ==>
1038        p ** n = n * f n + SIGMA (\d. d * f d) (divisors n DIFF {n}))���
1039    by (ntac 2 strip_tac >>
1040        qabbrev_tac `s = divisors n` >>
1041        qabbrev_tac `a = {n}` >>
1042        qabbrev_tac `b = s DIFF a` >>
1043        qabbrev_tac `m = n DIV 2` >>
1044        `b = less_divisors n` by rw[EXTENSION, Abbr`b`, Abbr`a`, Abbr`s`] >>
1045        `b SUBSET (natural m)` by metis_tac[less_divisors_subset_natural] >>
1046        strip_tac >- rw[] >>
1047        `a SUBSET s` by rw[divisors_has_last, SUBSET_DEF, Abbr`s`, Abbr`a`] >>
1048        `b SUBSET s` by rw[Abbr`b`] >>
1049        `s = b UNION a` by rw[UNION_DIFF, Abbr`b`] >>
1050        `DISJOINT b a`
1051          by (rw[DISJOINT_DEF, Abbr`b`, EXTENSION] >> metis_tac[]) >>
1052        `FINITE s` by rw[divisors_finite, Abbr`s`] >>
1053        `FINITE a /\ FINITE b` by metis_tac[SUBSET_FINITE] >>
1054        strip_tac >>
1055        `_ = SIGMA (\d. d * f d) (b UNION a)` by metis_tac[Abbr`s`] >>
1056        `_ = SIGMA (\d. d * f d) b + SIGMA (\d. d * f d) a`
1057          by rw[SUM_IMAGE_DISJOINT] >>
1058        `_ = SIGMA (\d. d * f d) b + n * f n` by rw[SUM_IMAGE_SING, Abbr`a`] >>
1059        rw[]) >>
1060  conj_asm1_tac >| [
1061    rpt strip_tac >>
1062    `p ** n = n * f n + SIGMA (\d. d * f d) (divisors n DIFF {n})` by rw[] >>
1063    decide_tac,
1064    rpt strip_tac >>
1065    qabbrev_tac `s = divisors n` >>
1066    qabbrev_tac `a = {n}` >>
1067    qabbrev_tac `b = s DIFF a` >>
1068    qabbrev_tac `m = n DIV 2` >>
1069    `b SUBSET (natural m) /\ (p ** n = n * f n + SIGMA (\d. d * f d) b)`
1070      by rw[Abbr`s`, Abbr`a`, Abbr`b`, Abbr`m`] >>
1071    `FINITE (natural m)` by rw[natural_finite] >>
1072    `SIGMA (\d. d * f d) b <= SIGMA (\d. d * f d) (natural m)`
1073      by rw[SUM_IMAGE_SUBSET_LE] >>
1074    `!d. d IN (natural m) ==> 0 < d` by rw[natural_element] >>
1075    `SIGMA (\d. d * f d) (natural m) <= SIGMA (\d. p ** d) (natural m)`
1076      by rw[SUM_IMAGE_MONO_LESS_EQ] >>
1077    `0 < p /\ (p - 1) <> 0` by decide_tac >>
1078    `(p - 1) * SIGMA (\d. p ** d) (natural m) = p * (p ** m - 1)`
1079      by rw[sigma_geometric_natural_eqn] >>
1080    `p * (p ** m - 1) < (p - 1) * (2 * p ** m)`
1081      by rw[perfect_power_special_inequality] >>
1082    `SIGMA (\d. d * f d) b < 2 * p ** m`
1083      by metis_tac[LE_MULT_LCANCEL, LESS_EQ_TRANS, LESS_EQ_LESS_TRANS,
1084                   LT_MULT_LCANCEL] >>
1085    `p ** n < n * f n + 2 * p ** m` by decide_tac >>
1086    `2 * p ** m <= p ** n` by rw[perfect_power_half_inequality_1, Abbr`m`] >>
1087    decide_tac
1088  ]
1089QED
1090
1091(* Theorem: 1 < p ==> !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
1092            (!n. 0 < n ==> n * (f n) <= p ** n) /\
1093            (!n. 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) < n * (f n)) *)
1094(* Proof:
1095   For the first goal: (!n. 0 < n ==> n * (f n) <= p ** n)
1096       True by sigma_eq_perfect_power_bounds_1.
1097   For the second goal: (!n. 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) < n * (f n))
1098       Let m = n DIV 2.
1099       Then p ** n - 2 * p ** m < n * (f n)                     by sigma_eq_perfect_power_bounds_1
1100        and (p ** m - 2) * p ** m <= p ** n - 2 * p ** m        by perfect_power_half_inequality_2
1101      Hence (p ** (n DIV 2) - 2) * p ** (n DIV 2) < n * (f n)   by LESS_EQ_LESS_TRANS
1102*)
1103val sigma_eq_perfect_power_bounds_2 = store_thm(
1104  "sigma_eq_perfect_power_bounds_2",
1105  ``!p. 1 < p ==> !f. (!n. 0 < n ==> (p ** n = SIGMA (\d. d * f d) (divisors n))) ==>
1106   (!n. 0 < n ==> n * (f n) <= p ** n) /\
1107   (!n. 0 < n ==> (p ** (n DIV 2) - 2) * p ** (n DIV 2) < n * (f n))``,
1108  rpt strip_tac >-
1109  rw[sigma_eq_perfect_power_bounds_1] >>
1110  qabbrev_tac `m = n DIV 2` >>
1111  `p ** n - 2 * p ** m < n * (f n)` by rw[sigma_eq_perfect_power_bounds_1, Abbr`m`] >>
1112  `(p ** m - 2) * p ** m <= p ** n - 2 * p ** m` by rw[perfect_power_half_inequality_2, Abbr`m`] >>
1113  decide_tac);
1114
1115(* This is a milestone theorem. *)
1116
1117(* ------------------------------------------------------------------------- *)
1118
1119(* export theory at end *)
1120val _ = export_theory();
1121
1122(*===========================================================================*)
1123