1(* ------------------------------------------------------------------------- *)
2(* AKS Parameter.                                                            *)
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 "computeParam";
12
13(* ------------------------------------------------------------------------- *)
14
15
16(* val _ = load "jcLib"; *)
17open jcLib;
18
19(* val _ = load "SatisfySimps"; (* for SatisfySimps.SATISFY_ss *) *)
20
21(* Get dependent theories local *)
22(* val _ = load "computeOrderTheory"; *)
23open computeOrderTheory;
24open helperFunctionTheory; (* for SQRT_LE *)
25open logPowerTheory; (* for ulog *)
26open ringTheory ringInstancesTheory; (* for ZN_coprime_order_alt *)
27open monoidOrderTheory;
28
29(* Get dependent theories in lib *)
30(* (* val _ = load "helperNumTheory"; -- in monoidTheory *) *)
31(* (* val _ = load "helperSetTheory"; -- in monoidTheory *) *)
32open helperNumTheory helperSetTheory;
33
34(* open dependent theories *)
35open pred_setTheory listTheory arithmeticTheory;
36(* (* val _ = load "dividesTheory"; -- in helperNumTheory *) *)
37(* (* val _ = load "gcdTheory"; -- in helperNumTheory *) *)
38open dividesTheory gcdTheory;
39
40(* (* val _ = load "GaussTheory"; *) *)
41open GaussTheory; (* for phi_pos *)
42(* (* val _ = load "EulerTheory"; *) *)
43open EulerTheory; (* for residue_def *)
44(* (* val _ = load "triangleTheory"; *) *)
45open triangleTheory; (* for list_lcm_pos *)
46
47
48(* ------------------------------------------------------------------------- *)
49(* AKS Parameter Documentation                                               *)
50(* ------------------------------------------------------------------------- *)
51(* Overloading:
52   param_search_result =
53       nice num
54     | good num
55     | bad
56*)
57(* Definitions and Theorems (# are exported):
58
59   Prime Candidates:
60   prime_candidates_def         |- !n m. prime_candidates n m =
61                        {k | prime k /\ ~(k divides n) /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))
62   prime_candidates_element     |- !n m k. k IN prime_candidates n m <=>
63                             prime k /\ ~(k divides n) /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))
64   prime_candidates_property    |- !n m k. k IN prime_candidates n m ==> prime k /\ coprime k n /\ m <= ordz k n
65   prime_candidates_not_empty   |- !n m. 1 < n /\ 0 < m ==> prime_candidates n m <> {}
66   ZN_order_big_enough          |- !n m. 1 < n /\ 0 < m ==> ?k. prime k /\ coprime k n /\ m <= ordz k n
67
68   Coprime Candidates:
69   coprime_candidates_def       |- !n m. coprime_candidates n m =
70                                   {k | coprime k n /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))}
71   coprime_candidates_element   |- !n m k. k IN coprime_candidates n m <=>
72                                        coprime k n /\ !j. 0 < j /\ j < m ==> ~(k divides n ** j - 1)
73   coprime_candidates_property  |- !n m. 1 < n /\ 1 < m ==> !k. k IN coprime_candidates n m ==>
74                                         coprime k n /\ m <= ordz k n
75   coprime_candidates_not_empty |- !n m. 1 < n /\ 0 < m ==> coprime_candidates n m <> {}
76   coprime_candidates_ne_0      |- !n m. 1 < n ==> 0 NOTIN coprime_candidates n m
77   coprime_candidates_ne_1      |- !n m. 1 < m ==> 1 NOTIN coprime_candidates n m
78   ZN_order_good_enough         |- !n m. 1 < n /\ 1 < m ==> ?k. 1 < k /\ coprime k n /\ m <= ordz k n
79
80
81   Smallest Candidate
82   least_prime_candidates_property    |- !n m. 1 < n /\ 0 < m ==> !h. prime h /\ coprime h n /\
83                                               h < MIN_SET (prime_candidates n m) ==>
84                                               ?j. 0 < j /\ j < m /\ h divides (n ** j - 1)
85   least_coprime_candidates_property  |- !n m. 1 < n /\ 1 < m ==> !h. 1 < h /\ coprime h n /\
86                                               h < MIN_SET (coprime_candidates n m) ==>
87                                               ?j. 0 < j /\ j < m /\ h divides (n ** j - 1):
88
89   Power Factors and Product:
90   power_factors_def        |- !n m. power_factors n m = {n ** j - 1 | 0 < j /\ j < m}
91   power_factors_element    |- !n m x. x IN power_factors n m <=> ?j. (x = n ** j - 1) /\ 0 < j /\ j < m
92   power_factors_alt        |- !n m. power_factors n m = IMAGE (\k. n ** k - 1) (residue m)
93   power_factors_finite     |- !n m. FINITE (power_factors n m)
94   power_factors_nonempty   |- !n m. 1 < m ==> power_factors n m <> {}
95   power_factors_sing       |- !m. 1 < m ==> (power_factors 1 m = {0})
96   power_factors_element_nonzero  |- !n m. 1 < n ==> 0 NOTIN power_factors n m
97   power_factors_image_suc  |- !n. 0 < n ==> !m. IMAGE SUC (power_factors n m) = IMAGE (\j. n ** j) (residue m)
98
99   Product of Power Factors:
100   product_factors_def       |- !n m. product_factors n m = PROD_SET (power_factors n m)
101   product_factors_zero      |- !m. 1 < m ==> (product_factors 1 m = 0)
102   proudct_factors_prime_divisor_property   |- !n p. 0 < n /\ prime p ==>
103                                !m. p divides (product_factors n m) ==> ordz p n < m
104   proudct_factors_coprime_divisor_property  |- !n p. 1 < p /\ coprime p n ==>
105                                !m. ordz p n < m ==> p divides (product_factors n m)
106   product_factors_divisors  |- !n m k. 1 < k /\ coprime k n /\ ordz k n < m ==> k divides (product_factors n m)
107   product_factors_pos       |- !n m. 1 < n ==> 0 < product_factors n m
108   product_factors_upper     |- !n m. 1 < n /\ 1 < m ==> product_factors n m < n ** (m ** 2 DIV 2)
109
110   Non-dividing coprimes of Product of Power Factors:
111   power_factors_coprime_nondivisor_def      |- !n m. power_factors_coprime_nondivisor n m =
112                                                 {x | coprime x n /\ ~(x divides (product_factors n m))}
113   power_factors_coprime_nondivisor_element  |- !n m x. x IN power_factors_coprime_nondivisor n m <=>
114                                                        coprime x n /\ ~(x divides (product_factors n m))
115   power_factors_coprime_nondivisor_property |- !n m. 1 < m ==>
116                                                !x. x IN power_factors_coprime_nondivisor n m ==> 1 < x
117   power_factors_coprime_nondivisor_nonempty |- !n m. 1 < n /\ 1 < m ==>
118                                                      power_factors_coprime_nondivisor n m <> {}
119   power_factors_coprime_nondivisor_order    |- !n m. 1 < m ==>
120                                                !k. k IN power_factors_coprime_nondivisor n m ==> m <= ordz k n
121
122   Common Multiples of Residue n:
123   residue_common_multiple_def      |- !n. residue_common_multiple n =
124                                           {m | 0 < m /\ !j. j IN residue n ==> j divides m}
125   residue_common_multiple_element  |- !n m. m IN residue_common_multiple n <=>
126                                             0 < m /\ !j. 1 < j /\ j < n ==> j divides m
127   residue_common_multiple_nonempty |- !n. residue_common_multiple n <> {}
128   residue_common_multiple_nonzero  |- !n. 0 NOTIN residue_common_multiple n
129   residue_common_multiple_has_multiple  |- !n m. m IN residue_common_multiple n ==>
130                                            !k. 0 < k ==> k * m IN residue_common_multiple n
131   residue_common_multiple_has_sub  |- !n a b. b < a /\ a IN residue_common_multiple n /\
132                                               b IN residue_common_multiple n ==>
133                                               a - b IN residue_common_multiple n
134   residue_common_multiple_element_1 |- !n m. 1 < n ==>
135                              (m IN residue_common_multiple n <=> 0 < m /\ !j. 1 < j /\ j < n ==> j divides m)
136
137   Least Common Multiple of Residue n:
138   residue_lcm_def         |- !n. residue_lcm n = MIN_SET (residue_common_multiple n)
139   residue_lcm_divides_common_multiple
140                           |- !n m. m IN residue_common_multiple n ==> (residue_lcm n) divides m
141   residue_lcm_property    |- !n c. (c = residue_lcm n) <=> c IN residue_common_multiple n /\
142                                    !m. m IN residue_common_multiple n ==> c divides m
143   residue_lcm_alt         |- !n. 1 < n ==> (residue_lcm n = list_lcm (leibniz_vertical (n - 2)))
144   residue_lcm_eq_list_lcm |- !n. 1 < n ==> (residue_lcm n = list_lcm [1 .. (n - 1)])
145
146   Bounds for Residue LCM:
147   product_factors_lower        |- !n m k. 1 < n /\ 1 < m /\ 1 < k /\
148                                    (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
149                                    residue_lcm k <= product_factors n m
150   residue_lcm_upper            |- !n m k. 1 < n /\ 1 < m /\ 1 < k /\
151                                    (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
152                                    residue_lcm k < n ** (m ** 2 DIV 2)
153   residue_lcm_lower            |- !n. 1 < n ==> 2 ** (n - 2) <= residue_lcm n
154   ZN_order_modulus_upper       |- !n m k. 1 < n /\ 1 < m /\ 1 < k /\
155                                    (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
156                                    k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2)
157   ZN_order_modulus_exists      |- !n m. 1 < n /\ 1 < m ==> ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\
158                                    (coprime k n ==> m <= ordz k n)
159   ZN_order_modulus_exists_alt  |- !n m. 1 < n /\ 1 < m ==> ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\
160                                    (gcd k n <> 1 \/ m <= ordz k n)
161   ZN_order_modulus_when_prime  |- !n m k. 0 < n /\ prime k /\ m <= ordz k n ==> ~(k divides (product_factors n m))
162
163   ZN bounds based on upper logarithm:
164   ZN_order_modulus_upper_1         |- !n m k. 1 < n /\ 1 < m /\ 1 < k /\
165                                       (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
166                                        k < 2 + ulog n * HALF (m ** 2)
167   ZN_order_modulus_exists_1        |- !n m. 1 < n /\ 1 < m ==>
168                                       ?k. 1 < k /\ k < 2 + ulog n * HALF (m ** 2) /\
169                                           (gcd k n <> 1 \/ m <= ordz k n)
170   ZN_order_modulus_upper_2         |- !n m k. 1 < n /\ 1 < m /\ 1 < k /\
171                                       (!j. 1 < j /\ j < k ==> ~(j divides n) /\ ordz j n < m) ==>
172                                        k < 2 + HALF (m ** 2) * ulog n
173   ZN_order_modulus_exists_2        |- !n m. 1 < n /\ 1 < m ==>
174                                       ?k. 1 < k /\ k < 2 + HALF (m ** 2) * ulog n /\
175                                           (k divides n \/ m <= ordz k n)
176   ZN_order_modulus_upper_3         |- !n m k c. 1 < n /\ 1 < m /\
177                                       (!k. 1 < k /\ k <= c ==> ~(k divides n) /\ ordz k n < m) ==>
178                                        c < 1 + HALF (m ** 2) * ulog n
179   ZN_order_modulus_exists_3        |- !n m c. 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * ulog n <= c ==>
180                                           ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)
181   ZN_order_modulus_exists_4        |- !n m c. 1 < n /\ 1 < m /\ (c = 1 + HALF (m ** 2) * ulog n) ==>
182                                           ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)
183   ZN_order_condition_property      |- !n k. 1 < n /\ coprime k n /\
184                                             TWICE (SUC (LOG2 n)) ** 2 <= ordz k n ==>
185                                             TWICE (SQRT (phi k)) * SUC (LOG2 n) <= phi k
186   ZN_order_condition_property_0    |- !n k. 1 < n /\ coprime k n ==>
187                                       !h. 0 < h /\ h ** 2 <= ordz k n ==>
188                                           SQRT (phi k) * h <= phi k
189   ZN_order_condition_property_1    |- !n k. 1 < n /\ coprime k n /\
190                                             SUC (LOG2 n) ** 2 <= ordz k n ==>
191                                             SQRT (phi k) * SUC (LOG2 n) <= phi k
192   ZN_order_condition_property_2    |- !n k. 1 < n /\ coprime k n /\
193                                             ulog n ** 2 <= ordz k n ==>
194                                             SQRT (phi k) * ulog n <= phi k
195   ZN_order_condition_property_3    |- !n k. 1 < n /\ coprime k n /\
196                                          TWICE (ulog n) ** 2 <= ordz k n ==>
197                                          TWICE (SQRT (phi k)) * ulog n <= phi k
198
199   AKS Parameter Search:
200   aks_param_search_def  |- !n m k c. aks_param_search n m k c =
201                                      if c < k then bad
202                                      else if k divides n then nice k
203                                      else if m <= k /\ m <= ordz_compute k n then good k
204                                      else aks_param_search n m (k + 1) c
205   aks_param_def         |- !n. aks_param n =
206                                    if n <= 2 then nice n
207                                    else aks_param_search n (SQ (ulog n)) 2 (1 + HALF (ulog n ** 5))
208   aks_param_alt         |- !n. aks_param n =
209                                    (let m = ulog n ;
210                                         c = 1 + HALF (m ** 5)
211                                      in if m <= 1 then nice n else aks_param_search n (SQ m) 2 c)
212   aks_param_search_alt  |- !n m k c. aks_param_search n m k c =
213                                      if c < k then bad
214                                      else if k divides n then nice k
215                                      else if m <= k /\ m <= ordz k n then good k
216                                      else aks_param_search n m (k + 1) c
217   aks_param_0           |- aks_param 0 = nice 0
218   aks_param_1           |- aks_param 1 = nice 1
219   aks_param_2           |- aks_param 2 = nice 2
220   aks_param_search_bad      |- !n m k c. 0 < k /\ (aks_param_search n m k c = bad) ==>
221                                !j. k <= j /\ j <= c ==> ~(j divides n) /\ ordz j n < m
222   aks_param_search_not_bad  |- !n m k j c. k <= j /\ j <= c /\ (j divides n \/ m <= ordz j n) ==>
223                                            (k = 0) \/ aks_param_search n m k c <> bad
224   aks_param_search_success  |- !n m c. 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * ulog n <= c ==>
225                                        aks_param_search n m 2 c <> bad
226   aks_param_exists          |- !n. aks_param n <> bad
227
228   aks_param_search_nice     |- !n m h c k. (aks_param_search n m h c = nice k) ==>
229                                            h <= k /\ k divides n /\ !j. h <= j /\ j < k ==> ~(j divides n)
230   aks_param_search_nice_le  |- !n m k c. 0 < n /\ (aks_param_search n m 2 c = nice k) ==> k <= n
231   aks_param_search_nice_bound
232                             |- !n m k c j. (aks_param_search n m k c = nice j) ==> j <= c
233   aks_param_search_nice_factor
234                             |- !n m k c j. (aks_param_search n m k c = nice j) ==> j divides n
235   aks_param_search_nice_coprime_all
236                             |- !n m k c. (aks_param_search n m 2 c = nice k) ==>
237                                !j. 1 < j /\ j < k ==> coprime j n
238   aks_param_search_good     |- !n m h c k. (aks_param_search n m h c = good k) ==>
239                                            h <= k /\ m <= ordz k n /\ !j. h <= j /\ j <= k ==> ~(j divides n)
240   aks_param_search_good_lt  |- !n m k c. 1 < n /\ (aks_param_search n m 2 c = good k) ==> k < n
241   aks_param_search_good_bound
242                             |- !n m k c j. (aks_param_search n m k c = good j) ==> j <= c
243   aks_param_search_good_coprime_all
244                             |- !n m k c. (aks_param_search n m 2 c = good k) ==>
245                                   1 < k /\ m <= ordz k n /\ !j. 1 < j /\ j <= k ==> coprime j n
246   aks_param_nice            |- !n k. 1 < n /\ (aks_param n = nice k) ==>
247                                      1 < k /\ k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n)
248   aks_param_nice_alt        |- !n k. (aks_param n = nice k) ==>
249                                      (if n <= 1 then k = n else 1 < k) /\ k divides n /\
250                                      !j. 1 < j /\ j < k ==> ~(j divides n)
251   aks_param_nice_le         |- !n k. (aks_param n = nice k) ==> k <= n
252   aks_param_nice_bound      |- !n k. 2 < n /\ (aks_param n = nice k) ==> k <= 1 + HALF (ulog n ** 5)
253   aks_param_nice_coprime_all|- !n k. (aks_param n = nice k) ==> !j. 1 < j /\ j < k ==> coprime j n
254   aks_param_nice_for_prime  |- !n k. 1 < n /\ (aks_param n = nice k) ==> (prime n <=> (k = n))
255   aks_param_good            |- !n k. (aks_param n = good k) ==>
256                                         1 < k /\ SQ (ulog n) <= ordz k n /\ !j. 1 < j /\ j <= k ==> ~(j divides n)
257   aks_param_good_lt         |- !n k. (aks_param n = good k) ==> k < n
258   aks_param_good_bound      |- !n k. (aks_param n = good k) ==> k <= 1 + HALF (ulog n ** 5)
259   aks_param_good_coprime_all|- !n k. (aks_param n = good k) ==> !j. 1 < j /\ j <= k ==> coprime j n
260   aks_param_good_range      |- !n k. (aks_param n = good k) ==> 1 < n /\ 1 < k /\ k < n
261   aks_param_good_coprime_order
262                             |- !n k. (aks_param n = good k) ==>
263                                       1 < k /\ k < n /\ coprime k n /\ ulog n ** 2 <= ordz k n
264
265   AKS parameter search simplified:
266   param_seek_def      |- !n m k c. param_seek m c n k =
267                                      if c <= k then bad
268                                      else if n MOD k = 0 then nice k
269                                      else if m <= ordz k n then good k
270                                      else param_seek m c n (k + 1)
271   param_def           |- !n. param n =
272                              if n <= 2 then nice n
273                              else param_seek (SQ (ulog n)) (2 + HALF (ulog n ** 5)) n 2
274   param_alt           |- !n. param n =
275                              (let m = ulog n ;
276                                   c = 2 + HALF (m ** 5)
277                                in if m <= 1 then nice n else param_seek (SQ m) c n 2)
278   param_seek_thm      |- !m c n k. 0 < k /\ k <= c ==>
279                                    (param_seek m c n k = aks_param_search n m k (c - 1))
280   param_eqn           |- !n. param n = aks_param n
281   param_0             |- param 0 = nice 0
282   param_1             |- param 1 = nice 1
283   param_2             |- param 2 = nice 2
284   param_nice_bound    |- !n k. 2 < n /\ (param n = nice k) ==> k <= 1 + HALF (ulog n ** 5)
285   param_good_bound    |- !n k. (param n = good k) ==> k <= 1 + HALF (ulog n ** 5)
286   param_exists        |- !n. param n <> bad
287   param_nice_for_prime|- !n k. 1 < n /\ (param n = nice k) ==> (prime n <=> (k = n))
288   param_good_range    |- !n k. (param n = good k) ==> 1 < n /\ 1 < k /\ k < n
289*)
290
291(* ------------------------------------------------------------------------- *)
292(* Helper Theorems                                                           *)
293(* ------------------------------------------------------------------------- *)
294
295(* ------------------------------------------------------------------------- *)
296(* AKS Parameter                                                             *)
297(* ------------------------------------------------------------------------- *)
298
299(* The AKS parameter k:
300   prime k /\ (2 * SUC (LOG2 n)) ** 2 <= ordz k n
301*)
302
303(* Problem: Given n and m, find the smallest k such that ordz k n > m *)
304
305(* Theory:
306   Note that n is supposed to be large, k is small, the max m is also small.
307     ordz k n
308   = exponent j such that n ** j = 1 (mod k)      by definition of order
309   = exponent j such that k divides (n ** j - 1)  by definition of mod k
310   So try k = 1, 2, 3, ... successively
311      If coprime k n, then n is invertible in (ZN k)  ... otherwise n is composite.
312      Being invertibe, j = ordz k n exists ... j divides phi(k) by Euler-Fermat
313      That is 0 < j < phi(k), j can be found by iteration.
314      If j <= m, reject: next k, else accept: k found.
315*)
316
317(* More Theory:
318   When found, k divides (n ** j - 1)  with j > m = 8, say.
319   Since k is the smallest, this means:
320         k does not divide (n ** 1 - 1), otherwise j = 1, but 1 < 8.
321    also k does not divide (n ** 2 - 1), otherwise j = 2, but 2 < 8.
322    also ...
323    also k does not divide (n ** m - 1), otherwise j = m, but m = 8.
324   Taken all together, k does not divide P = PROD (0 < j <= m) (n ** j - 1)
325   or k is the smallest value not a factor of P.
326   Of course, being the smallest value, k must be prime.
327*)
328
329(* Still More Theory:
330   When found, k divides (n ** j - 1)  with j > m.
331   Say k = 13, the smallest value not a factor of P = PROD (0 < j <= m) (n ** j - 1)
332   Since k is the smallest, this must be due to:
333         1 divides P, otherwise k = 1, but actually k = 13.
334         2 divides P, otherwise k = 2, but actually k = 13.
335         ....
336        12 divides P, otherwise k = 12, but actually k = 13.
337        and we know 13 is the smallest value not dividing P.
338   That means,
339        LCM (1..k) divides P = PROD (0 < j <= m) (n ** j - 1)
340   Since divisor LCM (1..k) increases with k,
341     but dividend P is fixed by n and m,
342   Hence a suitable k will eventually be found if searching starts from 1, 2, ...
343*)
344
345(* Examples:
346
347Let n = 97, m = 8.
348
349n^1 - 1 = (2^5)(3)
350n^2 - 1 = (2^6)(3)(7^2)
351n^3 - 1 = (2^5)(3^2)(3169)
352n^4 - 1 = (2^7)(3)(5)(7^2)(941)
353n^5 - 1 = (2^5)(3)(11)(31)(262321)
354n^6 - 1 = (2^6)(3^2)(7^2)(67)(139)(3169)
355n^7 - 1 = (2^5)(3)(43)(967)(20241187)
356n^8 - 1 = (2^8)(3)(5)(7^2)(233)(941)(189977)
357
358So pick k = 13, the smallest prime not appearing. order_13(97) = order_13(6) = 12 > 8.
359
360Let n = 98, m = 8.
361
362n^1 - 1 = (97)
363n^2 - 1 = (3^2)(11)(97)
364n^3 - 1 = (31)(97)(313)
365n^4 - 1 = (3^2)(5)(11)(17)(97)(113)
366n^5 - 1 = (41)(97)(241)(9431)
367n^6 - 1 = (3^3)(11)(31)(97)(313)(3169)
368n^7 - 1 = (97)(239)(1303)(2873879)
369n^8 - 1 = (3^2)(5)(11)(17)(97)(113)(401)(230017)
370
371The smallest prime not appearing is k = 2, but n = 98 has no order in (ZN 2).
372Need to pick k = 13, smallest coprime not appearing. order_13(98) = order_13(7) = 12 > 8.
373Need to pick k = 13, smallest prime not dividing 98 and not appearing. order_13(98) = order_13(7) = 12 > 8.
374Both are the same due to prime_not_divides_is_coprime: |- !n p. prime p /\ ~(p divides n) ==> coprime p
375*)
376
377(* ------------------------------------------------------------------------- *)
378(* Prime Candidates                                                          *)
379(* ------------------------------------------------------------------------- *)
380
381(* Define the set of candidate primes *)
382val prime_candidates_def = Define `
383  prime_candidates (n:num) (m:num) =
384    {k | prime k /\ ~(k divides n) /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))}
385`;
386
387(* Theorem: k IN (prime_candidates n m) <=>
388           prime k /\ ~(k divides n) /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1)) *)
389(* Proof: by prime_candidates_def *)
390val prime_candidates_element = store_thm(
391  "prime_candidates_element",
392  ``!n m k. k IN (prime_candidates n m) <=>
393           prime k /\ ~(k divides n) /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))``,
394  rw[prime_candidates_def]);
395
396(* Theorem: k IN prime_candidates n m ==> prime k /\ coprime k n /\ m <= ordz k n *)
397(* Proof:
398   Expand by prime_candidates_def, this is to show:
399   (1) prime k /\ ~(k divides n) ==> coprime k n
400       True by prime_not_divides_is_coprime.
401   (2) prime k /\ ~(k divides n) /\
402       !j. 0 < j /\ j <= m ==> ~(k divides (n ** j - 1)) ==> m <= ordz k n
403   Let h = ordz k n. To show: m <= h.
404   Note prime k /\ ~(k divides n) ==> coprime k n   by prime_not_divides_is_coprime
405   By contradiction. Suppose h < m.
406   Then prime k ==> 1 < k                           by ONE_LT_PRIME
407    and 1 < k /\ coprime k n
408    ==> 0 < h /\ (n ** h MOD k = 1)                 by ZN_coprime_order_alt, 1 < k, coprime k n
409   Thus 0 < h /\ h < m /\ (n ** h - 1) MOD k = 0    by MOD_EQ_DIFF, 0 < k
410   or   0 < h /\ h < m /\ k divides (n ** h - 1)    by DIVIDES_MOD_0, 0 < k
411   This is a contradiction with the implication when j = h.
412*)
413val prime_candidates_property = store_thm(
414  "prime_candidates_property",
415  ``!n m k. k IN prime_candidates n m ==> prime k /\ coprime k n /\ m <= ordz k n``,
416  rw[prime_candidates_def] >-
417  rw[prime_not_divides_is_coprime] >>
418  qabbrev_tac `h = ordz k n` >>
419  spose_not_then strip_assume_tac >>
420  `coprime k n` by rw[prime_not_divides_is_coprime] >>
421  `1 < k` by rw[ONE_LT_PRIME] >>
422  `0 < h /\ (n ** h MOD k = 1)` by rw[ZN_coprime_order_alt, Abbr`h`] >>
423  `0 < k /\ h < m` by decide_tac >>
424  metis_tac[MOD_EQ_DIFF, DIVIDES_MOD_0, ONE_MOD]);
425
426(* Theorem: 1 < n /\ 0 < m ==> prime_candidates n m <> {} *)
427(* Proof:
428   Expand by prime_candidates_def, this is to show:
429   ?x. prime x /\ ~(x divides n) /\ !j. j <= m /\ x divides (n ** j - 1) ==> (j = 0)
430   Let z = n ** m.
431       ?p. prime p /\ z < p           by prime_always_bigger
432   Take x = p. Need to show:
433   (1) prime p /\ z < p ==> ~(p divides n)
434       Since n <= z                   by X_LE_X_EXP, 0 < m
435       n < p, or ~(p <= n)
436       Hence ~(p divides n)           by DIVIDES_LE
437   (2) p divides (n ** j - 1) ==> j = 0
438       Note 0 < n ** j                by EXP_POS, 0 < n
439            n ** j < z                by EXP_BASE_LT_MONO, 1 < n
440         so n ** j - 1 < z
441        and n ** j - 1 < p
442       which means: ~(p <= n ** j - 1)
443      Given p divides (n ** j - 1),
444      this forces ~(0 < n ** j - 1)   by DIVIDES_LE
445         or n ** j = 1
446      Therefore j = 0                 by EXP_EQ_1, n <> 1
447*)
448val prime_candidates_not_empty = store_thm(
449  "prime_candidates_not_empty",
450  ``!n m. 1 < n /\ 0 < m ==> prime_candidates n m <> {}``,
451  rw[prime_candidates_def, EXTENSION] >>
452  `?x. prime x /\ ~(x divides n) /\ !j. j < m /\ x divides (n ** j - 1) ==> (j = 0)` suffices_by metis_tac[] >>
453  qabbrev_tac `z = n ** m` >>
454  `?p. prime p /\ z < p` by rw[prime_always_bigger] >>
455  qexists_tac `p` >>
456  `0 < n /\ n <> 1` by decide_tac >>
457  rw[] >| [
458    `n <= z` by rw[X_LE_X_EXP, Abbr`z`] >>
459    `~(p <= n)` by decide_tac >>
460    metis_tac[DIVIDES_LE],
461    `0 < n ** j` by rw[EXP_POS] >>
462    `n ** j < z` by rw[EXP_BASE_LT_MONO, Abbr`z`] >>
463    `~(p <= n ** j - 1)` by decide_tac >>
464    `~(0 < n ** j - 1)` by metis_tac[DIVIDES_LE] >>
465    `n ** j = 1` by decide_tac >>
466    metis_tac[EXP_EQ_1]
467  ]);
468
469(* Theorem: 1 < n /\ 0 < m ==> ?k. prime k /\ coprime k n /\ m <= ordz k n *)
470(* Proof:
471   Let s = prime_candidates n m
472   With 1 < n /\ 0 < m, s <> {}     by prime_candidates_not_empty
473     so MIN_SET s IN s              by MIN_SET_LEM
474   Hence true with k = MIN_SET s    by prime_candidates_property
475*)
476val ZN_order_big_enough = store_thm(
477  "ZN_order_big_enough",
478  ``!n m. 1 < n /\ 0 < m ==> ?k. prime k /\ coprime k n /\ m <= ordz k n``,
479  rpt strip_tac >>
480  qabbrev_tac `s = prime_candidates n m` >>
481  `s <> {}` by rw[prime_candidates_not_empty, Abbr`s`] >>
482  `MIN_SET s IN s` by rw[MIN_SET_LEM] >>
483  qexists_tac `MIN_SET s` >>
484  metis_tac[prime_candidates_property]);
485
486(* ------------------------------------------------------------------------- *)
487(* Coprime Candidates                                                        *)
488(* ------------------------------------------------------------------------- *)
489
490(* Define the set of candidate coprimes *)
491val coprime_candidates_def = Define `
492  coprime_candidates (n:num) (m:num) =
493    {k | coprime k n /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))}
494`;
495
496(* Theorem: k IN (coprime_candidates n m) <=> coprime k n /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1)) *)
497(* Proof: by coprime_candidates_def *)
498val coprime_candidates_element = store_thm(
499  "coprime_candidates_element",
500  ``!n m k. k IN (coprime_candidates n m) <=> coprime k n /\ !j. 0 < j /\ j < m ==> ~(k divides (n ** j - 1))``,
501  rw[coprime_candidates_def]);
502
503(* Theorem: 1 < n /\ 1 < m ==>
504   !k. k IN coprime_candidates n m ==> coprime k n /\ m <= ordz k n *)
505(* Proof:
506   Let h = ordz k n. By coprime_candidates_def, this is to show m <= h.
507   By contradiction, suppose h < m.
508   Note k <> 0                         by coprime_0L, n <> 1.
509   Since 0 < 1 /\ 1 < m, ~k divides (n ** m - 1)    by implication
510      so k <> 1                        by ONE_DIVIDES_ALL
511   Hence 0 < k /\ 1 < k
512   Now, 1 < k /\ coprime k n
513   gives 0 < h /\ (n ** h MOD k = 1)   by ZN_coprime_order_alt
514      or (n ** h - 1) MOD k = 0        by MOD_EQ_DIFF, k <> 0.
515   meaning k divides (n ** h - 1)      by DIVIDES_MOD_0
516   which contradicts the implication when j = h.
517*)
518val coprime_candidates_property = store_thm(
519  "coprime_candidates_property",
520  ``!n m. 1 < n /\ 1 < m ==>
521   !k. k IN coprime_candidates n m ==> coprime k n /\ m <= ordz k n``,
522  rw[coprime_candidates_def] >>
523  qabbrev_tac `h = ordz k n` >>
524  spose_not_then strip_assume_tac >>
525  `n <> 1 /\ h < m /\ 0 < 1` by decide_tac >>
526  `k <> 0` by metis_tac[coprime_0L] >>
527  `k <> 1` by metis_tac[ONE_DIVIDES_ALL] >>
528  `0 < k /\ 1 < k` by decide_tac >>
529  `0 < h /\ (n ** h MOD k = 1)` by rw[ZN_coprime_order_alt, Abbr`h`] >>
530  `(n ** h - 1) MOD k = 0` by rw[MOD_EQ_DIFF] >>
531  metis_tac[DIVIDES_MOD_0]);
532
533(* Theorem: 1 < n /\ 1 < m ==> coprime_candidates n m <> {} *)
534(* Proof:
535   By coprime_candidates_def, this is to show:
536      ?x. coprime x n /\ !j. j < m /\ x divides (n ** j - 1) ==> (j = 0)
537   Let z = n ** m.
538       ?p. prime p /\ z < p           by prime_always_bigger
539   Take x = p. Need to show:
540   (1) prime p /\ z < p ==> coprime p n
541       Since n <= z                   by X_LE_X_EXP, 0 < m
542          so n < p
543       Hence coprime p n              by prime_iff_coprime_all_lt
544   (2) p divides (n ** j - 1) ==> j = 0
545       Note 0 < n ** j                by EXP_POS, 0 < n
546            n ** j < z                by EXP_BASE_LT_MONO, 1 < n
547         so n ** j - 1 < z
548        and n ** j - 1 < p
549       which means: ~(p <= n ** j - 1)
550      Given p divides (n ** j - 1),
551      this forces ~(0 < n ** j - 1)   by DIVIDES_LE
552         or n ** j = 1
553      Therefore j = 0                 by EXP_EQ_1, n <> 1
554*)
555val coprime_candidates_not_empty = store_thm(
556  "coprime_candidates_not_empty",
557  ``!n m. 1 < n /\ 0 < m ==> coprime_candidates n m <> {}``,
558  rw[coprime_candidates_def, EXTENSION] >>
559  `?x. coprime x n /\ !j. j < m /\ x divides (n ** j - 1) ==> (j = 0)` suffices_by metis_tac[] >>
560  qabbrev_tac `z = n ** m` >>
561  `?p. prime p /\ z < p` by rw[prime_always_bigger] >>
562  qexists_tac `p` >>
563  `0 < n /\ n <> 1` by decide_tac >>
564  rw[] >| [
565    `n <= z` by rw[X_LE_X_EXP, Abbr`z`] >>
566    `n < p` by decide_tac >>
567    metis_tac[prime_iff_coprime_all_lt],
568    `0 < n ** j` by rw[EXP_POS] >>
569    `n ** j < z` by rw[EXP_BASE_LT_MONO, Abbr`z`] >>
570    `~(p <= n ** j - 1)` by decide_tac >>
571    `~(0 < n ** j - 1)` by metis_tac[DIVIDES_LE] >>
572    `n ** j = 1` by decide_tac >>
573    metis_tac[EXP_EQ_1]
574  ]);
575
576(* Theorem: 1 < n /\ 1 < m ==> 0 NOTIN coprime_candidates n m *)
577(* Proof:
578   By contradiction, suppose 0 IN coprime_candidates n m.
579   Then coprime 0 n     by coprime_candidates_element
580   But  gcd 0 n = n     by GCD_0L
581   and  1 < n means n <> 1, so this is a contradiction.
582*)
583val coprime_candidates_ne_0 = store_thm(
584  "coprime_candidates_ne_0",
585  ``!n m. 1 < n ==> 0 NOTIN coprime_candidates n m``,
586  spose_not_then strip_assume_tac >>
587  `n <> 1` by decide_tac >>
588  metis_tac[coprime_candidates_element, GCD_0L]);
589
590(* Theorem: 1 < m ==> 1 NOTIN coprime_candidates n m *)
591(* Proof:
592   By contradiction, suppose 1 IN coprime_candidates n m.
593   Then coprime 1 n is true   by GCD_1
594   But Take j = 1, as 0 < 1 and 1 < m.
595   and 1 divides (n ** j - 1)   by ONE_DIVIDES_ALL
596   This is a contradiction.
597*)
598val coprime_candidates_ne_1 = store_thm(
599  "coprime_candidates_ne_1",
600  ``!n m. 1 < m ==> 1 NOTIN coprime_candidates n m``,
601  metis_tac[coprime_candidates_element, ONE_DIVIDES_ALL, DECIDE``0 < 1``]);
602
603(* Theorem: 1 < n /\ 1 < m ==> ?k. 1 < k /\ coprime k n /\ m <= ordz k n *)
604(* Proof:
605   Note 1 < m ==> 0 < m             by arithmetic
606   Let s = coprime_candidates n m
607   With 1 < n /\ 0 < m, s <> {}     by coprime_candidates_not_empty
608     so MIN_SET s IN s              by MIN_SET_LEM
609    and MIN_SET s <> 0              by coprime_candidates_ne_0, 1 < n
610    and MIN_SET s <> 1              by coprime_candidates_ne_1, 1 < m
611     so 1 < MIN_SET s               by arithmetic
612   Hence true with k = MIN_SET s    by coprime_candidates_property, 1 < n, 1 < m
613*)
614val ZN_order_good_enough = store_thm(
615  "ZN_order_good_enough",
616  ``!n m. 1 < n /\ 1 < m ==> ?k. 1 < k /\ coprime k n /\ m <= ordz k n``,
617  rpt strip_tac >>
618  qabbrev_tac `s = coprime_candidates n m` >>
619  `s <> {}` by rw[coprime_candidates_not_empty, Abbr`s`] >>
620  `MIN_SET s IN s` by rw[MIN_SET_LEM] >>
621  `MIN_SET s <> 0` by metis_tac[coprime_candidates_ne_0] >>
622  `MIN_SET s <> 1` by metis_tac[coprime_candidates_ne_1] >>
623  `1 < MIN_SET s` by decide_tac >>
624  metis_tac[coprime_candidates_property]);
625
626(* ------------------------------------------------------------------------- *)
627(* Smallest Candidate                                                        *)
628(* ------------------------------------------------------------------------- *)
629
630(* Theorem: prime h /\ coprime h n /\ h < MIN_SET (prime_candidates n m) ==>
631            ?j. 0 < j /\ j < m /\ h divides (n ** j - 1) *)
632(* Proof:
633   By contradiction, assume !j. 0 < j /\ j < m ==> ~(h divides (n ** j - 1)).
634   Let k = MIN_SET (prime_candidates n m).
635   Since MIN_SET (prime_candidates n m) <> {}         by prime_candidates_not_empty
636         k IN (prime_candidates n m) and
637         !x. x IN (prime_candidates n m) ==> k <= x   by MIN_SET_LEM
638   But 1 < h                                          by ONE_LT_PRIME
639   and coprime h n ==> ~(h divides n)                 by coprime_not_divides
640   Hence h IN (prime_candidates n m)                  by prime_candidates_def
641   The implication gives k <= h, which contradicts h < k.
642*)
643val least_prime_candidates_property = store_thm(
644  "least_prime_candidates_property",
645  ``!n m.  1 < n /\ 0 < m ==>
646   !h. prime h /\ coprime h n /\ h < MIN_SET (prime_candidates n m) ==>
647    ?j. 0 < j /\ j < m /\ h divides (n ** j - 1)``,
648  rpt strip_tac >>
649  qabbrev_tac `k = MIN_SET (prime_candidates n m)` >>
650  `k IN (prime_candidates n m) /\ !x. x IN (prime_candidates n m) ==> k <= x` by rw[MIN_SET_LEM, prime_candidates_not_empty, Abbr`k`] >>
651  `prime k /\ coprime k n /\ m <= ordz k n` by metis_tac[prime_candidates_property] >>
652  spose_not_then strip_assume_tac >>
653  `1 < h` by rw[ONE_LT_PRIME] >>
654  `~(h divides n)` by rw[coprime_not_divides] >>
655  `h IN (prime_candidates n m)` by rw[prime_candidates_def] >>
656  `k <= h` by metis_tac[] >>
657  decide_tac);
658
659(* Theorem: 1 < h /\ coprime h n /\ h < MIN_SET (coprime_candidates n m) ==>
660            ?j. 0 < j /\ j < m /\ h divides (n ** j - 1) *)
661(* Proof:
662   By contradiction, assume !j. 0 < j /\ j < m ==> ~(h divides (n ** j - 1)).
663   Let k = MIN_SET (coprime_candidates n m).
664   Since MIN_SET (coprime_candidates n m) <> {}         by coprime_candidates_not_empty
665         k IN (coprime_candidates n m) and
666         !x. x IN (coprime_candidates n m) ==> k <= x   by MIN_SET_LEM
667   But 1 < h                                            by ONE_LT_PRIME
668   and coprime h n ==> ~(h divides n)                   by coprime_not_divides
669   Hence h IN (coprime_candidates n m)                  by coprime_candidates_def
670   The implication gives k <= h, which contradicts h < k.
671*)
672val least_coprime_candidates_property = store_thm(
673  "least_coprime_candidates_property",
674  ``!n m.  1 < n /\ 1 < m ==>
675   !h. 1 < h /\ coprime h n /\ h < MIN_SET (coprime_candidates n m) ==>
676    ?j. 0 < j /\ j < m /\ h divides (n ** j - 1)``,
677  rpt strip_tac >>
678  qabbrev_tac `k = MIN_SET (coprime_candidates n m)` >>
679  `0 < m` by decide_tac >>
680  `k IN (coprime_candidates n m) /\ !x. x IN (coprime_candidates n m) ==> k <= x` by rw[MIN_SET_LEM, coprime_candidates_not_empty, Abbr`k`] >>
681  `coprime k n /\ m <= ordz k n` by metis_tac[coprime_candidates_property] >>
682  spose_not_then strip_assume_tac >>
683  `~(h divides n)` by rw[coprime_not_divides] >>
684  `h IN (coprime_candidates n m)` by rw[coprime_candidates_def] >>
685  `k <= h` by metis_tac[] >>
686  decide_tac);
687
688(* ------------------------------------------------------------------------- *)
689(* Power Factors and Product                                                 *)
690(* ------------------------------------------------------------------------- *)
691
692(* Define set of factors *)
693val power_factors_def = Define`
694  power_factors n m = {n ** j - 1 | j | 0 < j /\ j < m}
695`;
696
697(* Theorem: x IN power_factors n m <=> ?j. (x = n ** j - 1) /\ 0 < j /\ j < m *)
698(* Proof: by power_factors_def. *)
699val power_factors_element = store_thm(
700  "power_factors_element",
701  ``!n m x. x IN power_factors n m <=> ?j. (x = n ** j - 1) /\ 0 < j /\ j < m``,
702  rw[power_factors_def]);
703
704(* Theorem: power_factors n m = IMAGE (\k. n ** k - 1) (residue m) *)
705(* Proof:
706   By power_factors_def and residue_def, this is to show:
707   {n ** j - 1 | 0 < j /\ j < m} = IMAGE (\k. n ** k - 1) {i | 0 < i /\ i < m}
708   which is true by IMAGE_DEF.
709*)
710val power_factors_alt = store_thm(
711  "power_factors_alt",
712  ``!n m. power_factors n m = IMAGE (\k. n ** k - 1) (residue m)``,
713  rw[power_factors_def, residue_def, IMAGE_DEF]);
714
715(* Theorem: FINITE (power_factors n m) *)
716(* Proof:
717   Since power_factors n m
718       = IMAGE (\k. n ** k - 1) (residue m)   by power_factors_alt
719     and FINITE (residue m)                   by residue_finite
720   Hence FINITE (power_factors n m)           by IMAGE_FINITE
721*)
722val power_factors_finite = store_thm(
723  "power_factors_finite",
724  ``!n m. FINITE (power_factors n m)``,
725  rw[power_factors_alt, residue_finite, IMAGE_FINITE]);
726
727(* Theorem: 1 < m ==> (power_factors n m) <> {} *)
728(* Proof:
729   By power_factors_alt, this is to show: 1 < m ==> residue m <> {}
730   This is true by residue_nonempty.
731*)
732val power_factors_nonempty = store_thm(
733  "power_factors_nonempty",
734  ``!n m. 1 < m ==> (power_factors n m) <> {}``,
735  rw[power_factors_alt] >>
736  rw[residue_nonempty]);
737
738(* Theorem: 1 < m ==> (power_factors 1 m = {0}) *)
739(* Proof:
740   By power_factors_def and residue_def, this is to show: 1 < m ==> ?j. 0 < j /\ j < m
741   Take j = 1, this is true.
742*)
743val power_factors_sing = store_thm(
744  "power_factors_sing",
745  ``!m. 1 < m ==> (power_factors 1 m = {0})``,
746  rw[power_factors_def, residue_def, EXTENSION, SPECIFICATION, EQ_IMP_THM] >>
747  metis_tac[DECIDE ``0 < 1``]);
748
749(* Theorem: 1 < n ==> 0 NOTIN power_factors n m *)
750(* Proof:
751   By power_factors_def, this is to show: ~(n ** j <= 1) \/ (j = 0) \/ ~(j < m)
752   By contradiction, this is to show: 1 < n /\ j <> 0 /\ j < m /\ n * j <= 1 ==> F
753   Note 0 < m             by j <> /\ j < m
754   Since 1 < n ** j       by ONE_LT_EXP, 1 < n, 0 < m
755   This contradicts n * j <= 1.
756*)
757val power_factors_element_nonzero = store_thm(
758  "power_factors_element_nonzero",
759  ``!n m. 1 < n ==> 0 NOTIN power_factors n m``,
760  rw[power_factors_def] >>
761  spose_not_then strip_assume_tac >>
762  `1 < n ** j` by rw[ONE_LT_EXP] >>
763  decide_tac);
764
765(* Theorem: 0 < n ==> !m. IMAGE SUC (power_factors n m) = IMAGE (\j. n ** j) (residue m) *)
766(* Proof:
767   Expanding by definitions, this is to show:
768   (1) k <> 0 /\ k < m ==> ?j. (SUC (n ** k - 1) = n ** j) /\ j < m /\ j <> 0
769       Take j = k, to show: SUC (n ** k - 1) = n ** k
770       Since 0 < n ** k      by EXP_POS, 0 < n
771         SUC (n ** k - 1)
772       = SUC (PRE (n ** k))  by PRE_SUB1
773       = n ** k              by SUC_PRE, 0 < n ** k
774   (2) j < m /\ j <> 0 ==> ?x'. (n ** j = SUC x') /\ ?k. (x' = n ** k - 1) /\ k < m /\ k <> 0
775       Take x' = PRE (n ** j)
776       Since 0 < n ** j                   by EXP_POS, 0 < n
777          so n ** j = SUC (PRE (n ** j))  by SUC_PRE
778       Take k = j, to show: PRE (n ** j) = n ** j - 1
779          which is true                   by PRE_SUB1
780*)
781val power_factors_image_suc = store_thm(
782  "power_factors_image_suc",
783  ``!n. 0 < n ==> !m. IMAGE SUC (power_factors n m) = IMAGE (\j. n ** j) (residue m)``,
784  rw[power_factors_alt, residue_def, IMAGE_DEF, EXTENSION, EQ_IMP_THM] >| [
785    qexists_tac `k` >>
786    rw[] >>
787    `0 < n ** k` by rw[] >>
788    metis_tac[PRE_SUB1, PRE_SUC_EQ],
789    qexists_tac `PRE (n ** j)` >>
790    `0 < n ** j` by rw[] >>
791    rw[GSYM SUC_PRE] >>
792    qexists_tac `j` >>
793    rw[PRE_SUB1]
794  ]);
795
796(* ------------------------------------------------------------------------- *)
797(* Product of Power Factors                                                  *)
798(* ------------------------------------------------------------------------- *)
799
800(* Define product of factors *)
801val product_factors_def = Define`
802  product_factors n m = PROD_SET (power_factors n m)
803`;
804
805(* Theorem: 1 < m ==> (product_factors 1 m = 0) *)
806(* Proof:
807   Since power_factors 1 m = {0}
808     product_factors 1 m
809   = PROD_SET (power_factors 1 m)  by product_factors_def
810   = PROD_SET {0}                  by power_factors_sing
811   = 0                             by PROD_SET_SING
812*)
813val product_factors_zero = store_thm(
814  "product_factors_zero",
815  ``!m. 1 < m ==> (product_factors 1 m = 0)``,
816  rw[product_factors_def, power_factors_sing, PROD_SET_SING]);
817
818(* Theorem: 0 < n /\ prime p ==> !m. p divides (product_factors n m) ==> (ordz p n) < m *)
819(* Proof:
820   Let s = power_factors n m
821    Then FINITE (power_factors n m)               by power_factors_finite
822      so ?b. b IN s /\ p divides b                by PROD_SET_EUCLID
823      or ?j. (b = n ** j - 1) /\ 0 < j /\ j < m   by power_factors_element
824    Now  p divides b means b MOD p = 0            by DIVIDES_MOD_0, 0 < p
825     or  (n ** j - 1) MOD p = 0                   by b = n ** j - 1
826     so  n ** j MOD p = 1                         by MOD_EQ, ONE_MOD
827    Since  Ring (ZN p)                            by ZN_ring, 0 < p
828       ==> Monoid (ZN p).prod                     by ring_mult_monoid
829      Let  x = n MOD p, x IN (ZN p).carrier       by MOD_LESS, ZN_property
830      and  ordz p n = ordz p x                    by ZN_order_mod, 0 < p
831    Hence  (ordz p x) j divides                   by monoid_order_condition
832       or  (ordz p n) <= j                        by DIVIDES_LE
833       or  (ordz p n) < m                         by j < m
834*)
835val proudct_factors_prime_divisor_property = store_thm(
836  "proudct_factors_prime_divisor_property",
837  ``!n p. 0 < n /\ prime p ==> !m. p divides (product_factors n m) ==> (ordz p n) < m``,
838  rw[product_factors_def] >>
839  qabbrev_tac `s = power_factors n m` >>
840  `FINITE s` by rw[power_factors_finite, Abbr`s`] >>
841  `?b. b IN s /\ p divides b` by rw[PROD_SET_EUCLID] >>
842  `?j. (b = n ** j - 1) /\ 0 < j /\ j < m` by metis_tac[power_factors_element] >>
843  `1 < p` by rw[ONE_LT_PRIME] >>
844  `0 < p` by decide_tac >>
845  `(n ** j - 1) MOD p = 0` by rw[GSYM DIVIDES_MOD_0] >>
846  `0 < n ** j` by rw[EXP_POS] >>
847  `1 <= n ** j` by decide_tac >>
848  `n ** j MOD p = 1` by metis_tac[MOD_EQ, ONE_MOD] >>
849  qabbrev_tac `x = n MOD p` >>
850  `x < p` by rw[MOD_LESS, Abbr`x`] >>
851  `x ** j MOD p = 1` by rw[Abbr`x`] >>
852  `ordz p n = ordz p x` by rw_tac std_ss[ZN_order_mod, Abbr`x`] >>
853  `Ring (ZN p)` by rw[ZN_ring] >>
854  `Monoid (ZN p).prod /\ x IN (ZN p).prod.carrier` by rw[ZN_property, ring_mult_monoid] >>
855  `p <> 1` by decide_tac >>
856  `(ZN p).prod.id = 1` by rw[ZN_property] >>
857  `(ordz p x) divides j` by metis_tac[monoid_order_condition, ZN_exp] >>
858  `ordz p n <= j` by metis_tac[DIVIDES_LE] >>
859  decide_tac);
860
861(* Theorem: 1 < p /\ coprime p n ==> !m. (ordz p n) < m ==> p divides (product_factors n m) *)
862(* Proof:
863   Let s = power_factors n m
864    Then FINITE s                        by power_factors_finite
865   Let j = ordz p n
866   Then 0 < j /\ (n ** j MOD p = 1)      by ZN_coprime_order_alt, 1 < p, coprime p n
867     or (n ** j - 1) MOD p = 0           by MOD_EQ_DIFF, ONE_MOD, need 0 < p
868    Let b = n ** j - 1, then b IN s      by power_factors_element, need j < m
869    and p divides b                      by DIVIDES_MOD_0
870   Hence  PROD_SET s
871        = PROD_SET (b INSERT s)          by ABSORPTION
872        = b * PROD_SET (s DELETE b)      by PROD_SET_THM, need FINITE s
873        = PROD_SET (s DELETE b) * b      by MULT_COMM
874   So   b divides (PROD_SET s)           by divides_def
875   or   p divides (PROD_SET s)           by DIVIDES_TRANS
876*)
877val proudct_factors_coprime_divisor_property = store_thm(
878  "proudct_factors_coprime_divisor_property",
879  ``!n p. 1 < p /\ coprime p n ==> !m. (ordz p n) < m ==> p divides (product_factors n m)``,
880  rw[product_factors_def] >>
881  qabbrev_tac `s = power_factors n m` >>
882  `FINITE s` by rw[power_factors_finite, Abbr`s`] >>
883  qabbrev_tac `j = ordz p n` >>
884  `0 < j /\ (n ** j MOD p = 1)` by rw[ZN_coprime_order_alt, Abbr`j`] >>
885  `0 < p` by decide_tac >>
886  `(n ** j - 1) MOD p = 0` by rw[MOD_EQ_DIFF, ONE_MOD] >>
887  qabbrev_tac `b = n ** j - 1` >>
888  `p divides b` by rw[DIVIDES_MOD_0, Abbr`b`] >>
889  `b IN s` by metis_tac[power_factors_element] >>
890  `b INSERT s = s` by rw[GSYM ABSORPTION] >>
891  `PROD_SET s = b * PROD_SET (s DELETE b)` by rw[GSYM PROD_SET_THM] >>
892  `b divides (PROD_SET s)` by metis_tac[divides_def, MULT_COMM] >>
893  metis_tac[DIVIDES_TRANS]);
894
895(* Theorem: 1 < k /\ coprime k n /\ ordz k n < m ==> k divides (product_factors n m) *)
896(* Proof:
897   By product_factors_def, this is to show:
898      ordz k n < m ==> k divides (PROD_SET (power_factors n m))
899   Let j = ordz k n
900   Then 0 < j /\ (n ** j MOD k = 1)           by ZN_coprime_order_alt, 1 < k
901     or (n ** j - 1) MOD k = 0                by MOD_EQ_DIFF, MOD_LESS, 0 < k
902     or k divides (n ** j - 1)                by DIVIDES_MOD_0, 0 < k
903   Since FINITE (power_factors n m            by power_factors_finite
904     and (n ** j - 1) IN (power_factors n m)  by power_factors_element
905   Hence  k divides (product_factors n m)     by PROD_SET_DIVISORS
906*)
907val product_factors_divisors = store_thm(
908  "product_factors_divisors",
909  ``!n m k. 1 < k /\ coprime k n /\ ordz k n < m ==> k divides (product_factors n m)``,
910  rw[product_factors_def] >>
911  qabbrev_tac `j = ordz k n` >>
912  `0 < j /\ (n ** j MOD k = 1)` by rw[ZN_coprime_order_alt, Abbr`j`] >>
913  `(n ** j - 1) MOD k = 0` by rw[MOD_EQ_DIFF, MOD_LESS] >>
914  `k divides (n ** j - 1)` by rw[DIVIDES_MOD_0] >>
915  metis_tac[PROD_SET_DIVISORS, power_factors_element, power_factors_finite]);
916
917(* Theorem: 1 < n ==> 0 < product_factors n m *)
918(* Proof:
919   By product_factors_def, this is to show: 0 < PROD_SET (power_factors n m)
920   Since FINITE (power_factors n m)             by power_factors_finite
921     and 0 NOTIN power_factors n m              by power_factors_element_nonzero
922   Hence 0 < PROD_SET (power_factors n m)       by PROD_SET_NONZERO
923*)
924val product_factors_pos = store_thm(
925  "product_factors_pos",
926  ``!n m. 1 < n ==> 0 < product_factors n m``,
927  rw[product_factors_def] >>
928  `FINITE (power_factors n m)` by rw[power_factors_finite] >>
929  `0 NOTIN power_factors n m` by metis_tac[power_factors_element_nonzero] >>
930  rw[PROD_SET_NONZERO]);
931
932(* ------------------------------------------------------------------------- *)
933(* Upper Bound of Product of Power Factors                                   *)
934(* ------------------------------------------------------------------------- *)
935
936(* Theorem: 1 < n /\ 1 < m ==> product_factors n m < n ** ((m ** 2) DIV 2) *)
937(* Proof:
938   By product_factors_def, this is to show: PROD_SET (power_factors n m) < n ** (m ** 2 DIV 2)
939   Let s = power_factors n m.
940    Then FINITE s                               by power_factors_finite
941     and s <> {}                                by power_factors_nonempty, 1 < m
942    also 0 NOTIN s                              by power_factors_element_nonzero, 1 < n
943   Hence PROD_SET s < PROD_SET (IMAGE SUC s)    by PROD_SET_LESS_SUC
944     PROD_SET (IMAGE SUC s)
945   = PROD_SET (IMAGE (\j. n ** j) (residue m))  by power_factors_image_suc, 0 < n
946   = PROD_SET (IMAGE (\j. n ** j) (count m))    by PROD_SET_IMAGE_EXP_NONZERO, 1 < n /\ 0 < m
947   = n ** SUM_SET (count m)                     by PROD_SET_IMAGE_EXP, 1 < n
948   = n ** (m * (m - 1) DIV 2)                   by SUM_SET_COUNT
949
950   Since m * (m - 1) < m * m = m ** 2           by EXP_2
951      or m * (m - 1) <= m ** 2                  by LESS_IMP_LESS_OR_EQ
952         (m * (m - 1)) DIV 2 <= (m ** 2) DIV 2  by DIV_LE_MONOTONE
953      so n ** ((m * (m - 1)) DIV 2) <= n ** ((m ** 2) DIV 2)   by EXP_BASE_LE_MONO, 1 < n
954   Therefore,
955      PROD_SET s < PROD_SET (IMAGE SUC s)
956                 = n ** (m * (m - 1) DIV 2)
957                 <= n ** ((m ** 2) DIV 2)
958   or PROD_SET s < n ** ((m ** 2) DIV 2)        by LESS_LESS_EQ_TRANS
959*)
960val product_factors_upper = store_thm(
961  "product_factors_upper",
962  ``!n m. 1 < n /\ 1 < m ==> product_factors n m < n ** ((m ** 2) DIV 2)``,
963  rw[product_factors_def] >>
964  qabbrev_tac `s = power_factors n m` >>
965  `FINITE s` by rw[power_factors_finite, Abbr`s`] >>
966  `s <> {}` by rw[power_factors_nonempty, Abbr`s`] >>
967  `0 NOTIN s` by metis_tac[power_factors_element_nonzero] >>
968  `PROD_SET s < PROD_SET (IMAGE SUC s)` by rw[PROD_SET_LESS_SUC] >>
969  `0 < n /\ 0 < m` by decide_tac >>
970  `PROD_SET (IMAGE SUC s) = PROD_SET (IMAGE (\j. n ** j) (residue m))` by rw[GSYM power_factors_image_suc] >>
971  `_ = PROD_SET (IMAGE (\j. n ** j) (count m))` by rw[PROD_SET_IMAGE_EXP_NONZERO] >>
972  `_ = n ** SUM_SET (count m)` by rw[PROD_SET_IMAGE_EXP] >>
973  `_ = n ** (m * (m - 1) DIV 2)` by rw[SUM_SET_COUNT] >>
974  `m * (m - 1) < m ** 2` by rw_tac std_ss[EXP_2] >>
975  `(m * (m - 1)) DIV 2 <= (m ** 2) DIV 2` by rw[DIV_LE_MONOTONE, LESS_IMP_LESS_OR_EQ] >>
976  `n ** ((m * (m - 1)) DIV 2) <= n ** ((m ** 2) DIV 2)` by rw[EXP_BASE_LE_MONO] >>
977  metis_tac[LESS_LESS_EQ_TRANS]);
978
979(* ------------------------------------------------------------------------- *)
980(* Non-dividing coprimes of Product of Power Factors                         *)
981(* ------------------------------------------------------------------------- *)
982
983(* Define non-dividing coprimes *)
984val power_factors_coprime_nondivisor_def = Define`
985    power_factors_coprime_nondivisor n m = {x | coprime x n /\ ~(x divides (product_factors n m))}
986`;
987
988(* Theorem: x IN power_factors_coprime_nondivisor n m <=>
989            coprime x n /\ ~(x divides (product_factors n m)) *)
990(* Proof: by power_factors_coprime_nondivisor_def. *)
991val power_factors_coprime_nondivisor_element = store_thm(
992  "power_factors_coprime_nondivisor_element",
993  ``!n m x. x IN power_factors_coprime_nondivisor n m <=>
994           coprime x n /\ ~(x divides (product_factors n m))``,
995  rw[power_factors_coprime_nondivisor_def]);
996
997(* Theorem: 1 < m ==> !x. x IN power_factors_coprime_nondivisor n m ==> 1 < x *)
998(* Proof:
999   Since ~(x divides (product_factors n m))  by power_factors_coprime_nondivisor_element
1000         x <> 1                            by ONE_DIVIDES_ALL
1001      If x = 0, coprime x n ==> n = 1      by coprime_0L
1002    But then product_factors 1 m = 0       by product_factors_zero, 1 < m.
1003    making x divides (product_factors 1 m) by ALL_DIVIDES_0
1004    Hence x <> 0. Overall, 1 < x.
1005*)
1006val power_factors_coprime_nondivisor_property = store_thm(
1007  "power_factors_coprime_nondivisor_property",
1008  ``!n m. 1 < m ==> !x. x IN power_factors_coprime_nondivisor n m ==> 1 < x``,
1009  rw[power_factors_coprime_nondivisor_element] >>
1010  `x <> 1` by metis_tac[ONE_DIVIDES_ALL] >>
1011  `x <> 0` by metis_tac[product_factors_zero, coprime_0L, ALL_DIVIDES_0] >>
1012  decide_tac);
1013
1014(* Theorem: 1 < n /\ 1 < m ==> power_factors_coprime_nondivisor n m <> {} *)
1015(* Proof
1016   By power_factors_coprime_nondivisor_def, this is to show:
1017      1 < n /\ 1 < m ==> ?x. coprime x n /\ ~(x divides (product_factors n m))
1018   Let z = product_factors n m
1019       k = (m ** 2) DIV 2
1020   First, get the bounds of (product_factors n m):
1021      z < n ** k                  by product_factors_upper, 1 < n, 1 < m
1022      0 < z                       by product_factors_pos, 1 < n
1023      0 < 1 < k                   by ONE_LT_HALF_SQ
1024   Now, get a large prime,
1025      ?p. prime p /\ n < p        by NEXT_LARGER_PRIME, PRIME_INDEX
1026   Take x = p ** k, then
1027         coprime p n              by prime_coprime_all_lt, 0 < n /\ n < p
1028      or coprime (p ** k) n       by coprime_exp
1029    Also n ** k < p ** k          by EXP_EXP_LT_MONO, n < p /\ 0 < k
1030      so z < p ** k               by above
1031      or ~(p ** k <= z)
1032      or ~(divides (p ** k) z)    by DIVIDES_LE
1033*)
1034val power_factors_coprime_nondivisor_nonempty = store_thm(
1035  "power_factors_coprime_nondivisor_nonempty",
1036  ``!n m. 1 < n /\ 1 < m ==> power_factors_coprime_nondivisor n m <> {}``,
1037  rw[power_factors_coprime_nondivisor_def, EXTENSION] >>
1038  qabbrev_tac `z = product_factors n m` >>
1039  qabbrev_tac `k = (m ** 2) DIV 2` >>
1040  `z < n ** k` by rw[product_factors_upper, Abbr`z`, Abbr`k`] >>
1041  `0 < z` by rw[product_factors_pos, Abbr`z`] >>
1042  `1 < k` by rw[ONE_LT_HALF_SQ, Abbr`k`] >>
1043  `0 < k /\ 0 < n` by decide_tac >>
1044  `?p. prime p /\ n < p` by metis_tac[NEXT_LARGER_PRIME, PRIME_INDEX] >>
1045  `coprime p n` by rw[prime_coprime_all_lt] >>
1046  `coprime (p ** k) n` by rw[coprime_exp] >>
1047  `n ** k < p ** k` by rw_tac std_ss[EXP_EXP_LT_MONO] >>
1048  `~(p ** k <= z)` by decide_tac >>
1049  metis_tac[DIVIDES_LE]);
1050
1051(* Theorem: Any power_factors_coprime_nondivisor has order at least m.
1052   1 < m ==> !k. k IN (power_factors_coprime_nondivisor n m) ==> m <= ordz k n *)
1053(* Proof
1054   Since coprime k n /\
1055        ~(k divides (product_factors n m)) by power_factors_coprime_nondivisor_element
1056     and 1 < k                             by power_factors_coprime_nondivisor_property, 1 < m
1057      so ~(ordz k n < m)                   by proudct_factors_coprime_divisor_property
1058      or m <= ordz k n
1059*)
1060val power_factors_coprime_nondivisor_order = store_thm(
1061  "power_factors_coprime_nondivisor_order",
1062  ``!n m. 1 < m ==> !k. k IN (power_factors_coprime_nondivisor n m) ==> m <= ordz k n``,
1063  rpt strip_tac >>
1064  `coprime k n /\ ~(k divides (product_factors n m))` by metis_tac[power_factors_coprime_nondivisor_element] >>
1065  `1 < k` by metis_tac[power_factors_coprime_nondivisor_property] >>
1066  `~(ordz k n < m)` by metis_tac[proudct_factors_coprime_divisor_property] >>
1067  decide_tac);
1068
1069(* ------------------------------------------------------------------------- *)
1070(* LCM of coprimes divides Product of Power Factors                          *)
1071(* ------------------------------------------------------------------------- *)
1072
1073(* ------------------------------------------------------------------------- *)
1074(* Common Multiples of Residue n                                             *)
1075(* ------------------------------------------------------------------------- *)
1076
1077(* Define Common Multiple of residue n *)
1078val residue_common_multiple_def = Define`
1079  residue_common_multiple n = {m | 0 < m /\ !j. j IN residue n ==> j divides m}
1080`;
1081
1082(* Theorem: m IN residue_common_multiple n <=> 0 < m /\ !j. 0 < j /\ j < n ==> j divides m *)
1083(* Proof: by residue_common_multiple_def. *)
1084val residue_common_multiple_element = store_thm(
1085  "residue_common_multiple_element",
1086  ``!n m. m IN residue_common_multiple n <=> 0 < m /\ !j. 0 < j /\ j < n ==> j divides m``,
1087  rw[residue_common_multiple_def, residue_def]);
1088
1089(* Theorem: residue_common_multiple n <> {} *)
1090(* Proof:
1091   By residue_common_multiple_def and residue_def, this is to show:
1092      ?x. x <> 0 /\ !j. ((j = 0) \/ ~(j < n)) \/ j divides x
1093   Let x = FACT n, this is to show:
1094   (1) FACT n <> 0, true by FACT_LESS.
1095   (2) ((j = 0) \/ ~(j < n)) \/ j divides (FACT n)
1096       By contradiction, to show: j <> 0 /\ j < n /\ ~(j divides (FACT n)) ==> F
1097       But  0 < j /\ j <= n ==> j divides (FACT n)  by LEQ_DIVIDES_FACT
1098       Hence true by contradiction.
1099*)
1100val residue_common_multiple_nonempty = store_thm(
1101  "residue_common_multiple_nonempty",
1102  ``!n. residue_common_multiple n <> {}``,
1103  rw[residue_common_multiple_def, residue_def, EXTENSION] >>
1104  qexists_tac `FACT n` >>
1105  rpt strip_tac >| [
1106    `0 < FACT n` by rw[FACT_LESS] >>
1107    decide_tac,
1108    spose_not_then strip_assume_tac >>
1109    `0 < j /\ j <= n` by decide_tac >>
1110    metis_tac[LEQ_DIVIDES_FACT]
1111  ]);
1112
1113(* Theorem: 0 NOTIN residue_common_multiple n *)
1114(* Proof: by residue_common_multiple_def *)
1115val residue_common_multiple_nonzero = store_thm(
1116  "residue_common_multiple_nonzero",
1117  ``!n. 0 NOTIN residue_common_multiple n``,
1118  rw[residue_common_multiple_def]);
1119
1120(* Theorem: m IN residue_common_multiple n ==>
1121            !k. 0 < k ==> k * m IN residue_common_multiple n *)
1122(* Proof:
1123   By residue_common_multiple_def, this is to show:
1124   (1) 0 < m /\ 0 < k ==> 0 < k * m
1125       True by ZERO_LESS_MULT.
1126   (2) j IN residue n /\ (!j. j IN residue n ==> j divides m) ==> j divides (k * m)
1127       Condition implies j divides m, hence true by DIVIDES_MULTIPLE.
1128*)
1129Theorem residue_common_multiple_has_multiple:
1130  !n m. m IN residue_common_multiple n ==>
1131        !k. 0 < k ==> k * m IN residue_common_multiple n
1132Proof
1133  rw[residue_common_multiple_def, DIVIDES_MULTIPLE]
1134QED
1135
1136(* Theorem: b < a /\ a IN residue_common_multiple n /\ b IN residue_common_multiple n ==>
1137           (a - b) IN residue_common_multiple n *)
1138(* Proof:
1139   By residue_common_multiple_def, this is to show:
1140   (1) b < a ==> 0 < a - b
1141       True by arithmetic.
1142   (2) j IN residue n /\ (!j. j IN residue n ==> j divides a) /\
1143                         (!j. j IN residue n ==> j divides b) ==> j divides (a - b)
1144       Condition implies j divides a /\ j divides b, hence true by DIVIDES_SUB.
1145*)
1146val residue_common_multiple_has_sub = store_thm(
1147  "residue_common_multiple_has_sub",
1148  ``!n a b. b < a /\ a IN residue_common_multiple n /\ b IN residue_common_multiple n ==>
1149           (a - b) IN residue_common_multiple n``,
1150  rw[residue_common_multiple_def, DIVIDES_SUB]);
1151
1152(* Theorem: m IN residue_common_multiple n <=> 0 < m /\ !j. 1 < j /\ j < n ==> j divides m *)
1153(* Proof:
1154   If part: m IN residue_common_multiple n ==> 0 < m /\ !j. 1 < j /\ j < n ==> j divides m
1155      True by residue_common_multiple_def, as that gives the range 0 < j /\ j < n.
1156   Only-if part: 0 < m /\ !j. 1 < j /\ j < n ==> j divides m ==> m IN residue_common_multiple n
1157   The implication    0 < m /\ !j. 1 < j /\ j < n ==> j divides m)
1158   can be extended to 0 < m /\ !j. 0 < j /\ j < n ==> j divides m)
1159   by putting j = 1, and 1 divides m is true       by ONE_DIVIDES_ALL, 1 < n.
1160   Thus m IN residue_common_multiple n.
1161*)
1162val residue_common_multiple_element_1 = store_thm(
1163  "residue_common_multiple_element_1",
1164  ``!n m. 1 < n ==> (m IN residue_common_multiple n <=> 0 < m /\ !j. 1 < j /\ j < n ==> j divides m)``,
1165  rw[residue_common_multiple_def, residue_def] >>
1166  rw[EQ_IMP_THM] >>
1167  Cases_on `j = 1` >> rw[]);
1168
1169(* ------------------------------------------------------------------------- *)
1170(* Least Common Multiple of Residue n                                        *)
1171(* ------------------------------------------------------------------------- *)
1172
1173(* Define LCM of 1..(n-1) = residue n *)
1174val residue_lcm_def = Define`
1175  residue_lcm n = MIN_SET (residue_common_multiple n)
1176`;
1177
1178(* Theorem: The residue_lcm divides any common multiple
1179            m IN residue_common_multiple n ==> (residue_lcm n) divides m *)
1180(* Proof:
1181   By residue_lcm_def, this is to show:
1182      m IN residue_common_multiple n ==> (MIN_SET (residue_common_multiple n)) divides m
1183   Let s = residue_common_multiple n, k = MIN_SET s.
1184    Then m IN s ==> s <> {}                 by MEMBER_NOT_EMPTY
1185   Hence k IN s /\ !x. x IN s ==> k <= x    by MIN_SET_LEM
1186     and k IN s ==> 0 < k /\ !j. 0 < j /\ j < n ==> j divides k  by residue_common_multiple_element
1187   Now, ?r q. (m = q * k + r) /\ r < k      by DA, division algorithm
1188   If r = 0,
1189      then m = q * k                        by ADD_0
1190        or k divides m                      by divides_def
1191   If r <> 0, 0 < r.
1192      If q = 0,
1193        then m = r                          by MULT, ADD
1194          or r IN s, hence k <= r           by above
1195        This contradicts r < k.
1196      If q <> 0, 0 < q.
1197        q * k < m /\ (r = m - q * k)        by 0 < r
1198        q * k IN s                          by residue_common_multiple_has_multiple, MULT_COMM, 0 < q
1199        so r IN s, hence k <= r             by residue_common_multiple_has_sub, q * k < m
1200        or k <= r
1201        This contradicts r < k.
1202*)
1203val residue_lcm_divides_common_multiple = store_thm(
1204  "residue_lcm_divides_common_multiple",
1205  ``!n m. m IN residue_common_multiple n ==> (residue_lcm n) divides m``,
1206  rw[residue_lcm_def] >>
1207  qabbrev_tac `s = residue_common_multiple n` >>
1208  qabbrev_tac `k = MIN_SET s` >>
1209  `s <> {}` by metis_tac[MEMBER_NOT_EMPTY] >>
1210  `k IN s /\ !x. x IN s ==> k <= x` by rw[MIN_SET_LEM, Abbr`k`] >>
1211  `0 < k /\ !j. 0 < j /\ j < n ==> j divides k` by metis_tac[residue_common_multiple_element] >>
1212  `?r q. (m = q * k + r) /\ r < k` by rw[DA] >>
1213  Cases_on `r = 0` >| [
1214    `m = q * k` by decide_tac >>
1215    metis_tac[divides_def],
1216    Cases_on `q = 0` >| [
1217      `m = r` by rw[] >>
1218      `k <= r` by metis_tac[] >>
1219      decide_tac,
1220      `0 < q /\ 0 < r /\ q * k < m /\ (r = m - q * k)` by decide_tac >>
1221      `q * k IN s` by rw[residue_common_multiple_has_multiple, MULT_COMM, Abbr`s`] >>
1222      `r IN s` by rw[residue_common_multiple_has_sub, Abbr`s`] >>
1223      `k <= r` by metis_tac[] >>
1224      decide_tac
1225    ]
1226  ]);
1227
1228(* ------------------------------------------------------------------------- *)
1229(* Relate Residue LCM to List LCM of Leibniz Triangle                        *)
1230(* ------------------------------------------------------------------------- *)
1231
1232(*
1233Idea:
1234If 1...(n-1) divides c
1235and !m. 1...(n-1) divides m ==> c divides m
1236Then c = residue_lcm n.
1237
1238Then prove that c = list_lcm (leibniz_vertical (n-1)).
1239*)
1240
1241(* Note: the following proof repeats a bit of the previous proof. *)
1242
1243(* Theorem: (c = residue_lcm n) <=>
1244            c IN (residue_common_multiple n) /\
1245            !m. m IN (residue_common_multiple n) ==> c divides m *)
1246(* Proof:
1247   Note residue_common_multiple n <> {}    by residue_common_multiple_nonempty
1248   By residue_lcm_def, this is to show:
1249   (1) MIN_SET (residue_common_multiple n) IN residue_common_multiple n
1250       True by MIN_SET_LEM.
1251   (2) m IN residue_common_multiple n ==> (MIN_SET (residue_common_multiple n)) divides m
1252       Let  c = MIN_SET (residue_common_multiple n)
1253       Then c IN residue_common_multiple n   by MIN_SET_LEM
1254         so 0 < c                            by residue_common_multiple_element
1255        and ?q r. m = q * c + r, with r < c  by DA, division algorithm, 0 < c.
1256       If r = 0,
1257         then m = q * c, so c divides m      by divides_def
1258       If r <> 0, 0 < r.
1259           and c <= m                        by MIN_SET_LEM
1260         hence q <> 0, otherwise m = r < c, a contradiction.
1261         Thus 0 < q,
1262         and q * c IN (residue_common_multiple n)  by residue_common_multiple_has_multiple
1263         Therefore, r = m - q * c, and q * c < m   by arithmetic
1264         so      r IN (residue_common_multiple n)` by residue_common_multiple_has_sub, q*c < m.
1265         giving  c <= r                      by MIN_SET_LEM
1266         which contradicts r < c.
1267   (3) !m. m IN residue_common_multiple n ==> c divides m ==> c = MIN_SET (residue_common_multiple n)
1268       Let z = MIN_SET (residue_common_multiple n)
1269       Then z IN residue_common_multiple n   by MIN_SET_LEM
1270        and z <= c                           by MIN_SET_LEM
1271        and 0 < z                            by residue_common_multiple_element
1272       With z IN residue_common_multiple n,
1273             divides c z                     by implication
1274       hence c <= z                          by DIVIDES_LE
1275       Therefore c = z                       by EQ_LESS_EQ
1276*)
1277val residue_lcm_property = store_thm(
1278  "residue_lcm_property",
1279  ``!n c. (c = residue_lcm n) <=>
1280      c IN (residue_common_multiple n) /\
1281      !m. m IN (residue_common_multiple n) ==> c divides m``,
1282  rpt strip_tac >>
1283  `residue_common_multiple n <> {}` by rw[residue_common_multiple_nonempty] >>
1284  rw[residue_lcm_def, EQ_IMP_THM] >| [
1285    rw[MIN_SET_LEM],
1286    qabbrev_tac `c = MIN_SET (residue_common_multiple n)` >>
1287    `c IN residue_common_multiple n` by rw[MIN_SET_LEM, Abbr`c`] >>
1288    `0 < c` by metis_tac[residue_common_multiple_element] >>
1289    `?r q. (m = q * c + r) /\ (r < c)` by rw[DA] >>
1290    Cases_on `r = 0` >| [
1291      `m = q * c` by decide_tac >>
1292      metis_tac[divides_def],
1293      `0 < r` by decide_tac >>
1294      `c <= m` by rw[MIN_SET_LEM, Abbr`c`] >>
1295      Cases_on `q = 0` >| [
1296        `m = r` by rw[] >>
1297        `m < c` by decide_tac >>
1298        decide_tac,
1299        `0 < q` by decide_tac >>
1300        `q * c IN (residue_common_multiple n)` by rw[residue_common_multiple_has_multiple] >>
1301        `(r = m - q * c) /\ q * c < m` by decide_tac >>
1302        `r IN (residue_common_multiple n)` by rw[residue_common_multiple_has_sub] >>
1303        `c <= r` by rw[MIN_SET_LEM, Abbr`c`] >>
1304        decide_tac
1305      ]
1306    ],
1307    qabbrev_tac `z = MIN_SET (residue_common_multiple n)` >>
1308    `z IN residue_common_multiple n` by rw[MIN_SET_LEM, Abbr`z`] >>
1309    `z <= c` by rw[MIN_SET_LEM, Abbr`z`] >>
1310    `0 < z` by metis_tac[residue_common_multiple_element] >>
1311    `c <= z` by rw[DIVIDES_LE] >>
1312    decide_tac
1313  ]);
1314
1315(*
1316> |- !n. residue_lcm n = MIN_SET (residue_common_multiple n)
1317> |- !n. residue_common_multiple n = {m | 0 < m /\ !j. j IN residue n ==> j divides m}
1318> |- !n. residue n = {i | 0 < i /\ i < n}
1319Hence 1 < n for residue_lcm n to be defined, as MIN_SET {} is meaningless.
1320residue_lcm n = MIN_SET { commomn multiples of 1..(n-1) inclusive}
1321list_lcm (leibniz_vertical n) = least (common multiples of 1 to (n+1) inclusive)
1322so residue_lcm n = list_lcm (leibniz_vertical (n-2)) >= 2 ** (n-2) = (2 ** n)/4
1323*)
1324
1325(* Theorem: 1 < n ==> (residue_lcm n = list_lcm (leibniz_vertical (n-2))) *)
1326(* Proof:
1327   Let c = list_lcm (leibniz_vertical (n-2)).
1328   Then 0 < c                              by list_lcm_pos, leibniz_vertical_pos
1329   We have !j. j < n <=> j <= SUC (n-2)    by 1 < n
1330        so !j. 0 < j /\ j < n <=> MEM j (leibniz_vertical (n-2))   by leibniz_vertical_mem
1331       and !j. 0 < j /\ j < n ==> j divides c                      by list_lcm_is_common_multiple
1332   Hence c IN (residue_common_multiple n)                          by residue_common_multiple_element
1333    Also !m. m IN (residue_common_multiple n)
1334         ==> 0 < m /\ (!j. MEM j (leibniz_vertical (n-2)) ==> j divides m) by residue_common_multiple_element
1335     and !m. m IN (residue_common_multiple n) ==> c divides m      by list_lcm_is_least_common_multiple
1336   Therefore c = residue_lcm n                                     by residue_lcm_property
1337*)
1338val residue_lcm_alt = store_thm(
1339  "residue_lcm_alt",
1340  ``!n. 1 < n ==> (residue_lcm n = list_lcm (leibniz_vertical (n-2)))``,
1341  rpt strip_tac >>
1342  qabbrev_tac `c = list_lcm (leibniz_vertical (n-2))` >>
1343  `0 < c` by rw[list_lcm_pos, leibniz_vertical_pos, Abbr`c`] >>
1344  `!j. j < n <=> j <= SUC (n-2)` by decide_tac >>
1345  `!j. 0 < j /\ j < n <=> MEM j (leibniz_vertical (n-2))` by rw[GSYM leibniz_vertical_mem] >>
1346  `!j. 0 < j /\ j < n ==> j divides c` by rw[list_lcm_is_common_multiple, Abbr`c`] >>
1347  metis_tac[residue_common_multiple_element, list_lcm_is_least_common_multiple, residue_lcm_property]);
1348
1349(* Theorem: 1 < n ==> (residue_lcm n = lcm_run (n - 1)) *)
1350(* Proof: by residue_lcm_alt *)
1351val residue_lcm_eq_list_lcm = store_thm(
1352  "residue_lcm_eq_list_lcm",
1353  ``!n. 1 < n ==> (residue_lcm n = list_lcm [1 .. (n - 1)])``,
1354  rw[residue_lcm_alt]);
1355
1356(* ------------------------------------------------------------------------- *)
1357(* Bounds for Residue LCM                                                    *)
1358(* ------------------------------------------------------------------------- *)
1359
1360(* Theorem: 1 < n /\ 1 < m /\ 1 < k /\
1361            (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1362            residue_lcm k <= (product_factors n m) *)
1363(* Proof:
1364   Since !j. 1 < j /\ j < k /\ coprime j n /\ ordz j n < m
1365     ==> j divides (product_factors n m)                 by product_factors_divisors
1366    and  0 < product_factors n m                         by product_factors_pos, 1 < n
1367     so  (product_factors n m) IN residue_common_multiple k     by residue_common_multiple_element_1, 1 < k
1368   Hence (residue_lcm k) divides (product_factors n m)   by residue_lcm_divides_common_multiple
1369      or (residue_lcm k) <= (product_factors n m)        by DIVIDES_LE, 0 < product_factors n m
1370*)
1371val product_factors_lower = store_thm(
1372  "product_factors_lower",
1373  ``!n m k. 1 < n /\ 1 < m /\ 1 < k /\
1374   (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==> residue_lcm k <= (product_factors n m)``,
1375  rpt strip_tac >>
1376  `!j. 1 < j /\ j < k /\ coprime j n /\ ordz j n < m ==> j divides (product_factors n m)` by rw[product_factors_divisors] >>
1377  `0 < product_factors n m` by rw[product_factors_pos] >>
1378  `(product_factors n m) IN residue_common_multiple k` by rw[residue_common_multiple_element_1] >>
1379  `(residue_lcm k) divides (product_factors n m)` by rw[residue_lcm_divides_common_multiple] >>
1380  rw[DIVIDES_LE]);
1381
1382(* Theorem: 1 < n /\ 1 < m /\ 1 < k /\
1383            (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1384            residue_lcm k < n ** ((m ** 2) DIV 2) *)
1385(* Proof:
1386   Since (residue_lcm k) <= (product_factors n m)        by product_factors_lower
1387     and product_factors n m < n ** ((m ** 2) DIV 2)     by product_factors_upper
1388   Therefore (residue_lcm k) < n ** ((m ** 2) DIV 2)     by LESS_EQ_LESS_TRANS
1389*)
1390val residue_lcm_upper = store_thm(
1391  "residue_lcm_upper",
1392  ``!n m k. 1 < n /\ 1 < m /\ 1 < k /\
1393   (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1394   residue_lcm k < n ** ((m ** 2) DIV 2)``,
1395  metis_tac[product_factors_upper, product_factors_lower, LESS_EQ_LESS_TRANS]);
1396
1397(* Theorem: 1 < n ==> 2 ** (n - 2) <= residue_lcm n *)
1398(* Proof:
1399     residue_lcm n
1400   = list_lcm (leibniz_vertical (n-2))   by residue_lcm_alt, 1 < n.
1401   >= 2 ** (n-2)                         by leibniz_vertical_lcm_lower
1402*)
1403val residue_lcm_lower = store_thm(
1404  "residue_lcm_lower",
1405  ``!n. 1 < n ==> 2 ** (n - 2) <= residue_lcm n``,
1406  rpt strip_tac >>
1407  `residue_lcm n = list_lcm (leibniz_vertical (n-2))` by rw[residue_lcm_alt] >>
1408  metis_tac[leibniz_vertical_lcm_lower]);
1409
1410(* Question: How big must k be in order that m <= ordz k n ?
1411   Searching from k = 1, k = 2, etc.
1412   Consider the j's in (residue k), i.e. 0 < j < k.
1413   Assume for all such j, no such candidate k is found.
1414   That is, because for these j's, ordz j n < m.
1415   or the following condition is satisfied:
1416   !j. 0 < j /\ j < k /\ coprime j n /\ ordz j n < m
1417   This gives: residue_lcm k < n ** ((m ** 2) DIV 2)    by residue_lcm_upper, independent of k.
1418   But always: 2 ** (k-2) <= residue_lcm k              by residue_lcm_lower, increases with k.
1419   Also: n ** ((m ** 2) DIV 2) < 2 ** ((SUC (LOG2 n)) * ((m ** 2) DIV 2))  by integer arithmetic
1420   Therefore   2 ** (k-2) < 2 ** ((SUC (LOG2 n)) * ((m ** 2) DIV 2))
1421   or                k-2 < (SUC (LOG2 n)) * ((m ** 2) DIV 2)
1422   or                  k < (SUC (LOG2 n)) * ((m ** 2) DIV 2) + 2 = f n m
1423   This means, the situation is possible only because k < f n m.
1424   So if you search from k = 1, k = 2, up to k = f n m,
1425   Then such situations are impossible, or there is some k such that m <= ordz k n
1426*)
1427
1428(* Theorem:  1 < n /\ 1 < m /\ 1 < k /\
1429             (!j. 1 < j /\ j < k /\ coprime j n /\ ordz j n < m) ==>
1430             k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) *)
1431(* Proof:
1432   Since residue_lcm k < n ** ((m ** 2) DIV 2)      by residue_lcm_upper, 1 < n, 1 < m
1433     and 2 ** (k-2) <= residue_lcm k                by residue_lcm_lower, 1 < k
1434   Step 1: show n ** ((m ** 2) DIV 2) < 2 ** SUC (LOG2 n)) ** ((m ** 2) DIV 2
1435   Since 1 < (m ** 2) DIV 2                         by ONE_LT_HALF_SQ, 1 < m
1436     and n < 2 ** SUC (LOG2 n)                      by LOG2_PROPERTY, 0 < n
1437      so n ** ((m ** 2) DIV 2)
1438         < (2 ** SUC (LOG2 n)) ** ((m ** 2) DIV 2)  by EXP_EXP_LT_MONO
1439         = 2 ** ((SUC (LOG2 n)) * ((m ** 2) DIV 2)) by EXP_EXP_MULT
1440   Step 2: conclusion
1441   Hence 2 ** (k-2) < 2 ** (SUC (LOG2 n) * (m ** 2 DIV 2))
1442      or k - 2 < SUC (LOG2 n) * (m ** 2 DIV 2)      by EXP_BASE_LT_MONO, 1 < 2
1443      or k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2)      by arithmetic
1444*)
1445val ZN_order_modulus_upper = store_thm(
1446  "ZN_order_modulus_upper",
1447  ``!n m k. 1 < n /\ 1 < m /\ 1 < k /\
1448    (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1449    k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2)``,
1450  rpt strip_tac >>
1451  `residue_lcm k < n ** ((m ** 2) DIV 2)` by rw[residue_lcm_upper] >>
1452  `2 ** (k-2) <= residue_lcm k` by rw[residue_lcm_lower] >>
1453  `1 < (m ** 2) DIV 2` by rw[ONE_LT_HALF_SQ] >>
1454  `0 < (m ** 2) DIV 2 /\ 0 < n /\ 0 < k` by decide_tac >>
1455  `n < 2 ** SUC (LOG2 n)` by rw[LOG2_PROPERTY] >>
1456  `n ** ((m ** 2) DIV 2) < (2 ** SUC (LOG2 n)) ** ((m ** 2) DIV 2)` by rw[EXP_EXP_LT_MONO] >>
1457  `(2 ** SUC (LOG2 n)) ** ((m ** 2) DIV 2) = 2 ** ((SUC (LOG2 n)) * ((m ** 2) DIV 2))` by rw[EXP_EXP_MULT] >>
1458  `2 ** (k - 2) < 2 ** (SUC (LOG2 n) * (m ** 2 DIV 2))` by decide_tac >>
1459  `k - 2 < SUC (LOG2 n) * (m ** 2 DIV 2)` by metis_tac[EXP_BASE_LT_MONO, DECIDE``1 < 2``] >>
1460  decide_tac);
1461
1462(* Theorem: ZN modulus for order of n to be at least m exists.
1463            1 < n /\ 1 < m ==>
1464            ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\
1465            (coprime k n ==> m <= ordz k n) *)
1466(* Proof:
1467   Let h = 2 + SUC (LOG2 n) * (m ** 2 DIV 2)
1468   By contradiction, assume !k. 1 < k /\ k < h /\ ~(coprime k n ==> m <= ordz k n)).
1469    Since  1 < h       by arithmetic
1470      and ~(p => q) == ~(~p or q) == p and ~q,
1471    where ~(m <= ordz k n)) means ordz k n < m  by NOT_LESS_EQUAL
1472       so !k. 1 < k /\ k < h ==> coprime k n /\ ordz k n < m
1473   giving  h < h       by ZN_order_modulus_upper, 1 < h.
1474   which is false      by prim_recTheory.LESS_REFL
1475*)
1476val ZN_order_modulus_exists = store_thm(
1477  "ZN_order_modulus_exists",
1478  ``!n m. 1 < n /\ 1 < m ==>
1479    ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\ (coprime k n ==> m <= ordz k n)``,
1480  rpt strip_tac >>
1481  `1 < 2 + SUC (LOG2 n) * (m ** 2 DIV 2)` by decide_tac >>
1482  metis_tac[ZN_order_modulus_upper, NOT_LESS_EQUAL, DECIDE``!n. ~(n < n)``]);
1483
1484(* Theorem: ZN modulus for order of n to be at least m exists.
1485            1 < n /\ 1 < m ==>
1486            ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\
1487            (~ coprime k n \/ m <= ordz k n) *)
1488(* Proof: by ZN_order_modulus_exists, and p ==> q == ~p or q. *)
1489val ZN_order_modulus_exists_alt = store_thm(
1490  "ZN_order_modulus_exists_alt",
1491  ``!n m. 1 < n /\ 1 < m ==>
1492    ?k. 1 < k /\ k < 2 + SUC (LOG2 n) * (m ** 2 DIV 2) /\ (~ coprime k n \/ m <= ordz k n)``,
1493  metis_tac[ZN_order_modulus_exists]);
1494
1495(* Theorem: 0 < n /\ prime k /\ m <= ordz k n ==> ~(k divides (product_factors n m)) *)
1496(* Proof:
1497   Since p divides (product_factors n m) ==> ordz p n < m   by proudct_factors_prime_divisor_property
1498   So given m <= ordz k n, p divides (product_factors n m) is false.
1499*)
1500val ZN_order_modulus_when_prime = store_thm(
1501  "ZN_order_modulus_when_prime",
1502  ``!n m k. 0 < n /\ prime k /\ m <= ordz k n ==> ~(k divides (product_factors n m))``,
1503  metis_tac[proudct_factors_prime_divisor_property, DECIDE``!n m. ~(n < m) <=> m <= n``]);
1504
1505(* ------------------------------------------------------------------------- *)
1506(* ZN bounds based on upper logarithm                                        *)
1507(* ------------------------------------------------------------------------- *)
1508
1509(* Theorem: 1 < n /\ 1 < m /\ 1 < k /\
1510            (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1511            k < 2 + (ulog n) * HALF (m ** 2) *)
1512(* Proof:
1513   Let h = HALF (m ** 2)
1514   Note residue_lcm k < n ** h               by residue_lcm_upper, 1 < n, 1 < m
1515    and 2 ** (k - 2) <= residue_lcm k        by residue_lcm_lower, 1 < k
1516
1517   Step 1: show n ** h <= 2 ** ((ulog n) * h)
1518   Note      n <= 2 ** ulog n                by ulog_property, 0 < n
1519        n ** h <= (2 ** ulog n) ** h         by EXP_EXP_LE_MONO
1520                = 2 ** ((ulog n) * h))       by EXP_EXP_MULT
1521
1522   Step 2: conclusion
1523   Hence 2 ** (k - 2) < 2 ** ((ulog n) * h)
1524      or        k - 2 < (ulog n) * h         by EXP_BASE_LT_MONO, 1 < 2
1525      or            k < 2 + (ulog n) * h     by arithmetic
1526*)
1527val ZN_order_modulus_upper_1 = store_thm(
1528  "ZN_order_modulus_upper_1",
1529  ``!n m k. 1 < n /\ 1 < m /\ 1 < k /\
1530    (!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m) ==>
1531    k < 2 + (ulog n) * HALF (m ** 2)``,
1532  rpt strip_tac >>
1533  qabbrev_tac `h = HALF (m ** 2)` >>
1534  `residue_lcm k < n ** h` by rw[residue_lcm_upper, Abbr`h`] >>
1535  `2 ** (k - 2) <= residue_lcm k` by rw[residue_lcm_lower] >>
1536  `n <= 2 ** (ulog n)` by rw[ulog_property] >>
1537  `n ** h <= (2 ** ulog n) ** h` by rw[EXP_EXP_LE_MONO] >>
1538  `(2 ** ulog n) ** h = 2 ** ((ulog n) * h)` by rw[EXP_EXP_MULT] >>
1539  `2 ** (k - 2) < 2 ** ((ulog n) * h)` by decide_tac >>
1540  `k - 2 < (ulog n) * h` by metis_tac[EXP_BASE_LT_MONO, DECIDE``1 < 2``] >>
1541  decide_tac);
1542
1543(* Theorem: ZN modulus for order of n to be at least m exists if coprime.
1544            1 < n /\ 1 < m ==>
1545            ?k. 1 < k /\ k < 2 + (ulog n) * HALF (m ** 2) /\ (~ coprime k n \/ m <= ordz k n) *)
1546(* Proof:
1547   Let h = 2 + (ulog n) * HALF (m ** 2).
1548   By contradiction, assume !k. 1 < k /\ k < h /\ ~(coprime k n ==> m <= ordz k n)).
1549    Since  1 < h             by arithmetic
1550      and ~(p => q) == ~(~p or q) == p and ~q,
1551    where ~(m <= ordz k n)) means ordz k n < m  by NOT_LESS_EQUAL
1552       so !k. 1 < k /\ k < h ==> coprime k n /\ ordz k n < m
1553   giving  h < h             by ZN_order_modulus_upper_1, 1 < h.
1554   which is false            by prim_recTheory.LESS_REFL
1555*)
1556val ZN_order_modulus_exists_1 = store_thm(
1557  "ZN_order_modulus_exists_1",
1558  ``!n m. 1 < n /\ 1 < m ==>
1559    ?k. 1 < k /\ k < 2 + (ulog n) * HALF (m ** 2) /\ (~ coprime k n \/ m <= ordz k n)``,
1560  rpt strip_tac >>
1561  `1 < 2 + (ulog n) * HALF (m ** 2)` by decide_tac >>
1562  metis_tac[ZN_order_modulus_upper_1, NOT_LESS_EQUAL, DECIDE``!n. ~(n < n)``]);
1563
1564(* Theorem: 1 < n /\ 1 < m /\ 1 < k /\
1565            (!j. 1 < j /\ j < k ==> ~(j divides n) /\ ordz j n < m) ==>
1566            k < 2 + HALF (m ** 2) * (ulog n) *)
1567(* Proof: by coprime_by_le_not_divides, ZN_order_modulus_upper_1 *)
1568val ZN_order_modulus_upper_2 = store_thm(
1569  "ZN_order_modulus_upper_2",
1570  ``!n m k. 1 < n /\ 1 < m /\ 1 < k /\
1571    (!j. 1 < j /\ j < k ==> ~(j divides n) /\ ordz j n < m) ==>
1572    k < 2 + HALF (m ** 2) * (ulog n)``,
1573  rpt strip_tac >>
1574  `!j. 1 < j /\ j < k ==> coprime j n /\ ordz j n < m` by rw[coprime_by_le_not_divides] >>
1575  metis_tac[ZN_order_modulus_upper_1, MULT_COMM]);
1576
1577(* Theorem: ZN modulus for order of n to be at least m exists if coprime.
1578            1 < n /\ 1 < m ==>
1579            ?k. 1 < k /\ k < 2 + HALF (m ** 2) * (ulog n) /\ (k divides n \/ m <= ordz k n) *)
1580(* Proof:
1581   Let h = 2 + HALF (m ** 2) * (ulog n).
1582   By contradiction, assume !k. 1 < k /\ k < h ==> ~(k divides n) /\ ~(m <= ordz k n).
1583    Since  1 < h             by arithmetic
1584      and ~(m <= ordz k n)) means ordz k n < m  by NOT_LESS_EQUAL
1585       so !k. 1 < k /\ k < h ==> ~(k divides n) /\ ordz k n < m
1586   giving  h < h             by ZN_order_modulus_upper_2, 1 < h.
1587   which is false            by prim_recTheory.LESS_REFL
1588*)
1589val ZN_order_modulus_exists_2 = store_thm(
1590  "ZN_order_modulus_exists_2",
1591  ``!n m. 1 < n /\ 1 < m ==>
1592    ?k. 1 < k /\ k < 2 + HALF (m ** 2) * (ulog n) /\ (k divides n \/ m <= ordz k n)``,
1593  rpt strip_tac >>
1594  `1 < 2 + HALF (m ** 2) * (ulog n)` by decide_tac >>
1595  metis_tac[ZN_order_modulus_upper_2, NOT_LESS_EQUAL, DECIDE``!n. ~(n < n)``]);
1596
1597(* Theorem: 1 < n /\ 1 < m /\
1598            (!k. 1 < k /\ k <= c ==> ~(k divides n) /\ ordz k n < m) ==> c < 1 + HALF (m ** 2) * (ulog n) *)
1599(* Proof: by ZN_order_modulus_upper_2 *)
1600val ZN_order_modulus_upper_3 = store_thm(
1601  "ZN_order_modulus_upper_3",
1602  ``!n m k c. 1 < n /\ 1 < m /\
1603    (!k. 1 < k /\ k <= c ==> ~(k divides n) /\ ordz k n < m) ==> c < 1 + HALF (m ** 2) * (ulog n)``,
1604  rpt strip_tac >>
1605  `!j. 1 < j /\ j < c + 1 ==> ~(j divides n) /\ ordz j n < m` by rw[] >>
1606  (Cases_on `c = 0` >> simp[]) >>
1607  `1 < c + 1` by rw[] >>
1608  `c + 1 < 2 + HALF (m ** 2) * ulog n` by rw[ZN_order_modulus_upper_2] >>
1609  decide_tac);
1610
1611(* Theorem: 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * (ulog n) <= c ==>
1612            ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n) *)
1613(* Proof: by ZN_order_modulus_upper_3 *)
1614val ZN_order_modulus_exists_3 = store_thm(
1615  "ZN_order_modulus_exists_3",
1616  ``!n m c. 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * (ulog n) <= c ==>
1617    ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)``,
1618  spose_not_then strip_assume_tac >>
1619  `c < 1 + HALF (m ** 2) * ulog n` by metis_tac[ZN_order_modulus_upper_3, NOT_LESS_EQUAL] >>
1620  decide_tac);
1621
1622(* Theorem: 1 < n /\ 1 < m /\ (c = 1 + HALF (m ** 2) * (ulog n)) ==>
1623            ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n) *)
1624(* Proof: by ZN_order_modulus_exists_3 *)
1625val ZN_order_modulus_exists_4 = store_thm(
1626  "ZN_order_modulus_exists_4",
1627  ``!n m c. 1 < n /\ 1 < m /\ (c = 1 + HALF (m ** 2) * (ulog n)) ==>
1628    ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)``,
1629  rw[ZN_order_modulus_exists_3]);
1630
1631(* The following is the generic form -- cannot move to ringInstances due to SQRT_LE, SQRT_OF_SQ, etc. *)
1632
1633(* Theorem: 1 < n /\ coprime k n /\ (2 * SUC (LOG2 n)) ** 2 <= ordz k n ==>
1634            2 * (SQRT (phi k)) * SUC (LOG2 n) <= phi k *)
1635(* Proof:
1636   Let j = ordz k n
1637       s = 2 * (SQRT (phi k)) * SUC (LOG2 n)
1638   Note 0 < k                            by coprime_gt_1, 1 < n
1639
1640   Step 1: show SQRT j <= SQRT (phi k)
1641   Since j divides (phi k)               by ZN_coprime_order_divides_phi, 0 < k, coprime k n
1642     and 0 < phi k                       by phi_pos, 0 < k
1643      so j <= phi k                      by DIVIDES_LE, 0 < phi k
1644    Thus SQRT j <= SQRT (phi k)          by SQRT_LE
1645
1646   Step 2: show s <= SQRT (phi k) * SQRT (phi k)
1647   Since 0 < LOG2 n                          by LOG2_POS, 1 < n
1648      so 0 < 2 * SUC (LOG2 n)
1649    Now, 2 * (SUC (LOG2 n)) ** 2 <= j        by given
1650      so 2 * SUC (LOG2 n) <= SQRT j          by SQRT_OF_SQ, SQRT_LE, 0 < 2 * SUC (LOG2 n)
1651      or 2 * SUC (LOG2 n) <= SQRT (phi k)    from Step 1
1652   Given s = 2 * (SQRT (phi k)) * SUC (LOG2 n)
1653           = 2 * SUC (LOG2 n) * SQRT (phi k) by arithmetic
1654         s <= SQRT (phi k) * SQRT (phi k)    by LESS_MONO_MULT
1655
1656   Step 3: conclude
1657   Since SQ (SQRT (phi k) <= phi k           by SQ_SQRT_LE
1658   Hence s <= phi k                          from Step 2
1659*)
1660val ZN_order_condition_property = store_thm(
1661  "ZN_order_condition_property",
1662  ``!n k. 1 < n /\ coprime k n /\ (2 * SUC (LOG2 n)) ** 2 <= ordz k n ==>
1663         2 * (SQRT (phi k)) * SUC (LOG2 n) <= phi k``,
1664  rpt strip_tac >>
1665  qabbrev_tac `j = ordz k n` >>
1666  qabbrev_tac `s = 2 * (SQRT (phi k)) * SUC (LOG2 n)` >>
1667  `0 < k` by metis_tac[coprime_gt_1] >>
1668  `j <= phi k` by
1669  (`j divides (phi k)` by rw[ZN_coprime_order_divides_phi, Abbr`j`] >>
1670  `0 < phi k` by rw[phi_pos] >>
1671  metis_tac[DIVIDES_LE]) >>
1672  `SQRT j <= SQRT (phi k)` by rw[SQRT_LE] >>
1673  `s <= SQRT (phi k) * SQRT (phi k)` by
1674    (`0 < LOG2 n` by rw[] >>
1675  `0 < 2 * SUC (LOG2 n)` by decide_tac >>
1676  `2 * SUC (LOG2 n) <= SQRT j` by metis_tac[SQRT_OF_SQ, SQRT_LE] >>
1677  `2 * SUC (LOG2 n) <= SQRT (phi k)` by decide_tac >>
1678  `s = 2 * SUC (LOG2 n) * SQRT (phi k)` by rw_tac arith_ss[Abbr`s`] >>
1679  rw[LESS_MONO_MULT]) >>
1680  `SQ (SQRT (phi k)) <= phi k` by rw[SQ_SQRT_LE] >>
1681  decide_tac);
1682
1683(* Theorem: 1 < n /\ coprime k n ==>
1684            !h. 0 < h /\ h ** 2 <= ordz k n ==> (SQRT (phi k)) * h <= phi k *)
1685(* Proof:
1686   Idea:               h ** 2 <= ordz k n <= phi k     by ZN_coprime_order_divides_phi
1687   Taking square roots,     h <= SQRT (phi k)
1688   Thus      SQRT (phi k) * h <= SQRT (phi k) * SQRT (phi k)
1689     or      SQRT (phi k) * h <= phi k
1690
1691   Let t = phi k, s = SQRT t * h.
1692   The goal becomes: s <= t.
1693
1694   Note 0 < k                      by coprime_gt_1, 1 < n
1695    and 0 < t                      by phi_pos, 0 < k
1696
1697   Note (ordz k n) divides t       by ZN_coprime_order_divides_phi, 0 < k /\ coprime k n
1698     so       ordz k n <= t        by DIVIDES_LE, 0 < t
1699     or         h ** 2 <= t        by given, h ** 2 <= ordz k n
1700    ==>  SQRT (h ** 2) <= SQRT t   by SQRT_LE
1701     or              h <= SQRT t   by SQRT_OF_SQ
1702    ==>     SQRT t * h <= SQRT t * SQRT t    by LESS_MONO_MULT
1703     or              s <= t        by SQ_SQRT_LE, notation
1704*)
1705val ZN_order_condition_property_0 = store_thm(
1706  "ZN_order_condition_property_0",
1707  ``!n k. 1 < n /\ coprime k n ==>
1708   !h. 0 < h /\ h ** 2 <= ordz k n ==> (SQRT (phi k)) * h <= phi k``,
1709  rpt strip_tac >>
1710  qabbrev_tac `t = phi k` >>
1711  qabbrev_tac `s = SQRT t * h` >>
1712  `0 < k` by metis_tac[coprime_gt_1] >>
1713  `0 < t` by rw[phi_pos, Abbr`t`] >>
1714  `(ordz k n) divides t` by rw[ZN_coprime_order_divides_phi, Abbr`t`] >>
1715  `ordz k n <= t` by rw[DIVIDES_LE] >>
1716  `h ** 2 <= t` by decide_tac >>
1717  `SQRT (h ** 2) <= SQRT t` by rw[SQRT_LE] >>
1718  `h <= SQRT t` by metis_tac[GSYM SQRT_LE, SQRT_OF_SQ] >>
1719  `s <= SQRT t * SQRT t` by rw[LESS_MONO_MULT, Abbr`s`] >>
1720  `SQ (SQRT t) <= t` by rw[SQ_SQRT_LE] >>
1721  decide_tac);
1722
1723(* Theorem: 1 < n /\ coprime k n /\ (SUC (LOG2 n)) ** 2 <= ordz k n ==>
1724         SQRT (phi k) * SUC (LOG2 n) <= phi k *)
1725(* Proof: by ZN_order_condition_property_0 *)
1726val ZN_order_condition_property_1 = store_thm(
1727  "ZN_order_condition_property_1",
1728  ``!n k. 1 < n /\ coprime k n /\ (SUC (LOG2 n)) ** 2 <= ordz k n ==>
1729         SQRT (phi k) * SUC (LOG2 n) <= phi k``,
1730  rpt strip_tac >>
1731  `0 < SUC (LOG2 n)` by decide_tac >>
1732  metis_tac[ZN_order_condition_property_0]);
1733
1734(* Theorem: 1 < n /\ coprime k n /\ (ulog n) ** 2 <= ordz k n ==>
1735            SQRT (phi k) * (ulog n) <= phi k *)
1736(* Proof:
1737   Note 0 < ulog n                         by ulog_pos, 1 < n
1738   Thus SQRT (phi k) * (ulog n) <= phi k   by ZN_order_condition_property_0
1739*)
1740val ZN_order_condition_property_2 = store_thm(
1741  "ZN_order_condition_property_2",
1742  ``!n k. 1 < n /\ coprime k n /\ (ulog n) ** 2 <= ordz k n ==>
1743         SQRT (phi k) * (ulog n) <= phi k``,
1744  metis_tac[ZN_order_condition_property_0, ulog_pos]);
1745
1746(* Theorem: 1 < n /\ coprime k n /\ (2 * ulog n) ** 2 <= ordz k n ==>
1747            2 * SQRT (phi k) * (ulog n) <= phi k *)
1748(* Proof:
1749   Note 0 < ulog n                             by ulog_pos, 1 < n
1750     so 0 < 2 * ulog n                         by arithmetic
1751   Thus SQRT (phi k) * (2 * ulog n) <= phi k   by ZN_order_condition_property_0
1752     or   2 * SQRT (phi k) * ulog n <= phi k
1753*)
1754val ZN_order_condition_property_3 = store_thm(
1755  "ZN_order_condition_property_3",
1756  ``!n k. 1 < n /\ coprime k n /\ (2 * ulog n) ** 2 <= ordz k n ==>
1757             2 * SQRT (phi k) * ulog n <= phi k``,
1758  rpt strip_tac >>
1759  `0 < ulog n` by rw[ulog_pos] >>
1760  `0 < 2 * ulog n` by decide_tac >>
1761  `2 * SQRT (phi k) * ulog n = SQRT (phi k) * (2 * ulog n)` by decide_tac >>
1762  metis_tac[ZN_order_condition_property_0]);
1763
1764(* ------------------------------------------------------------------------- *)
1765(* AKS Parameter Search                                                      *)
1766(* ------------------------------------------------------------------------- *)
1767
1768(* Possible outcomes of search loop *)
1769val _ = Datatype `param_search_result =
1770       nice num    (* found an nice of type :num telling prime or composite *)
1771     | good num    (* found a good of type :num good for AKS *)
1772     | bad         (* when search exceeds the upper bound *)
1773`;
1774
1775(* Note: pseudocode to search for parameter k.
1776
1777Input: n, m
1778Output: 1 < k that either (gcd k n <> 1) or (m <= ordz k n)
1779
1780k <- 2
1781while k < (2 + (ulog n) * HALF (m ** 2)) {
1782    if (gcd k n <> 1) exit
1783    if (m <= ordz k n) exit
1784    k <- k + 1
1785}
1786
1787This version using divides:
1788k <- 2
1789while k < (2 + (ulog n) * HALF (m ** 2)) {
1790    if (k divides n) exit
1791    if (m <= ordz k n) exit
1792    k <- k + 1
1793}
1794The following theorem guarantees that either exit must occur within the while-loop.
1795
1796ZN_order_modulus_exists_1
1797|- !n m. 1 < n /\ 1 < m ==>
1798   ?k. 1 < k /\ k < 2 + ulog n * HALF (m ** 2) /\ (gcd k n <> 1 \/ m <= ordz k n)
1799
1800Usually m = (ulog n) ** 2, and m = 1 iff (ulog n = 1), iff n = 2.
1801In this case, there is no suitable parameter k,
1802so just make (gcd k n <> 1), e.g. take k = n.
1803
1804Discussion
1805----------
1806The test for order: if (m <= ordz k n) exit
1807can be modified to: if (m <= k) and (m <= ordz k n) exit
1808that is, skipping computation of (ordz k n) when k < m.
1809
1810This is due to ZN_order_lt |- !k n m. 0 < k /\ k < m ==> ordz k n < m
1811As the proof shows,
1812     (ordz k n <= phi k)   by ZN_order_upper, (ordz k n) divides (phi k)
1813and     (phi k <  k)       by phi_lt, 1 < k
1814 so   ordz k n < k
1815and if k < m, then ordz k n < m, no hope to have (m <= ordz k n).
1816
1817This modification is used in the next version, with two loops.
1818*)
1819
1820(* Pseudo-Code:
1821   AKSparam n =
1822       m <- SQ (ulog n)               // the order to exceed
1823       c <- 2 + HALF ((ulog n) ** 5)  // the upper limit
1824       k <- 2                         // search from 2 onwards
1825       loop:
1826          if (m <= k) goto next                  // to stage two
1827          if (n MOD k == 0) return 'nice' k      // a factor k
1828          k <- k + 1
1829          goto loop
1830       next:
1831          if (c < k) return 'bad'                // not found
1832          if (n MOD k == 0) return 'nice' k      // a factor k
1833          // Here k must be coprime to n, by no factor from 2 to k
1834          if (m <= order(k, n)) return 'good' k  // suitable k
1835          k <- k + 1
1836          goto next
1837
1838   'bad': flag <- F, k <- 0
1839   'nice': flag <- F, k
1840   'good': flag <- T, k
1841*)
1842
1843(* Define the search for AKS parameter: given n and m, returns k, with cutoff at c. *)
1844val aks_param_search_def = tDefine "aks_param_search"  `
1845    aks_param_search n m k c =
1846      if c < k then bad (* unreachable *)
1847      else if k divides n then nice k    (* if (k divides n) exit *)
1848      else if m <= k /\ m <= ordz_compute k n then good k  (* if (m <= ordz k n) exit *)
1849      else aks_param_search n m (k + 1) c  (* k <- k + 1 *)
1850`(WF_REL_TAC `measure (\(n, m, k, c). c + 1 - k)`);
1851(* Check discussion in pseudo-code above for the good check: m <= k /\ m <= ordz_compute k n *)
1852
1853(* Define the AKS parameter: given n, return k, for better estimate:
1854   In theory: c = 1 + (ulog n) * (HALF (SQ m))
1855   In practice:       (ulog n) * (HALF (SQ m)) <= HALF ((ulog n) * (SQ m))
1856                                                = HALF ((ulog n) ** 5)
1857   The last one is easier for the estimate: O((ulog n)^5)
1858   HALF_MULT  |- !m n. n * HALF m <= HALF (n * m)
1859*)
1860(*  Originally aks_param_search has
1861              if m <= 1 then nice n : special since HALF (m ** 2) will be 0
1862    SQ (ulog n) <= 1, ulog 0 = 0, ulog 1 = 0, ulog 2 = 1, ulog 3 = 2, n <= 2.
1863    The check is moved to here.
1864*)
1865val aks_param_def = Define `
1866    aks_param n = if n <= 2 then nice n
1867                  else aks_param_search n (SQ (ulog n)) 2 (1 + HALF ((ulog n) ** 5))
1868`;
1869
1870(*
1871> EVAL ``aks_param 31``;
1872val it = |- aks_param 31 = good 29: thm
1873> EVAL ``aks_param 17``;
1874val it = |- aks_param 17 = nice 17: thm
1875> EVAL ``aks_param 95477``;
1876val it = |- aks_param 95477 = good 293: thm
1877> > EVAL ``MAP aks_param [0; 1; 2]``;
1878val it = |- MAP aks_param [0; 1; 2] = [nice 0; nice 1; nice 2]: thm
1879*)
1880
1881(* Theorem: aks_param n =
1882            let m = ulog n;
1883                c = 1 + HALF (m ** 5)
1884             in if m <= 1 then nice n else aks_param_search n (SQ m) 2 c *)
1885(* Proof: by aks_param_def, ulog_le_1 *)
1886val aks_param_alt = store_thm(
1887  "aks_param_alt",
1888  ``!n. aks_param n =
1889       let m = ulog n;
1890           c = 1 + HALF (m ** 5)
1891        in if m <= 1 then nice n else aks_param_search n (SQ m) 2 c``,
1892  rw[aks_param_def, ulog_le_1]);
1893
1894(* Theorem: aks_param_search n m k c =
1895            if c < k then bad
1896            else if k divides n then nice k
1897            else if m <= k /\ m <= ordz k n then good k
1898            else aks_param_search n m (k + 1) c *)
1899(* Proof: by aks_param_search_def, ordz_compute_eqn *)
1900val aks_param_search_alt = store_thm(
1901  "aks_param_search_alt",
1902  ``!n m k c. aks_param_search n m k c =
1903             if c < k then bad
1904             else if k divides n then nice k
1905             else if m <= k /\ m <= ordz k n then good k
1906             else aks_param_search n m (k + 1) c``,
1907  rw[Once aks_param_search_def, ordz_compute_eqn]);
1908
1909(* Theorem: aks_param 0 = nice 0 *)
1910(* Proof: by aks_param_def, 0 <= 2. *)
1911val aks_param_0 = store_thm(
1912  "aks_param_0",
1913  ``aks_param 0 = nice 0``,
1914  rw[aks_param_def]);
1915
1916(* Theorem: aks_param 1 = nice 1 *)
1917(* Proof: by aks_param_def, 1 <= 2. *)
1918val aks_param_1 = store_thm(
1919  "aks_param_1",
1920  ``aks_param 1 = nice 1``,
1921  rw[aks_param_def]);
1922
1923(* Theorem: aks_param 2 = nice 2 *)
1924(* Proof: by aks_param_def, 2 <= 2. *)
1925val aks_param_2 = store_thm(
1926  "aks_param_2",
1927  ``aks_param 2 = nice 2``,
1928  rw[aks_param_def]);
1929
1930(* Theorem: 0 < k /\ (aks_param_search n m k c = bad) ==>
1931            !j. k <= j /\ j <= c ==> ~(j divides n) /\ (ordz j n < m) *)
1932(* Proof:
1933   By aks_param_search_alt, and aks_param_search_ind induction, this is to show:
1934   (1) k <= j /\ j <= c /\ ~(k divides n) ==> ~(j divides n)
1935       If k = j, this is trivially true.
1936       If k <> j, then k < j, true by induction hypothesis, k + 1 <= j
1937   (2) k <= j /\ j <= c /\ ~(m <= k /\ m <= ordz k n) ==> ~(m <= ordz j n)
1938       If k = j,
1939          If m <= k, this is trivially true.
1940          If ~(m <= k), then k < m,
1941             By contradiction, suppose m <= ordz k n.
1942             Then ordz k n < k         by ZN_order_lt, 0 < k, k < m
1943             Hence there is a contradiction.
1944       If k <> j, then k < j, true by induction hypothesis, k + 1 <= j
1945*)
1946val aks_param_search_bad = store_thm(
1947  "aks_param_search_bad",
1948  ``!n m k c. 0 < k /\ (aks_param_search n m k c = bad) ==>
1949             !j. k <= j /\ j <= c ==> ~(j divides n) /\ (ordz j n < m)``,
1950  ho_match_mp_tac (theorem "aks_param_search_ind") >>
1951  (rpt gen_tac >> strip_tac) >>
1952  simp[Once aks_param_search_alt] >>
1953  rw[] >| [
1954    Cases_on `k = j` >-
1955    rw[] >>
1956    fs[ordz_compute_eqn],
1957    Cases_on `k = j` >| [
1958      Cases_on `m <= k` >-
1959      rw[] >>
1960      spose_not_then strip_assume_tac >>
1961      `ordz k n < m` by fs[ZN_order_lt] >>
1962      rw[],
1963      fs[ordz_compute_eqn]
1964    ]
1965  ]);
1966
1967(* Theorem: k <= j /\ j <= c /\ (j divides n \/ m <= ordz j n) ==>
1968            (k = 0) \/ (aks_param_search n m k c <> bad) *)
1969(* Proof: by aks_param_search_bad *)
1970val aks_param_search_not_bad = store_thm(
1971  "aks_param_search_not_bad",
1972  ``!n m k j c. k <= j /\ j <= c /\ (j divides n \/ m <= ordz j n) ==>
1973    (k = 0) \/ (aks_param_search n m k c <> bad)``,
1974  metis_tac[aks_param_search_bad, NOT_ZERO_LT_ZERO, NOT_LESS_EQUAL]);
1975
1976(* Theorem: 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * (ulog n) <= c ==> aks_param_search n m 2 c <> bad *)
1977(* Proof:
1978   Note ?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)
1979                                          by ZN_order_modulus_exists_3, 1 < m
1980    now 1 < k ==> 2 <= k,
1981    and k divides n ==> nice k <> bad     by aks_param_search_not_bad
1982     or m <= ordz k n ==> good k <> bad   by aks_param_search_not_bad
1983*)
1984val aks_param_search_success = store_thm(
1985  "aks_param_search_success",
1986  ``!n m c. 1 < n /\ 1 < m /\ 1 + HALF (m ** 2) * (ulog n) <= c ==>
1987           aks_param_search n m 2 c <> bad``,
1988  rpt strip_tac >>
1989  `?k. 1 < k /\ k <= c /\ (k divides n \/ m <= ordz k n)` by rw[ZN_order_modulus_exists_3] >| [
1990    `2 <= k /\ 2 <> 0` by decide_tac >>
1991    metis_tac[aks_param_search_not_bad],
1992    `2 <= k /\ 2 <> 0` by decide_tac >>
1993    metis_tac[aks_param_search_not_bad]
1994  ]);
1995
1996(* Theorem: (aks_param_search n m k c = bad) /\ d <= c ==> (aks_param_search n m k d = bad) *)
1997(* Proof: by aks_param_search_def and its induction *)
1998val aks_param_search_bad_le = store_thm(
1999  "aks_param_search_bad_le",
2000  ``!n m k c d. (aks_param_search n m k c = bad) /\ d <= c ==>
2001                (aks_param_search n m k d = bad)``,
2002  ho_match_mp_tac (theorem "aks_param_search_ind") >>
2003  (rpt gen_tac >> strip_tac) >>
2004  simp[Once aks_param_search_def] >>
2005  rw[] >-
2006  rw[Once aks_param_search_def] >>
2007  rw[Once aks_param_search_def]);
2008
2009(* Theorem: aks_param n <> bad *)
2010(* Proof:
2011   By aks_param_def, this is to show:
2012      aks_param_search n (SQ (ulog n)) 2 (1 + HALF (ulog n ** 5)) <> bad
2013   Let m = SQ (ulog n), k = 1 + HALF (m ** 2) * ulog n, c = 1 + HALF (ulog n ** 5).
2014
2015   If n <= 2, aks_param n = nice <> bad       by aks_param_def
2016   Otherwise 2 < n, so 1 < m                  by ulog_sq_gt_1
2017
2018   Claim: k <= c
2019   Proof: k = 1 + HALF ((SQ (ulog n)) ** 2) * ulog n     by notation
2020            = 1 + HALF (((ulog n) ** 2) ** 2) * ulog n   by EXP_2
2021            = 1 + ulog n * HALF ((ulog n) ** 4)          by EXP_EXP_MULT, MULT_COMM
2022            <= 1 + HALF (ulog n * (ulog n) ** 4)         by HALF_MULT
2023             = 1 + HALF ((ulog n) ** 5)                  by EXP
2024             = c                                         by notation
2025
2026   The result follows             by aks_param_search_success
2027*)
2028val aks_param_exists = store_thm(
2029  "aks_param_exists",
2030  ``!n. aks_param n <> bad``,
2031  rw_tac std_ss[aks_param_def] >>
2032  qabbrev_tac `m = SQ (ulog n)` >>
2033  qabbrev_tac `k = 1 + HALF (m ** 2) * ulog n` >>
2034  qabbrev_tac `c = 1 + HALF (ulog n ** 5)` >>
2035  `1 < n /\ 2 < n` by decide_tac >>
2036  `1 < m` by fs[ulog_sq_gt_1, Abbr`m`] >>
2037  `k <= c` by
2038  (rw[Abbr`k`, Abbr`c`, Abbr`m`] >>
2039  `(ulog n ** 2) ** 2 = (ulog n) ** 4` by rw[GSYM EXP_EXP_MULT] >>
2040  metis_tac[HALF_MULT, EXP, DECIDE``SUC 4 = 5``]) >>
2041  metis_tac[aks_param_search_success]);
2042
2043(* Theorem: (aks_param_search n m k_ c = nice k) ==>
2044            (h <= k /\ k divides n /\ (!j. h <= j /\ j < k ==> ~(j divides n))) *)
2045(* Proof:
2046   By aks_param_search_def, and aks_param_search_ind induction, this is to show:
2047   (1) k divides n ==> k divides n, true trivially.
2048   (2) ~(c < h) /\ ~(h divides n) /\ ~(m <= h /\ m <= ordz_compute h n) /\
2049       aks_param_search n m (h + 1) c = nice k /\ ~(h divides n) ==> h <= k
2050       Note !j. k + 1 <= j /\ j < k ==> ~(j divides n)    by induction hypothesis
2051       Take j = k + 1, then h + 1 <= k + 1, or h <= k.
2052   (3) h <= j /\ j < k /\ ~(h divides n) ==> ~(j divides n)
2053       If h = j, this is trivially true.
2054       If h <> j, then h < j, true by induction hypothesis, h + 1 <= j
2055*)
2056val aks_param_search_nice = store_thm(
2057  "aks_param_search_nice",
2058  ``!n m h c k. (aks_param_search n m h c = nice k) ==>
2059        (h <= k /\ k divides n /\ (!j. h <= j /\ j < k ==> ~(j divides n)))``,
2060  ho_match_mp_tac (theorem "aks_param_search_ind") >>
2061  (rpt gen_tac >> strip_tac) >>
2062  simp[Once aks_param_search_def] >>
2063  rw[] >-
2064  rw[] >-
2065  fs[] >>
2066  (Cases_on `h = j` >> fs[]));
2067
2068(* Theorem: 0 < n /\ (aks_param_search n m 2 c = nice k) ==> k <= n *)
2069(* Proof:
2070   By contradiction, suppose ~(k <= n), or n < k.
2071   Then k divides n      by aks_param_search_nice
2072   But k divides n ==> k <= n      by DIVIDES_LE, 0 < n
2073   This is a contradiction.
2074*)
2075val aks_param_search_nice_le = store_thm(
2076  "aks_param_search_nice_le",
2077  ``!n m k c. 0 < n /\ (aks_param_search n m 2 c = nice k) ==> k <= n``,
2078  metis_tac[aks_param_search_nice, DIVIDES_LE]);
2079
2080(* Theorem: (aks_param_search n m k c = nice j) ==> j <= c *)
2081(* Proof:
2082   By induction from aks_param_search_def, noting cases.
2083   The only case to give (nice j) is:
2084       ~(c < k) /\ k divides n
2085       giving (nice k), so j = k, and k <= c        by ~(c < k)
2086   Otherwise,
2087           aks_param_search n m k c = nice j
2088       <=> aks_param_search n m (k + 1) c = nice j  by aks_param_search_def
2089       ==> j <= c                                   by induction hypothesis
2090*)
2091val aks_param_search_nice_bound = store_thm(
2092  "aks_param_search_nice_bound",
2093  ``!n m k c j. (aks_param_search n m k c = nice j) ==> j <= c``,
2094  ho_match_mp_tac (theorem "aks_param_search_ind") >>
2095  rw[] >>
2096  Cases_on `c < k` >-
2097  fs[Once aks_param_search_def] >>
2098  Cases_on `k divides n` >-
2099  fs[Once aks_param_search_def] >>
2100  Cases_on `m <= k /\ m <= ordz_compute k n` >-
2101  fs[Once aks_param_search_def] >>
2102  fs[Once aks_param_search_def]);
2103
2104(* Theorem: (aks_param_search n m k c = nice j) ==> j divides n *)
2105(* Proof: by aks_param_search_nice *)
2106val aks_param_search_nice_factor = store_thm(
2107  "aks_param_search_nice_factor",
2108  ``!n m k c j. (aks_param_search n m k c = nice j) ==> j divides n``,
2109  metis_tac[aks_param_search_nice, NOT_LESS]);
2110
2111(* Theorem: (aks_param_search n m 2 c = nice k) ==> !j. 1 < j /\ j < k ==> coprime j n *)
2112(* Proof:
2113   By contradiction, suppose ?j. 1 < j /\ j < k, but ~(coprime j n).
2114   Note 2 <= k /\ k divides n /\
2115        (!j. 2 <= j /\ j < k ==> ~(j divides n))    by aks_param_search_nice
2116
2117   Let d = gcd j n, then d <> 1.
2118   Note j <> 0                       by 1 < j
2119     so d <> 0                       by GCD_EQ_0
2120   Thus 2 <= d                       by d <> 0, d <> 1
2121    Now d divides j /\ d divides n   by GCD_PROPERTY
2122    ==> d <= j                       by DIVIDES_LE, 0 < j
2123     or d < k                        by LESS_EQ_LESS_TRANS
2124    ==> ~(d divides n)               by implication
2125    This contradicts d divides n.
2126*)
2127val aks_param_search_nice_coprime_all = store_thm(
2128  "aks_param_search_nice_coprime_all",
2129  ``!n m k c. (aks_param_search n m 2 c = nice k) ==> !j. 1 < j /\ j < k ==> coprime j n``,
2130  spose_not_then strip_assume_tac >>
2131  `2 <= k /\ k divides n /\ !j. 2 <= j /\ j < k ==> ~(j divides n)` by metis_tac[aks_param_search_nice] >>
2132  `2 <= j` by decide_tac >>
2133  qabbrev_tac `d = gcd j n` >>
2134  `j <> 0 /\ 0 < j` by decide_tac >>
2135  `d <> 0` by metis_tac[GCD_EQ_0] >>
2136  `2 <= d` by decide_tac >>
2137  `d divides j /\ d divides n` by metis_tac[GCD_PROPERTY] >>
2138  `d <= j` by rw[DIVIDES_LE] >>
2139  metis_tac[LESS_EQ_LESS_TRANS]);
2140
2141(* Theorem: (aks_param_search n m h c = good k) ==>
2142            h <= k /\ m <= ordz k n /\ !j. h <= j /\ j <= k ==> ~(j divides n) *)
2143(* Proof:
2144   By aks_param_search_alt, and aks_param_search_ind induction, this is to show:
2145   (1) h <= j /\ j <= h /\ ~(h divides n) ==> ~(j divides n)
2146       Note j = h, hence true                   by h <= j /\ j <= h
2147   (2) ~(c < h) /\ ~(m <= h /\ m <= ordz h n) /\ (aks_param_search n m (h + 1) c = good k) /\ ~(h divides n) ==> h <= k
2148       Note !j. h + 1 <= j /\ m <= ordz j n    by induction hypothesis, ordz_compute_eqn
2149       Take j = k + 1, then h + 1 <= k + 1, or h <= k.
2150   (3) ~(c < h) /\ ~(m <= h /\ m <= ordz h n) /\ (aks_param_search n m (h + 1) c = good k) /\ ~(h divides n) ==> m <= ordz k n, true trivially.
2151       True                                     by induction hypothesis, ordz_compute_eqn
2152   (4) h <= j /\ j <= k /\ ~(h divides n) ==> ~(j divides n)
2153       If h = j, this is trivially true.
2154       If h <> j, then h < j, true by induction hypothesis, h + 1 <= j, ordz_compute_eqn
2155*)
2156val aks_param_search_good = store_thm(
2157  "aks_param_search_good",
2158  ``!n m h c k. (aks_param_search n m h c = good k) ==>
2159    h <= k /\ m <= ordz k n /\ !j. h <= j /\ j <= k ==> ~(j divides n)``,
2160  ho_match_mp_tac (theorem "aks_param_search_ind") >>
2161  (rpt gen_tac >> strip_tac) >>
2162  simp[Once aks_param_search_alt] >>
2163  rw[] >| [
2164    `j = h` by decide_tac >>
2165    rw[],
2166    fs[ordz_compute_eqn],
2167    fs[ordz_compute_eqn],
2168    (Cases_on `h = j` >> fs[ordz_compute_eqn])
2169  ]);
2170
2171(* Theorem: 1 < n /\ (aks_param_search n m 2 c = good k) ==> k < n *)
2172(* Proof:
2173   By contradiction, suppose ~(k < n), or n <= k.
2174   Note n <> 1 /\ 2 <= n                          by 1 < n
2175    and !j. 2 <= j /\ j <= k ==> ~(j divides n)   by aks_param_search_good
2176   Putting j = n, this gives  ~(n divides n)      by implication
2177   This contradicts n divides n                   by DIVIDES_REFL
2178*)
2179val aks_param_search_good_lt = store_thm(
2180  "aks_param_search_good_lt",
2181  ``!n m k c. 1 < n /\ (aks_param_search n m 2 c = good k) ==> k < n``,
2182  spose_not_then strip_assume_tac >>
2183  `n <> 1 /\ 2 <= n /\ n <= k` by decide_tac >>
2184  metis_tac[aks_param_search_good, DIVIDES_REFL]);
2185
2186(* Theorem: (aks_param_search n m k c = good j) ==> j <= c *)
2187(* Proof:
2188   By induction from aks_param_search_def, noting cases.
2189   The only case to give (good j) is:
2190       ~(c < k) /\ ~(k divides n) /\ m <= k /\ m <= ordz_compute k n
2191       giving (good k), so j = k, and k <= c        by ~(c < k)
2192   Otherwise,
2193           aks_param_search n m k c = good j
2194       <=> aks_param_search n m (k + 1) c = good j  by aks_param_search_def
2195       ==> j <= c                                   by induction hypothesis
2196*)
2197val aks_param_search_good_bound = store_thm(
2198  "aks_param_search_good_bound",
2199  ``!n m k c j. (aks_param_search n m k c = good j) ==> j <= c``,
2200  ho_match_mp_tac (theorem "aks_param_search_ind") >>
2201  rw[] >>
2202  Cases_on `c < k` >-
2203  fs[Once aks_param_search_def] >>
2204  Cases_on `k divides n` >-
2205  fs[Once aks_param_search_def] >>
2206  Cases_on `m <= k /\ m <= ordz_compute k n` >-
2207  fs[Once aks_param_search_def] >>
2208  fs[Once aks_param_search_def]);
2209
2210(* Theorem: (aks_param_search n m 2 c = good k) ==>
2211            1 < k /\ m <= ordz k n /\ !j. 1 < j /\ j <= k ==> coprime j n *)
2212(* Proof
2213   Note 2 <= k /\ m <= ordz k n /\ !j. 2 <= j /\ j <= k ==> ~(j divides n)  by aks_param_search_good
2214   Thus to show:
2215   (1) 1 < k, true                       by 2 <= k.
2216   (2) m <= ordz k n, true               trivially.
2217   (3) 1 < j /\ j <= k ==> coprime j n
2218       Let d = gcd j n, then d <> 1.
2219       Note j <> 0                       by 1 < j
2220         so d <> 0                       by GCD_EQ_0
2221       Thus 2 <= d                       by d <> 0, d <> 1
2222        Now d divides j /\ d divides n   by GCD_PROPERTY
2223        ==> d <= j                       by DIVIDES_LE, 0 < j
2224         or d <= k                       by LESS_EQ_TRANS
2225        ==> ~(d divides n)               by implication
2226       This contradicts d divides n.
2227*)
2228val aks_param_search_good_coprime_all = store_thm(
2229  "aks_param_search_good_coprime_all",
2230  ``!n m k c. (aks_param_search n m 2 c = good k) ==>
2231    1 < k /\ m <= ordz k n /\ !j. 1 < j /\ j <= k ==> coprime j n``,
2232  rpt strip_tac >-
2233  metis_tac[aks_param_search_good, DECIDE``!n. 1 < n <=> 2 <= n``] >-
2234  metis_tac[aks_param_search_good] >>
2235  spose_not_then strip_assume_tac >>
2236  `2 <= k /\ !j. 2 <= j /\ j <= k ==> ~(j divides n)` by metis_tac[aks_param_search_good] >>
2237  qabbrev_tac `d = gcd j n` >>
2238  `j <> 0 /\ 0 < j` by decide_tac >>
2239  `d <> 0` by metis_tac[GCD_EQ_0] >>
2240  `2 <= d` by decide_tac >>
2241  `d divides j /\ d divides n` by metis_tac[GCD_PROPERTY] >>
2242  `d <= j` by rw[DIVIDES_LE] >>
2243  metis_tac[LESS_EQ_TRANS]);
2244
2245(* Theorem: 1 < n /\ (aks_param n = nice k) ==> 1 < k /\ k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n) *)
2246(* Proof:
2247   If n = 2,
2248      Then aks_param 2 = nice 2     by aks_param_2
2249        or k = 2                    by nice k = nice 2
2250        so 2 divides 2              by DIVIDES_REFL
2251       and !j. 1 < j /\ j < 2 ==> ~(j divides n), as no j will satisfy the premise.
2252   If n <> 2,
2253      Then 2 < n                    by n <> 2, 1 < n
2254      Let m = SQ (ulog n), c = 1 + HALF ((ulog n) ** 5).
2255      Note aks_param_search n m 2 c
2256         = aks_param n              by aks_param_def, 2 < n
2257         = nice k                   by given
2258       ==> 2 <= k /\ k divides n /\ !j. 2 <= j /\ j < k ==> ~(j divides n)
2259                                    by aks_param_search_nice
2260        or 1 < k /\ k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n)
2261                                    by arithmetic
2262*)
2263val aks_param_nice = store_thm(
2264  "aks_param_nice",
2265  ``!n k. 1 < n /\ (aks_param n = nice k) ==>
2266         1 < k /\ k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n)``,
2267  ntac 3 strip_tac >>
2268  Cases_on `n = 2` >| [
2269    `nice 2 = nice k` by metis_tac[aks_param_2] >>
2270    `k = 2` by rw[] >>
2271    rw[],
2272    qabbrev_tac `m = SQ (ulog n)` >>
2273    qabbrev_tac `c = 1 + HALF ((ulog n) ** 5)` >>
2274    `~(n <= 2)` by decide_tac >>
2275    `aks_param_search n m 2 c = nice k` by metis_tac[aks_param_def] >>
2276    `2 <= k /\ k divides n /\ !j. 2 <= j /\ j < k ==> ~(j divides n)` by metis_tac[aks_param_search_nice] >>
2277    rw[]
2278  ]);
2279
2280(* Theorem: (aks_param n = nice k) ==>
2281            (if n <= 1 then (k = n) else 1 < k) /\
2282            k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n) *)
2283(* Proof:
2284   If n = 0,
2285      Then aks_param 0 = nice 0     by aks_param_0
2286       and 0 divides 0, hence true.
2287   If n = 1,
2288      Then aks_param 1 = nice 1     by aks_param_1
2289       and 1 divides 1, hence true.
2290   Otherwise 1 < n, true            by aks_param_nice
2291*)
2292val aks_param_nice_alt = store_thm(
2293  "aks_param_nice_alt",
2294  ``!n k. (aks_param n = nice k) ==>
2295         (if n <= 1 then (k = n) else 1 < k) /\
2296         k divides n /\ !j. 1 < j /\ j < k ==> ~(j divides n)``,
2297  ntac 3 strip_tac >>
2298  (Cases_on `n <= 1` >> simp[]) >| [
2299    `(n = 0) \/ (n = 1)` by decide_tac >-
2300    fs[aks_param_0] >>
2301    fs[aks_param_1],
2302    metis_tac[aks_param_nice, NOT_LESS_EQUAL]
2303  ]);
2304
2305(* Theorem: (aks_param n = nice k) ==> k <= n *)
2306(* Proof:
2307   Let m = SQ (ulog n), c = 2 + HALF ((ulog n) ** 5).
2308   If n <= 2,
2309      aks_param n = nice n       by aks_param_def
2310      so k = n, k <= n is true trivially.
2311   Otherwise 2 < n,
2312   Note aks_param_search n m 2 c
2313      = aks_param n              by aks_param_def, 2 < n
2314      = nice k                   by given
2315   Thus k <= n                   by aks_param_search_nice_le, 0 < n
2316*)
2317val aks_param_nice_le = store_thm(
2318  "aks_param_nice_le",
2319  ``!n k. (aks_param n = nice k) ==> k <= n``,
2320  rw[aks_param_def] >>
2321  `0 < n` by decide_tac >>
2322  metis_tac[aks_param_search_nice_le]);
2323
2324(* Theorem: 2 < n /\ (aks_param n = nice k) ==> k <= 1 + HALF (ulog n ** 5) *)
2325(* Proof:
2326   Let c = 1 + HALF (ulog n ** 5).
2327   Note aks_param n = nice k
2328    <=> aks_param_search n (SQ (ulog n)) 2 c = nice k    by aks_param_def, 2 < n
2329    ==> k <= c                                           by aks_param_search_nice_bound
2330*)
2331val aks_param_nice_bound = store_thm(
2332  "aks_param_nice_bound",
2333  ``!n k. 2 < n /\ (aks_param n = nice k) ==> k <= 1 + HALF (ulog n ** 5)``,
2334  metis_tac[aks_param_def, aks_param_search_nice_bound, NOT_LESS]);
2335
2336(* Theorem: (aks_param n = nice k) ==> !j. 1 < j /\ j < k ==> coprime j n *)
2337(* Proof:
2338   If n <= 2,
2339      Then aks_param n = nice n        by aks_param_def
2340      so k = n,
2341      and n = 0, or n = 1, or n = 2    by n <= 2
2342      and there is no j such that 1 < j /\ j < k.
2343   If ~(n <= 2),
2344      Then 2 < n                        by n <> 2, 1 < n
2345      Let m = SQ (ulog n), c = 1 + HALF ((ulog n) ** 5).
2346      Then aks_param n = nice k
2347       <=> aks_param_search n m 2 c = nice k
2348                                        by aks_param_def, 2 < n
2349      The result follows                by aks_param_search_nice_coprime_all
2350*)
2351val aks_param_nice_coprime_all = store_thm(
2352  "aks_param_nice_coprime_all",
2353  ``!n k. (aks_param n = nice k) ==> !j. 1 < j /\ j < k ==> coprime j n``,
2354  rw_tac std_ss[aks_param_def] >-
2355  decide_tac >>
2356  metis_tac[aks_param_search_nice_coprime_all]);
2357
2358(* Theorem: 1 < n /\ (aks_param n = nice k) ==> (prime n <=> (k = n)) *)
2359(* Proof:
2360   If part: prime n ==> k = n
2361      Note 1 < k /\ k divides n              by aks_param_nice
2362        or n <> 1 /\ k <> 1                  by 1 < n, 1 < k
2363      Thus k = n                             by prime_def, k divides n, prime n.
2364   Only-if part: k = n ==> prime n
2365      Note !j. 1 < j /\ j < n ==> ~(j divides n)  by aks_param_nice, k = n
2366       ==> prime n                                by prime_iff_no_proper_factor
2367*)
2368val aks_param_nice_for_prime = store_thm(
2369  "aks_param_nice_for_prime",
2370  ``!n k. 1 < n /\ (aks_param n = nice k) ==> (prime n <=> (k = n))``,
2371  rw_tac std_ss[EQ_IMP_THM] >-
2372  metis_tac[aks_param_nice, prime_def, LESS_NOT_EQ] >>
2373  rw[aks_param_nice, prime_iff_no_proper_factor]);
2374
2375(* Theorem: (aks_param n = good k) ==>
2376            1 < k /\ SQ (ulog n) <= ordz k n /\ !j. 1 < j /\ j <= k ==> ~(j divides n) *)
2377(* Proof:
2378   If n <= 2, aks_param n = nice k <> good k.
2379   Otherwise,
2380   Let m = SQ (ulog n), c = 1 + HALF ((ulog n) ** 5).
2381   Note aks_param_search n m 2 c
2382      = aks_param n              by aks_param_def, 2 < n
2383      = good k                   by given
2384   Thus 2 <= k /\ m <= ordz k n /\ !j. 2 <= j /\ j <= k ==> ~(j divides n)
2385                                 by aks_param_search_good
2386     or 1 < k /\ m <= ordz k n /\ !j. 1 < j /\ j <= k ==> ~(j divides n)
2387                                 by arithmetic
2388*)
2389val aks_param_good = store_thm(
2390  "aks_param_good",
2391  ``!n k. (aks_param n = good k) ==>
2392         1 < k /\ SQ (ulog n) <= ordz k n /\ !j. 1 < j /\ j <= k ==> ~(j divides n)``,
2393  ntac 3 strip_tac >>
2394  qabbrev_tac `m = SQ (ulog n)` >>
2395  qabbrev_tac `c = 1 + HALF ((ulog n) ** 5)` >>
2396  Cases_on `n <= 2` >-
2397  fs[aks_param_def] >>
2398  `aks_param_search n m 2 c = good k` by metis_tac[aks_param_def] >>
2399  `2 <= k /\ m <= ordz k n /\ !j. 2 <= j /\ j <= k ==> ~(j divides n)` by metis_tac[aks_param_search_good] >>
2400  rw[]);
2401
2402(* Theorem: (aks_param n = good k) ==> k < n *)
2403(* Proof: by aks_param_def, aks_param_search_good_lt *)
2404val aks_param_good_lt = store_thm(
2405  "aks_param_good_lt",
2406  ``!n k. (aks_param n = good k) ==> k < n``,
2407  rw[aks_param_def] >>
2408  `1 < n` by decide_tac >>
2409  metis_tac[aks_param_search_good_lt]);
2410
2411(* Theorem: (aks_param n = good k) ==> k <= 1 + HALF (ulog n ** 5) *)
2412(* Proof:
2413   If n <= 2, aks_param n = nice n <> good k   by aks_param_def
2414   Otherwise,
2415   Let c = 1 + HALF (ulog n ** 5).
2416       aks_param n = good k
2417   <=> aks_param_search n (SQ (ulog n)) 2 c = good k    by aks_param_def, 2 < n
2418   ==> k <= c                                           by aks_param_search_good_bound
2419*)
2420val aks_param_good_bound = store_thm(
2421  "aks_param_good_bound",
2422  ``!n k. (aks_param n = good k) ==> k <= 1 + HALF (ulog n ** 5)``,
2423  rw_tac std_ss[aks_param_def] >>
2424  metis_tac[aks_param_search_good_bound]);
2425
2426(* Theorem: (aks_param n = good k) ==> 1 < k /\ SQ (ulog n) <= ordz k n /\ !j. 1 < j /\ j <= k ==> coprime j n *)
2427(* Proof: by aks_param_def, aks_param_search_good_coprime_all *)
2428val aks_param_good_coprime_all = store_thm(
2429  "aks_param_good_coprime_all",
2430  ``!n k. (aks_param n = good k) ==> !j. 1 < j /\ j <= k ==> coprime j n``,
2431  rw_tac std_ss[aks_param_def] >>
2432  metis_tac[aks_param_search_good_coprime_all]);
2433
2434(* Theorem: (aks_param n = good k) ==> 1 < n /\ 1 < k /\ k < n *)
2435(* Proof:
2436   Claim: 1 < n
2437   Proof: By contradiction, suppose ~(1 < n).
2438          Then n = 0 \/ n = 1                  by arithmetic
2439           But aks_param n = nice n <> good k  by aks_param_0, aks_param_1
2440
2441   Thus 1 < n                           by Claim
2442     so 1 < k                           by aks_param_good
2443    and k < n                           by aks_param_good_lt, 1 < n
2444*)
2445val aks_param_good_range = store_thm(
2446  "aks_param_good_range",
2447  ``!n k. (aks_param n = good k) ==> 1 < n /\ 1 < k /\ k < n``,
2448  ntac 3 strip_tac >>
2449  `1 < n` by
2450  (spose_not_then strip_assume_tac >>
2451  `aks_param n = nice n` by metis_tac[aks_param_0, aks_param_1, DECIDE``~(1 < n) <=> (n = 0) \/ (n = 1)``] >>
2452  `!n k. good n <> nice k` by rw[] >>
2453  metis_tac[]) >>
2454  metis_tac[aks_param_good, aks_param_good_lt]);
2455
2456(* Theorem: (aks_param n = good k) ==> 1 < k /\ k < n /\ coprime k n /\ (ulog n) ** 2 <= ordz k n *)
2457(* Proof: by aks_param_good, aks_param_good_range, aks_param_good_coprime_all *)
2458val aks_param_good_coprime_order = store_thm(
2459  "aks_param_good_coprime_order",
2460  ``!n k. (aks_param n = good k) ==> 1 < k /\ k < n /\ coprime k n /\ (ulog n) ** 2 <= ordz k n``,
2461  metis_tac[aks_param_good_range, aks_param_good, EXP_2, aks_param_good_coprime_all, LESS_EQ_REFL]);
2462
2463(* ------------------------------------------------------------------------- *)
2464(* AKS parameter search simplified                                           *)
2465(* ------------------------------------------------------------------------- *)
2466(*
2467val aks_param_search_def = tDefine "aks_param_search"  `
2468    aks_param_search n m k c =
2469           if c < k then bad (* unreachable *)
2470      else if k divides n then nice k    (* if (k divides n) exit *)
2471      else if m <= k /\ m <= ordz_compute k n then good k  (* if (m <= ordz k n) exit *)
2472      else aks_param_search n m (k + 1) c  (* k <- k + 1 *)
2473`(WF_REL_TAC `measure (\(n, m, k, c). c + 1 - k)`);
2474val aks_param_def = Define `
2475    aks_param n = if n <= 2 then nice n
2476                  else aks_param_search n (SQ (ulog n)) 2 (1 + HALF ((ulog n) ** 5))
2477`;
2478*)
2479(* Note: the simplified version is better for analysis:
2480   (1) Have k = 0 check, so that k divides n is replaced by (n MOD k = 0).
2481   (2) Have c <= k, which replaces the calling c by incrementing 1.
2482   (3) Skip m <= k, which is an optimisation by ZN_order_le: ordz k n <= k when 0 < k.
2483*)
2484(* Define the parameter seek loop *)
2485(*
2486val param_seek_def = tDefine "param_seek" `
2487  param_seek m c n k =
2488       if k = 0 then bad
2489  else if c <= k then bad
2490  else if n MOD k = 0 then nice k (* same as k divides n when k <> 0 *)
2491  else if m <= ordz k n then good k
2492  else param_seek m c n (k + 1)
2493`(WF_REL_TAC `measure (\(m,c,n,k). c - k)`);
2494*)
2495(* Skip k = 0 check, as caller uses k = 2 *)
2496val param_seek_def = tDefine "param_seek" `
2497  param_seek m c n k =
2498       if c <= k then bad
2499  else if n MOD k = 0 then nice k (* same as k divides n when k <> 0 *)
2500  else if m <= ordz k n then good k
2501  else param_seek m c n (k + 1)
2502`(WF_REL_TAC `measure (\(m,c,n,k). c - k)`);
2503
2504(* Define the caller to parameter seek loop *)
2505val param_def = Define`
2506    param n = if n <= 2 then nice n
2507                  else param_seek (SQ (ulog n)) (2 + HALF ((ulog n) ** 5)) n 2
2508`;
2509(* Note: 2 + HALF ((ulog n) ** 5) is one more than aks_param *)
2510
2511(*
2512> EVAL ``param 31``;
2513val it = |- param 31 = good 29: thm
2514*)
2515
2516(* Theorem: param n = let m = ulog n;
2517                             c = 2 + HALF (m ** 5)
2518                          in if m <= 1 then nice n else param_seek (SQ m) c n 2 *)
2519(* Proof: by param_def, ulog_le_1 *)
2520val param_alt = store_thm(
2521  "param_alt",
2522  ``!n. param n =
2523          let m = ulog n;
2524              c = 2 + HALF (m ** 5)
2525           in if m <= 1 then nice n else param_seek (SQ m) c n 2``,
2526  rw[param_def, ulog_le_1]);
2527
2528(* Theorem: 0 < k /\ k <= c ==> (param_seek m c n k = aks_param_search n m k c) *)
2529(* Proof:
2530   Note k divides n <=> (n MOD k = 0)    by DIVIDES_MOD_0, 0 < k.
2531   By induction from param_seek_def, to show successive cases:
2532   (1) c <= k, means c - 1 < k           by 0 < c
2533         param_seek m c n k
2534       = bad                             by param_seek_def
2535       = aks_param_search n m k (c - 1)  by aks_param_search_def, c - 1 < k.
2536   (2) n MOD k = 0
2537       Note k divides n                  by DIVIDES_MOD_0, 0 < k
2538         param_seek m c n k
2539       = nice k                          by param_seek_def
2540       = aks_param_search n m k (c - 1)  by aks_param_search_def, k divides n
2541   (3) m <= ordz k n
2542       Note ~(k divides n)               by DIVIDES_MOD_0, 0 < k
2543       Note ordz k n <= k                by ZN_order_le, 0 < k
2544         so m <= k                       by LESS_EQ_TRANS
2545         param_seek m c n k
2546       = good k                          by param_seek_def
2547       = aks_param_search n m k (c - 1)  by aks_param_search_alt, m <= k /\ m <= ordz k n
2548   (4) Otherwise,
2549       Note k <= c - 1                   by ~(c <= k), or k < c
2550         or ~(c - 1 < k)
2551         param_seek m c n k
2552       = param_seek m c n (k + 1)              by param_seek_def
2553       = aks_param_search n m (k + 1) (c - 1)  by induction hypothesis
2554       = aks_param_search n m k (c - 1)        by aks_param_search_alt, ~(c - 1 < k)
2555*)
2556val param_seek_thm = store_thm(
2557  "param_seek_thm",
2558  ``!m c n k. 0 < k /\ k <= c ==> (param_seek m c n k = aks_param_search n m k (c - 1))``,
2559  ho_match_mp_tac (theorem "param_seek_ind") >>
2560  rw[] >>
2561  `k divides n <=> (n MOD k = 0)` by rw[DIVIDES_MOD_0] >>
2562  rw[Once param_seek_def] >-
2563  rw[Once aks_param_search_def] >-
2564  rw[Once aks_param_search_alt] >-
2565 (`ordz k n <= k` by rw[ZN_order_le] >>
2566  `m <= k` by decide_tac >>
2567  rw[Once aks_param_search_alt]) >>
2568  `~(c - 1 < k)` by decide_tac >>
2569  metis_tac[aks_param_search_alt]);
2570
2571(* Theorem: param n = aks_param n *)
2572(* Proof:
2573   If n <= 2, param n = nice n = aks_param n               by param_def, aks_param_def
2574   Otherwise,
2575     param n
2576   = param_seek (SQ (ulog n)) (2 + HALF (ulog n ** 5)) n 2    by param_def
2577   = aks_param_search n m 2 (1 + HALF (ulog n ** 5)           by param_seek_thm, 0 < 2
2578   = aks_param n                                              by aks_param_def
2579*)
2580val param_eqn = store_thm(
2581  "param_eqn",
2582  ``!n. param n = aks_param n``,
2583  rw[param_def, param_seek_thm, aks_param_def]);
2584
2585(* Obtain theorems *)
2586
2587val param_0 = save_thm("param_0", aks_param_0 |> SIMP_RULE bool_ss [GSYM param_eqn]);
2588(* val param_0 = |- param 0 = nice 0: thm *)
2589
2590val param_1 = save_thm("param_1", aks_param_1 |> SIMP_RULE bool_ss [GSYM param_eqn]);
2591(* val param_1 = |- param 1 = nice 1: thm *)
2592
2593val param_2 = save_thm("param_2", aks_param_2 |> SIMP_RULE bool_ss [GSYM param_eqn]);
2594(* val param_2 = |- param 2 = nice 2: thm *)
2595
2596val param_nice_bound = save_thm("param_nice_bound",
2597    aks_param_nice_bound |> SIMP_RULE bool_ss [GSYM param_eqn]);
2598(*
2599val param_nice_bound =
2600   |- !n k. 2 < n /\ (param n = nice k) ==> k <= 1 + HALF (ulog n ** 5): thm
2601*)
2602
2603val param_good_bound = save_thm("param_good_bound",
2604    aks_param_good_bound |> SIMP_RULE bool_ss [GSYM param_eqn]);
2605(*
2606val param_good_bound =
2607   |- !n k. (param n = good k) ==> k <= 1 + HALF (ulog n ** 5): thm
2608*)
2609
2610val param_exists = save_thm("param_exists",
2611    aks_param_exists |> SIMP_RULE bool_ss [GSYM param_eqn]);
2612(*
2613val param_exists = |- !n. param n <> bad: thm
2614*)
2615
2616val param_nice_for_prime = save_thm("param_nice_for_prime",
2617    aks_param_nice_for_prime |> SIMP_RULE bool_ss [GSYM param_eqn]);
2618(*
2619val param_nice_for_prime = |- !n k. 1 < n /\ (param n = nice k) ==> (prime n <=> (k = n)): thm
2620*)
2621
2622(* The following is absent:
2623val param_good_for_prime = save_thm("param_good_for_prime",
2624    aks_param_good_for_prime |> SIMP_RULE bool_ss [GSYM param_eqn]);
2625giving
2626val param_good_for_prime =
2627   |- !n k. (param n = good k) /\ power_free n ==>
2628          (prime n <=> poly_intro_checks n k (SQRT (phi k) * ulog n)): thm
2629
2630This is because:
2631(1) aks_param_good_for_prime is in AKSclean, due to poly_intro_checks.
2632(2) (param n) corresponds to poly_intro_checks n k k, not the range indicated.
2633*)
2634
2635val param_good_range = save_thm("param_good_range",
2636    aks_param_good_range |> SIMP_RULE bool_ss [GSYM param_eqn]);
2637(*
2638val param_good_range =
2639   |- !n k. (param n = good k) ==> 1 < n /\ 1 < k /\ k < n: thm
2640*)
2641
2642
2643(* ------------------------------------------------------------------------- *)
2644
2645(* export theory at end *)
2646val _ = export_theory();
2647
2648(*===========================================================================*)
2649