51173
|
1 |
(* Title: HOL/Number_Theory/Eratosthenes.thy
|
|
2 |
Author: Florian Haftmann, TU Muenchen
|
|
3 |
*)
|
|
4 |
|
|
5 |
header {* The sieve of Eratosthenes *}
|
|
6 |
|
|
7 |
theory Eratosthenes
|
|
8 |
imports Primes
|
|
9 |
begin
|
|
10 |
|
|
11 |
subsection {* Preliminary: strict divisibility *}
|
|
12 |
|
|
13 |
context dvd
|
|
14 |
begin
|
|
15 |
|
|
16 |
abbreviation dvd_strict :: "'a \<Rightarrow> 'a \<Rightarrow> bool" (infixl "dvd'_strict" 50)
|
|
17 |
where
|
|
18 |
"b dvd_strict a \<equiv> b dvd a \<and> \<not> a dvd b"
|
|
19 |
|
|
20 |
end
|
|
21 |
|
|
22 |
subsection {* Main corpus *}
|
|
23 |
|
|
24 |
text {* The sieve is modelled as a list of booleans, where @{const False} means \emph{marked out}. *}
|
|
25 |
|
|
26 |
type_synonym marks = "bool list"
|
|
27 |
|
|
28 |
definition numbers_of_marks :: "nat \<Rightarrow> marks \<Rightarrow> nat set"
|
|
29 |
where
|
|
30 |
"numbers_of_marks n bs = fst ` {x \<in> set (enumerate n bs). snd x}"
|
|
31 |
|
|
32 |
lemma numbers_of_marks_simps [simp, code]:
|
|
33 |
"numbers_of_marks n [] = {}"
|
|
34 |
"numbers_of_marks n (True # bs) = insert n (numbers_of_marks (Suc n) bs)"
|
|
35 |
"numbers_of_marks n (False # bs) = numbers_of_marks (Suc n) bs"
|
|
36 |
by (auto simp add: numbers_of_marks_def intro!: image_eqI)
|
|
37 |
|
|
38 |
lemma numbers_of_marks_Suc:
|
|
39 |
"numbers_of_marks (Suc n) bs = Suc ` numbers_of_marks n bs"
|
|
40 |
by (auto simp add: numbers_of_marks_def enumerate_Suc_eq image_iff Bex_def)
|
|
41 |
|
|
42 |
lemma numbers_of_marks_replicate_False [simp]:
|
|
43 |
"numbers_of_marks n (replicate m False) = {}"
|
|
44 |
by (auto simp add: numbers_of_marks_def enumerate_replicate_eq)
|
|
45 |
|
|
46 |
lemma numbers_of_marks_replicate_True [simp]:
|
|
47 |
"numbers_of_marks n (replicate m True) = {n..<n+m}"
|
|
48 |
by (auto simp add: numbers_of_marks_def enumerate_replicate_eq image_def)
|
|
49 |
|
|
50 |
lemma in_numbers_of_marks_eq:
|
|
51 |
"m \<in> numbers_of_marks n bs \<longleftrightarrow> m \<in> {n..<n + length bs} \<and> bs ! (m - n)"
|
|
52 |
by (simp add: numbers_of_marks_def in_set_enumerate_eq image_iff add_commute)
|
|
53 |
|
|
54 |
|
|
55 |
text {* Marking out multiples in a sieve *}
|
|
56 |
|
|
57 |
definition mark_out :: "nat \<Rightarrow> marks \<Rightarrow> marks"
|
|
58 |
where
|
|
59 |
"mark_out n bs = map (\<lambda>(q, b). b \<and> \<not> Suc n dvd Suc (Suc q)) (enumerate n bs)"
|
|
60 |
|
|
61 |
lemma mark_out_Nil [simp]:
|
|
62 |
"mark_out n [] = []"
|
|
63 |
by (simp add: mark_out_def)
|
|
64 |
|
|
65 |
lemma length_mark_out [simp]:
|
|
66 |
"length (mark_out n bs) = length bs"
|
|
67 |
by (simp add: mark_out_def)
|
|
68 |
|
|
69 |
lemma numbers_of_marks_mark_out:
|
|
70 |
"numbers_of_marks n (mark_out m bs) = {q \<in> numbers_of_marks n bs. \<not> Suc m dvd Suc q - n}"
|
|
71 |
by (auto simp add: numbers_of_marks_def mark_out_def in_set_enumerate_eq image_iff
|
|
72 |
nth_enumerate_eq less_dvd_minus)
|
|
73 |
|
|
74 |
|
|
75 |
text {* Auxiliary operation for efficient implementation *}
|
|
76 |
|
|
77 |
definition mark_out_aux :: "nat \<Rightarrow> nat \<Rightarrow> marks \<Rightarrow> marks"
|
|
78 |
where
|
|
79 |
"mark_out_aux n m bs =
|
|
80 |
map (\<lambda>(q, b). b \<and> (q < m + n \<or> \<not> Suc n dvd Suc (Suc q) + (n - m mod Suc n))) (enumerate n bs)"
|
|
81 |
|
|
82 |
lemma mark_out_code [code]:
|
|
83 |
"mark_out n bs = mark_out_aux n n bs"
|
|
84 |
proof -
|
|
85 |
{ fix a
|
|
86 |
assume A: "Suc n dvd Suc (Suc a)"
|
|
87 |
and B: "a < n + n"
|
|
88 |
and C: "n \<le> a"
|
|
89 |
have False
|
|
90 |
proof (cases "n = 0")
|
|
91 |
case True with A B C show False by simp
|
|
92 |
next
|
|
93 |
def m \<equiv> "Suc n" then have "m > 0" by simp
|
|
94 |
case False then have "n > 0" by simp
|
|
95 |
from A obtain q where q: "Suc (Suc a) = Suc n * q" by (rule dvdE)
|
|
96 |
have "q > 0"
|
|
97 |
proof (rule ccontr)
|
|
98 |
assume "\<not> q > 0"
|
|
99 |
with q show False by simp
|
|
100 |
qed
|
|
101 |
with `n > 0` have "Suc n * q \<ge> 2" by (auto simp add: gr0_conv_Suc)
|
|
102 |
with q have a: "a = Suc n * q - 2" by simp
|
|
103 |
with B have "q + n * q < n + n + 2"
|
|
104 |
by auto
|
|
105 |
then have "m * q < m * 2" by (simp add: m_def)
|
|
106 |
with `m > 0` have "q < 2" by simp
|
|
107 |
with `q > 0` have "q = 1" by simp
|
|
108 |
with a have "a = n - 1" by simp
|
|
109 |
with `n > 0` C show False by simp
|
|
110 |
qed
|
|
111 |
} note aux = this
|
|
112 |
show ?thesis
|
|
113 |
by (auto simp add: mark_out_def mark_out_aux_def in_set_enumerate_eq intro: aux)
|
|
114 |
qed
|
|
115 |
|
|
116 |
lemma mark_out_aux_simps [simp, code]:
|
|
117 |
"mark_out_aux n m [] = []" (is ?thesis1)
|
|
118 |
"mark_out_aux n 0 (b # bs) = False # mark_out_aux n n bs" (is ?thesis2)
|
|
119 |
"mark_out_aux n (Suc m) (b # bs) = b # mark_out_aux n m bs" (is ?thesis3)
|
|
120 |
proof -
|
|
121 |
show ?thesis1
|
|
122 |
by (simp add: mark_out_aux_def)
|
|
123 |
show ?thesis2
|
|
124 |
by (auto simp add: mark_out_code [symmetric] mark_out_aux_def mark_out_def
|
|
125 |
enumerate_Suc_eq in_set_enumerate_eq less_dvd_minus)
|
|
126 |
{ def v \<equiv> "Suc m" and w \<equiv> "Suc n"
|
|
127 |
fix q
|
|
128 |
assume "m + n \<le> q"
|
|
129 |
then obtain r where q: "q = m + n + r" by (auto simp add: le_iff_add)
|
|
130 |
{ fix u
|
|
131 |
from w_def have "u mod w < w" by simp
|
|
132 |
then have "u + (w - u mod w) = w + (u - u mod w)"
|
|
133 |
by simp
|
|
134 |
then have "u + (w - u mod w) = w + u div w * w"
|
|
135 |
by (simp add: div_mod_equality' [symmetric])
|
|
136 |
}
|
|
137 |
then have "w dvd v + w + r + (w - v mod w) \<longleftrightarrow> w dvd m + w + r + (w - m mod w)"
|
|
138 |
by (simp add: add_assoc add_left_commute [of m] add_left_commute [of v]
|
|
139 |
dvd_plus_eq_left dvd_plus_eq_right)
|
|
140 |
moreover from q have "Suc q = m + w + r" by (simp add: w_def)
|
|
141 |
moreover from q have "Suc (Suc q) = v + w + r" by (simp add: v_def w_def)
|
|
142 |
ultimately have "w dvd Suc (Suc (q + (w - v mod w))) \<longleftrightarrow> w dvd Suc (q + (w - m mod w))"
|
|
143 |
by (simp only: add_Suc [symmetric])
|
|
144 |
then have "Suc n dvd Suc (Suc (Suc (q + n) - Suc m mod Suc n)) \<longleftrightarrow>
|
|
145 |
Suc n dvd Suc (Suc (q + n - m mod Suc n))"
|
|
146 |
by (simp add: v_def w_def Suc_diff_le trans_le_add2)
|
|
147 |
}
|
|
148 |
then show ?thesis3
|
|
149 |
by (auto simp add: mark_out_aux_def
|
|
150 |
enumerate_Suc_eq in_set_enumerate_eq not_less)
|
|
151 |
qed
|
|
152 |
|
|
153 |
|
|
154 |
text {* Main entry point to sieve *}
|
|
155 |
|
|
156 |
fun sieve :: "nat \<Rightarrow> marks \<Rightarrow> marks"
|
|
157 |
where
|
|
158 |
"sieve n [] = []"
|
|
159 |
| "sieve n (False # bs) = False # sieve (Suc n) bs"
|
|
160 |
| "sieve n (True # bs) = True # sieve (Suc n) (mark_out n bs)"
|
|
161 |
|
|
162 |
text {*
|
|
163 |
There are the following possible optimisations here:
|
|
164 |
|
|
165 |
\begin{itemize}
|
|
166 |
|
|
167 |
\item @{const sieve} can abort as soon as @{term n} is too big to let
|
|
168 |
@{const mark_out} have any effect.
|
|
169 |
|
|
170 |
\item Search for further primes can be given up as soon as the search
|
|
171 |
position exceeds the square root of the maximum candidate.
|
|
172 |
|
|
173 |
\end{itemize}
|
|
174 |
|
|
175 |
This is left as an constructive exercise to the reader.
|
|
176 |
*}
|
|
177 |
|
|
178 |
lemma numbers_of_marks_sieve:
|
|
179 |
"numbers_of_marks (Suc n) (sieve n bs) =
|
|
180 |
{q \<in> numbers_of_marks (Suc n) bs. \<forall>m \<in> numbers_of_marks (Suc n) bs. \<not> m dvd_strict q}"
|
|
181 |
proof (induct n bs rule: sieve.induct)
|
|
182 |
case 1 show ?case by simp
|
|
183 |
next
|
|
184 |
case 2 then show ?case by simp
|
|
185 |
next
|
|
186 |
case (3 n bs)
|
|
187 |
have aux: "\<And>M n. n \<in> Suc ` M \<longleftrightarrow> n > 0 \<and> n - 1 \<in> M"
|
|
188 |
proof
|
|
189 |
fix M and n
|
|
190 |
assume "n \<in> Suc ` M" then show "n > 0 \<and> n - 1 \<in> M" by auto
|
|
191 |
next
|
|
192 |
fix M and n :: nat
|
|
193 |
assume "n > 0 \<and> n - 1 \<in> M"
|
|
194 |
then have "n > 0" and "n - 1 \<in> M" by auto
|
|
195 |
then have "Suc (n - 1) \<in> Suc ` M" by blast
|
|
196 |
with `n > 0` show "n \<in> Suc ` M" by simp
|
|
197 |
qed
|
|
198 |
{ fix m :: nat
|
|
199 |
assume "Suc (Suc n) \<le> m" and "m dvd Suc n"
|
|
200 |
from `m dvd Suc n` obtain q where "Suc n = m * q" ..
|
|
201 |
with `Suc (Suc n) \<le> m` have "Suc (m * q) \<le> m" by simp
|
|
202 |
then have "m * q < m" by arith
|
|
203 |
then have "q = 0" by simp
|
|
204 |
with `Suc n = m * q` have False by simp
|
|
205 |
} note aux1 = this
|
|
206 |
{ fix m q :: nat
|
|
207 |
assume "\<forall>q>0. 1 < q \<longrightarrow> Suc n < q \<longrightarrow> q \<le> Suc (n + length bs)
|
|
208 |
\<longrightarrow> bs ! (q - Suc (Suc n)) \<longrightarrow> \<not> Suc n dvd q \<longrightarrow> q dvd m \<longrightarrow> m dvd q"
|
|
209 |
then have *: "\<And>q. Suc n < q \<Longrightarrow> q \<le> Suc (n + length bs)
|
|
210 |
\<Longrightarrow> bs ! (q - Suc (Suc n)) \<Longrightarrow> \<not> Suc n dvd q \<Longrightarrow> q dvd m \<Longrightarrow> m dvd q"
|
|
211 |
by auto
|
|
212 |
assume "\<not> Suc n dvd m" and "q dvd m"
|
|
213 |
then have "\<not> Suc n dvd q" by (auto elim: dvdE)
|
|
214 |
moreover assume "Suc n < q" and "q \<le> Suc (n + length bs)"
|
|
215 |
and "bs ! (q - Suc (Suc n))"
|
|
216 |
moreover note `q dvd m`
|
|
217 |
ultimately have "m dvd q" by (auto intro: *)
|
|
218 |
} note aux2 = this
|
|
219 |
from 3 show ?case
|
|
220 |
apply (simp_all add: numbers_of_marks_mark_out numbers_of_marks_Suc Compr_image_eq inj_image_eq_iff
|
|
221 |
in_numbers_of_marks_eq Ball_def imp_conjL aux)
|
|
222 |
apply safe
|
|
223 |
apply (simp_all add: less_diff_conv2 le_diff_conv2 dvd_minus_self not_less)
|
|
224 |
apply (clarsimp dest!: aux1)
|
|
225 |
apply (simp add: Suc_le_eq less_Suc_eq_le)
|
|
226 |
apply (rule aux2) apply (clarsimp dest!: aux1)+
|
|
227 |
done
|
|
228 |
qed
|
|
229 |
|
|
230 |
|
|
231 |
text {* Relation the sieve algorithm to actual primes *}
|
|
232 |
|
|
233 |
definition primes_upto :: "nat \<Rightarrow> nat set"
|
|
234 |
where
|
|
235 |
"primes_upto n = {m. m \<le> n \<and> prime m}"
|
|
236 |
|
|
237 |
lemma in_primes_upto:
|
|
238 |
"m \<in> primes_upto n \<longleftrightarrow> m \<le> n \<and> prime m"
|
|
239 |
by (simp add: primes_upto_def)
|
|
240 |
|
|
241 |
lemma primes_upto_sieve [code]:
|
|
242 |
"primes_upto n = numbers_of_marks 2 (sieve 1 (replicate (n - 1) True))"
|
|
243 |
proof (cases "n > 1")
|
|
244 |
case False then have "n = 0 \<or> n = 1" by arith
|
|
245 |
then show ?thesis
|
|
246 |
by (auto simp add: numbers_of_marks_sieve One_nat_def numeral_2_eq_2 primes_upto_def dest: prime_gt_Suc_0_nat)
|
|
247 |
next
|
|
248 |
{ fix m q
|
|
249 |
assume "Suc (Suc 0) \<le> q"
|
|
250 |
and "q < Suc n"
|
|
251 |
and "m dvd q"
|
|
252 |
then have "m < Suc n" by (auto dest: dvd_imp_le)
|
|
253 |
assume *: "\<forall>m\<in>{Suc (Suc 0)..<Suc n}. m dvd q \<longrightarrow> q dvd m"
|
|
254 |
and "m dvd q" and "m \<noteq> 1"
|
|
255 |
have "m = q" proof (cases "m = 0")
|
|
256 |
case True with `m dvd q` show ?thesis by simp
|
|
257 |
next
|
|
258 |
case False with `m \<noteq> 1` have "Suc (Suc 0) \<le> m" by arith
|
|
259 |
with `m < Suc n` * `m dvd q` have "q dvd m" by simp
|
|
260 |
with `m dvd q` show ?thesis by (simp add: dvd.eq_iff)
|
|
261 |
qed
|
|
262 |
}
|
|
263 |
then have aux: "\<And>m q. Suc (Suc 0) \<le> q \<Longrightarrow>
|
|
264 |
q < Suc n \<Longrightarrow>
|
|
265 |
m dvd q \<Longrightarrow>
|
|
266 |
\<forall>m\<in>{Suc (Suc 0)..<Suc n}. m dvd q \<longrightarrow> q dvd m \<Longrightarrow>
|
|
267 |
m dvd q \<Longrightarrow> m \<noteq> q \<Longrightarrow> m = 1" by auto
|
|
268 |
case True then show ?thesis
|
|
269 |
apply (auto simp add: numbers_of_marks_sieve One_nat_def numeral_2_eq_2 primes_upto_def dest: prime_gt_Suc_0_nat)
|
|
270 |
apply (simp add: prime_nat_def dvd_def)
|
|
271 |
apply (auto simp add: prime_nat_def aux)
|
|
272 |
done
|
|
273 |
qed
|
|
274 |
|
|
275 |
end
|
|
276 |
|