author  haftmann 
Sun, 09 Nov 2014 10:03:18 +0100  
changeset 58954  18750e86d5b8 
parent 58889  5b7a9633cfa8 
child 60526  fad653acf58f 
permissions  rwrr 
51173  1 
(* Title: HOL/Number_Theory/Eratosthenes.thy 
2 
Author: Florian Haftmann, TU Muenchen 

3 
*) 

4 

58889  5 
section {* The sieve of Eratosthenes *} 
51173  6 

7 
theory Eratosthenes 

52379  8 
imports Main Primes 
51173  9 
begin 
10 

52379  11 

51173  12 
subsection {* Preliminary: strict divisibility *} 
13 

14 
context dvd 

15 
begin 

16 

17 
abbreviation dvd_strict :: "'a \<Rightarrow> 'a \<Rightarrow> bool" (infixl "dvd'_strict" 50) 

18 
where 

19 
"b dvd_strict a \<equiv> b dvd a \<and> \<not> a dvd b" 

20 

21 
end 

22 

23 
subsection {* Main corpus *} 

24 

25 
text {* The sieve is modelled as a list of booleans, where @{const False} means \emph{marked out}. *} 

26 

27 
type_synonym marks = "bool list" 

28 

29 
definition numbers_of_marks :: "nat \<Rightarrow> marks \<Rightarrow> nat set" 

30 
where 

31 
"numbers_of_marks n bs = fst ` {x \<in> set (enumerate n bs). snd x}" 

32 

33 
lemma numbers_of_marks_simps [simp, code]: 

34 
"numbers_of_marks n [] = {}" 

35 
"numbers_of_marks n (True # bs) = insert n (numbers_of_marks (Suc n) bs)" 

36 
"numbers_of_marks n (False # bs) = numbers_of_marks (Suc n) bs" 

37 
by (auto simp add: numbers_of_marks_def intro!: image_eqI) 

38 

39 
lemma numbers_of_marks_Suc: 

40 
"numbers_of_marks (Suc n) bs = Suc ` numbers_of_marks n bs" 

41 
by (auto simp add: numbers_of_marks_def enumerate_Suc_eq image_iff Bex_def) 

42 

43 
lemma numbers_of_marks_replicate_False [simp]: 

44 
"numbers_of_marks n (replicate m False) = {}" 

45 
by (auto simp add: numbers_of_marks_def enumerate_replicate_eq) 

46 

47 
lemma numbers_of_marks_replicate_True [simp]: 

48 
"numbers_of_marks n (replicate m True) = {n..<n+m}" 

49 
by (auto simp add: numbers_of_marks_def enumerate_replicate_eq image_def) 

50 

51 
lemma in_numbers_of_marks_eq: 

52 
"m \<in> numbers_of_marks n bs \<longleftrightarrow> m \<in> {n..<n + length bs} \<and> bs ! (m  n)" 

57512
cc97b347b301
reduced name variants for assoc and commute on plus and mult
haftmann
parents:
55337
diff
changeset

53 
by (simp add: numbers_of_marks_def in_set_enumerate_eq image_iff add.commute) 
51173  54 

52379  55 
lemma sorted_list_of_set_numbers_of_marks: 
56 
"sorted_list_of_set (numbers_of_marks n bs) = map fst (filter snd (enumerate n bs))" 

57 
by (auto simp add: numbers_of_marks_def distinct_map 

58 
intro!: sorted_filter distinct_filter inj_onI sorted_distinct_set_unique) 

59 

51173  60 

61 
text {* Marking out multiples in a sieve *} 

62 

63 
definition mark_out :: "nat \<Rightarrow> marks \<Rightarrow> marks" 

64 
where 

65 
"mark_out n bs = map (\<lambda>(q, b). b \<and> \<not> Suc n dvd Suc (Suc q)) (enumerate n bs)" 

66 

67 
lemma mark_out_Nil [simp]: 

68 
"mark_out n [] = []" 

69 
by (simp add: mark_out_def) 

70 

71 
lemma length_mark_out [simp]: 

72 
"length (mark_out n bs) = length bs" 

73 
by (simp add: mark_out_def) 

74 

75 
lemma numbers_of_marks_mark_out: 

76 
"numbers_of_marks n (mark_out m bs) = {q \<in> numbers_of_marks n bs. \<not> Suc m dvd Suc q  n}" 

77 
by (auto simp add: numbers_of_marks_def mark_out_def in_set_enumerate_eq image_iff 

54222  78 
nth_enumerate_eq less_eq_dvd_minus) 
51173  79 

80 

81 
text {* Auxiliary operation for efficient implementation *} 

82 

83 
definition mark_out_aux :: "nat \<Rightarrow> nat \<Rightarrow> marks \<Rightarrow> marks" 

84 
where 

85 
"mark_out_aux n m bs = 

86 
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)" 

87 

88 
lemma mark_out_code [code]: 

89 
"mark_out n bs = mark_out_aux n n bs" 

90 
proof  

91 
{ fix a 

92 
assume A: "Suc n dvd Suc (Suc a)" 

93 
and B: "a < n + n" 

94 
and C: "n \<le> a" 

95 
have False 

96 
proof (cases "n = 0") 

97 
case True with A B C show False by simp 

98 
next 

99 
def m \<equiv> "Suc n" then have "m > 0" by simp 

100 
case False then have "n > 0" by simp 

101 
from A obtain q where q: "Suc (Suc a) = Suc n * q" by (rule dvdE) 

102 
have "q > 0" 

103 
proof (rule ccontr) 

104 
assume "\<not> q > 0" 

105 
with q show False by simp 

106 
qed 

107 
with `n > 0` have "Suc n * q \<ge> 2" by (auto simp add: gr0_conv_Suc) 

108 
with q have a: "a = Suc n * q  2" by simp 

109 
with B have "q + n * q < n + n + 2" 

110 
by auto 

111 
then have "m * q < m * 2" by (simp add: m_def) 

112 
with `m > 0` have "q < 2" by simp 

113 
with `q > 0` have "q = 1" by simp 

114 
with a have "a = n  1" by simp 

115 
with `n > 0` C show False by simp 

116 
qed 

117 
} note aux = this 

118 
show ?thesis 

119 
by (auto simp add: mark_out_def mark_out_aux_def in_set_enumerate_eq intro: aux) 

120 
qed 

121 

122 
lemma mark_out_aux_simps [simp, code]: 

123 
"mark_out_aux n m [] = []" (is ?thesis1) 

124 
"mark_out_aux n 0 (b # bs) = False # mark_out_aux n n bs" (is ?thesis2) 

125 
"mark_out_aux n (Suc m) (b # bs) = b # mark_out_aux n m bs" (is ?thesis3) 

126 
proof  

127 
show ?thesis1 

128 
by (simp add: mark_out_aux_def) 

129 
show ?thesis2 

130 
by (auto simp add: mark_out_code [symmetric] mark_out_aux_def mark_out_def 

54222  131 
enumerate_Suc_eq in_set_enumerate_eq less_eq_dvd_minus) 
51173  132 
{ def v \<equiv> "Suc m" and w \<equiv> "Suc n" 
133 
fix q 

134 
assume "m + n \<le> q" 

135 
then obtain r where q: "q = m + n + r" by (auto simp add: le_iff_add) 

136 
{ fix u 

137 
from w_def have "u mod w < w" by simp 

138 
then have "u + (w  u mod w) = w + (u  u mod w)" 

139 
by simp 

140 
then have "u + (w  u mod w) = w + u div w * w" 

141 
by (simp add: div_mod_equality' [symmetric]) 

142 
} 

143 
then have "w dvd v + w + r + (w  v mod w) \<longleftrightarrow> w dvd m + w + r + (w  m mod w)" 

57512
cc97b347b301
reduced name variants for assoc and commute on plus and mult
haftmann
parents:
55337
diff
changeset

144 
by (simp add: add.assoc add.left_commute [of m] add.left_commute [of v] 
58649
a62065b5e1e2
generalized and consolidated some theorems concerning divisibility
haftmann
parents:
57512
diff
changeset

145 
dvd_add_left_iff dvd_add_right_iff) 
51173  146 
moreover from q have "Suc q = m + w + r" by (simp add: w_def) 
147 
moreover from q have "Suc (Suc q) = v + w + r" by (simp add: v_def w_def) 

148 
ultimately have "w dvd Suc (Suc (q + (w  v mod w))) \<longleftrightarrow> w dvd Suc (q + (w  m mod w))" 

149 
by (simp only: add_Suc [symmetric]) 

150 
then have "Suc n dvd Suc (Suc (Suc (q + n)  Suc m mod Suc n)) \<longleftrightarrow> 

151 
Suc n dvd Suc (Suc (q + n  m mod Suc n))" 

152 
by (simp add: v_def w_def Suc_diff_le trans_le_add2) 

153 
} 

154 
then show ?thesis3 

155 
by (auto simp add: mark_out_aux_def 

156 
enumerate_Suc_eq in_set_enumerate_eq not_less) 

157 
qed 

158 

159 

160 
text {* Main entry point to sieve *} 

161 

162 
fun sieve :: "nat \<Rightarrow> marks \<Rightarrow> marks" 

163 
where 

164 
"sieve n [] = []" 

165 
 "sieve n (False # bs) = False # sieve (Suc n) bs" 

166 
 "sieve n (True # bs) = True # sieve (Suc n) (mark_out n bs)" 

167 

168 
text {* 

169 
There are the following possible optimisations here: 

170 

171 
\begin{itemize} 

172 

173 
\item @{const sieve} can abort as soon as @{term n} is too big to let 

174 
@{const mark_out} have any effect. 

175 

176 
\item Search for further primes can be given up as soon as the search 

177 
position exceeds the square root of the maximum candidate. 

178 

179 
\end{itemize} 

180 

181 
This is left as an constructive exercise to the reader. 

182 
*} 

183 

184 
lemma numbers_of_marks_sieve: 

185 
"numbers_of_marks (Suc n) (sieve n bs) = 

186 
{q \<in> numbers_of_marks (Suc n) bs. \<forall>m \<in> numbers_of_marks (Suc n) bs. \<not> m dvd_strict q}" 

187 
proof (induct n bs rule: sieve.induct) 

188 
case 1 show ?case by simp 

189 
next 

190 
case 2 then show ?case by simp 

191 
next 

192 
case (3 n bs) 

193 
have aux: "\<And>M n. n \<in> Suc ` M \<longleftrightarrow> n > 0 \<and> n  1 \<in> M" 

194 
proof 

195 
fix M and n 

196 
assume "n \<in> Suc ` M" then show "n > 0 \<and> n  1 \<in> M" by auto 

197 
next 

198 
fix M and n :: nat 

199 
assume "n > 0 \<and> n  1 \<in> M" 

200 
then have "n > 0" and "n  1 \<in> M" by auto 

201 
then have "Suc (n  1) \<in> Suc ` M" by blast 

202 
with `n > 0` show "n \<in> Suc ` M" by simp 

203 
qed 

204 
{ fix m :: nat 

205 
assume "Suc (Suc n) \<le> m" and "m dvd Suc n" 

206 
from `m dvd Suc n` obtain q where "Suc n = m * q" .. 

207 
with `Suc (Suc n) \<le> m` have "Suc (m * q) \<le> m" by simp 

208 
then have "m * q < m" by arith 

209 
then have "q = 0" by simp 

210 
with `Suc n = m * q` have False by simp 

211 
} note aux1 = this 

212 
{ fix m q :: nat 

213 
assume "\<forall>q>0. 1 < q \<longrightarrow> Suc n < q \<longrightarrow> q \<le> Suc (n + length bs) 

214 
\<longrightarrow> bs ! (q  Suc (Suc n)) \<longrightarrow> \<not> Suc n dvd q \<longrightarrow> q dvd m \<longrightarrow> m dvd q" 

215 
then have *: "\<And>q. Suc n < q \<Longrightarrow> q \<le> Suc (n + length bs) 

216 
\<Longrightarrow> bs ! (q  Suc (Suc n)) \<Longrightarrow> \<not> Suc n dvd q \<Longrightarrow> q dvd m \<Longrightarrow> m dvd q" 

217 
by auto 

218 
assume "\<not> Suc n dvd m" and "q dvd m" 

219 
then have "\<not> Suc n dvd q" by (auto elim: dvdE) 

220 
moreover assume "Suc n < q" and "q \<le> Suc (n + length bs)" 

221 
and "bs ! (q  Suc (Suc n))" 

222 
moreover note `q dvd m` 

223 
ultimately have "m dvd q" by (auto intro: *) 

224 
} note aux2 = this 

225 
from 3 show ?case 

226 
apply (simp_all add: numbers_of_marks_mark_out numbers_of_marks_Suc Compr_image_eq inj_image_eq_iff 

227 
in_numbers_of_marks_eq Ball_def imp_conjL aux) 

228 
apply safe 

229 
apply (simp_all add: less_diff_conv2 le_diff_conv2 dvd_minus_self not_less) 

230 
apply (clarsimp dest!: aux1) 

231 
apply (simp add: Suc_le_eq less_Suc_eq_le) 

232 
apply (rule aux2) apply (clarsimp dest!: aux1)+ 

233 
done 

234 
qed 

235 

236 

52379  237 
text {* Relation of the sieve algorithm to actual primes *} 
51173  238 

52379  239 
definition primes_upto :: "nat \<Rightarrow> nat list" 
51173  240 
where 
52379  241 
"primes_upto n = sorted_list_of_set {m. m \<le> n \<and> prime m}" 
51173  242 

52379  243 
lemma set_primes_upto: 
244 
"set (primes_upto n) = {m. m \<le> n \<and> prime m}" 

51173  245 
by (simp add: primes_upto_def) 
246 

52379  247 
lemma sorted_primes_upto [iff]: 
248 
"sorted (primes_upto n)" 

249 
by (simp add: primes_upto_def) 

250 

251 
lemma distinct_primes_upto [iff]: 

252 
"distinct (primes_upto n)" 

253 
by (simp add: primes_upto_def) 

254 

255 
lemma set_primes_upto_sieve: 

256 
"set (primes_upto n) = numbers_of_marks 2 (sieve 1 (replicate (n  1) True))" 

51173  257 
proof (cases "n > 1") 
258 
case False then have "n = 0 \<or> n = 1" by arith 

259 
then show ?thesis 

55130
70db8d380d62
Restored Suc rather than +1, and using Library/Binimial
paulson <lp15@cam.ac.uk>
parents:
54222
diff
changeset

260 
by (auto simp add: numbers_of_marks_sieve numeral_2_eq_2 set_primes_upto dest: prime_gt_Suc_0_nat) 
51173  261 
next 
262 
{ fix m q 

263 
assume "Suc (Suc 0) \<le> q" 

264 
and "q < Suc n" 

265 
and "m dvd q" 

266 
then have "m < Suc n" by (auto dest: dvd_imp_le) 

267 
assume *: "\<forall>m\<in>{Suc (Suc 0)..<Suc n}. m dvd q \<longrightarrow> q dvd m" 

268 
and "m dvd q" and "m \<noteq> 1" 

269 
have "m = q" proof (cases "m = 0") 

270 
case True with `m dvd q` show ?thesis by simp 

271 
next 

272 
case False with `m \<noteq> 1` have "Suc (Suc 0) \<le> m" by arith 

273 
with `m < Suc n` * `m dvd q` have "q dvd m" by simp 

274 
with `m dvd q` show ?thesis by (simp add: dvd.eq_iff) 

275 
qed 

276 
} 

277 
then have aux: "\<And>m q. Suc (Suc 0) \<le> q \<Longrightarrow> 

278 
q < Suc n \<Longrightarrow> 

279 
m dvd q \<Longrightarrow> 

280 
\<forall>m\<in>{Suc (Suc 0)..<Suc n}. m dvd q \<longrightarrow> q dvd m \<Longrightarrow> 

281 
m dvd q \<Longrightarrow> m \<noteq> q \<Longrightarrow> m = 1" by auto 

282 
case True then show ?thesis 

55337
5d45fb978d5a
Number_Theory no longer introduces One_nat_def as a simprule. Tidied some proofs.
paulson <lp15@cam.ac.uk>
parents:
55130
diff
changeset

283 
apply (auto simp add: One_nat_def numbers_of_marks_sieve numeral_2_eq_2 set_primes_upto 
5d45fb978d5a
Number_Theory no longer introduces One_nat_def as a simprule. Tidied some proofs.
paulson <lp15@cam.ac.uk>
parents:
55130
diff
changeset

284 
dest: prime_gt_Suc_0_nat) 
55130
70db8d380d62
Restored Suc rather than +1, and using Library/Binimial
paulson <lp15@cam.ac.uk>
parents:
54222
diff
changeset

285 
apply (metis One_nat_def Suc_le_eq less_not_refl prime_nat_def) 
70db8d380d62
Restored Suc rather than +1, and using Library/Binimial
paulson <lp15@cam.ac.uk>
parents:
54222
diff
changeset

286 
apply (metis One_nat_def Suc_le_eq aux prime_nat_def) 
51173  287 
done 
288 
qed 

289 

52379  290 
lemma primes_upto_sieve [code]: 
291 
"primes_upto n = map fst (filter snd (enumerate 2 (sieve 1 (replicate (n  1) True))))" 

292 
proof  

293 
have "primes_upto n = sorted_list_of_set (numbers_of_marks 2 (sieve 1 (replicate (n  1) True)))" 

294 
apply (rule sorted_distinct_set_unique) 

295 
apply (simp_all only: set_primes_upto_sieve numbers_of_marks_def) 

296 
apply auto 

297 
done 

298 
then show ?thesis by (simp add: sorted_list_of_set_numbers_of_marks) 

299 
qed 

300 

301 
lemma prime_in_primes_upto: 

302 
"prime n \<longleftrightarrow> n \<in> set (primes_upto n)" 

303 
by (simp add: set_primes_upto) 

304 

305 

306 
subsection {* Application: smallest prime beyond a certain number *} 

307 

308 
definition smallest_prime_beyond :: "nat \<Rightarrow> nat" 

309 
where 

310 
"smallest_prime_beyond n = (LEAST p. prime p \<and> p \<ge> n)" 

311 

312 
lemma 

313 
prime_smallest_prime_beyond [iff]: "prime (smallest_prime_beyond n)" (is ?P) 

314 
and smallest_prime_beyond_le [iff]: "smallest_prime_beyond n \<ge> n" (is ?Q) 

315 
proof  

316 
let ?least = "LEAST p. prime p \<and> p \<ge> n" 

317 
from primes_infinite obtain q where "prime q \<and> q \<ge> n" 

318 
by (metis finite_nat_set_iff_bounded_le mem_Collect_eq nat_le_linear) 

319 
then have "prime ?least \<and> ?least \<ge> n" by (rule LeastI) 

320 
then show ?P and ?Q by (simp_all add: smallest_prime_beyond_def) 

321 
qed 

322 

323 
lemma smallest_prime_beyond_smallest: 

324 
"prime p \<Longrightarrow> p \<ge> n \<Longrightarrow> smallest_prime_beyond n \<le> p" 

325 
by (simp only: smallest_prime_beyond_def) (auto intro: Least_le) 

326 

327 
lemma smallest_prime_beyond_eq: 

328 
"prime p \<Longrightarrow> p \<ge> n \<Longrightarrow> (\<And>q. prime q \<Longrightarrow> q \<ge> n \<Longrightarrow> q \<ge> p) \<Longrightarrow> smallest_prime_beyond n = p" 

329 
by (simp only: smallest_prime_beyond_def) (auto intro: Least_equality) 

330 

331 
definition smallest_prime_between :: "nat \<Rightarrow> nat \<Rightarrow> nat option" 

332 
where 

333 
"smallest_prime_between m n = 

334 
(if (\<exists>p. prime p \<and> m \<le> p \<and> p \<le> n) then Some (smallest_prime_beyond m) else None)" 

335 

336 
lemma smallest_prime_between_None: 

337 
"smallest_prime_between m n = None \<longleftrightarrow> (\<forall>q. m \<le> q \<and> q \<le> n \<longrightarrow> \<not> prime q)" 

338 
by (auto simp add: smallest_prime_between_def) 

339 

340 
lemma smallest_prime_betwen_Some: 

341 
"smallest_prime_between m n = Some p \<longleftrightarrow> smallest_prime_beyond m = p \<and> p \<le> n" 

342 
by (auto simp add: smallest_prime_between_def dest: smallest_prime_beyond_smallest [of _ m]) 

343 

344 
lemma [code]: 

345 
"smallest_prime_between m n = List.find (\<lambda>p. p \<ge> m) (primes_upto n)" 

346 
proof  

347 
{ fix p 

348 
def A \<equiv> "{p. p \<le> n \<and> prime p \<and> m \<le> p}" 

349 
assume assms: "m \<le> p" "prime p" "p \<le> n" 

350 
then have "smallest_prime_beyond m \<le> p" by (auto intro: smallest_prime_beyond_smallest) 

351 
from this `p \<le> n` have *: "smallest_prime_beyond m \<le> n" by (rule order_trans) 

352 
from assms have ex: "\<exists>p\<le>n. prime p \<and> m \<le> p" by auto 

353 
then have "finite A" by (auto simp add: A_def) 

354 
with * have "Min A = smallest_prime_beyond m" 

355 
by (auto simp add: A_def intro: Min_eqI smallest_prime_beyond_smallest) 

356 
with ex sorted_primes_upto have "List.find (\<lambda>p. p \<ge> m) (primes_upto n) = Some (smallest_prime_beyond m)" 

357 
by (auto simp add: set_primes_upto sorted_find_Min A_def) 

358 
} then show ?thesis 

359 
by (auto simp add: smallest_prime_between_def find_None_iff set_primes_upto intro!: sym [of _ None]) 

360 
qed 

361 

362 
definition smallest_prime_beyond_aux :: "nat \<Rightarrow> nat \<Rightarrow> nat" 

363 
where 

364 
"smallest_prime_beyond_aux k n = smallest_prime_beyond n" 

365 

366 
lemma [code]: 

367 
"smallest_prime_beyond_aux k n = 

368 
(case smallest_prime_between n (k * n) 

369 
of Some p \<Rightarrow> p 

370 
 None \<Rightarrow> smallest_prime_beyond_aux (Suc k) n)" 

371 
by (simp add: smallest_prime_beyond_aux_def smallest_prime_betwen_Some split: option.split) 

372 

373 
lemma [code]: 

374 
"smallest_prime_beyond n = smallest_prime_beyond_aux 2 n" 

375 
by (simp add: smallest_prime_beyond_aux_def) 

376 

51173  377 
end 
378 