Theory Sieve_Primes
theory Sieve_Primes
imports
"HOL-Computational_Algebra.Primes"
"../Num_Class"
"../HOLCF_Prelude"
begin
section ‹The Sieve of Eratosthenes›
declare [[coercion int]]
declare [[coercion_enabled]]
text ‹
This example proves that the well-known Haskell two-liner that lazily
calculates the list of all primes does indeed do so. This proof is using
coinduction.
›
text ‹
We need to hide some constants again since we imported something from HOL
not via @{theory "HOLCF-Prelude.HOLCF_Main"}.
›
no_notation
Rings.divide (infixl ‹div› 70) and
Rings.modulo (infixl ‹mod› 70)
no_notation
Set.member (‹(:)›) and
Set.member (‹(‹notation=‹infix :››_/ : _)› [51, 51] 50)
text ‹This is the implementation. We also need a modulus operator.›
fixrec sieve :: "[Integer] → [Integer]" where
"sieve⋅(p : xs) = p : (sieve⋅(filter⋅(Λ x. neg⋅(eq⋅(mod⋅x⋅p)⋅0))⋅xs))"
fixrec primes :: "[Integer]" where
"primes = sieve⋅[2..]"
text ‹Simplification rules for modI:›
definition MkI' :: "int ⇒ Integer" where
"MkI' x = MkI⋅x"
lemma MkI'_simps [simp]:
shows "MkI' 0 = 0" and "MkI' 1 = 1" and "MkI' (numeral k) = numeral k"
unfolding MkI'_def zero_Integer_def one_Integer_def numeral_Integer_eq
by rule+
lemma modI_numeral_numeral [simp]:
"mod⋅(numeral i)⋅(numeral j) = MkI' (Rings.modulo (numeral i) (numeral j))"
unfolding numeral_Integer_eq mod_Integer_def MkI'_def by simp
text ‹Some lemmas demonstrating evaluation of our list:›
lemma "primes !! 0 = 2"
unfolding primes.simps
apply (simp only: enumFrom_intsFrom_conv)
apply (subst intsFrom.simps)
apply simp
done
lemma "primes !! 1 = 3"
unfolding primes.simps
apply (simp only: enumFrom_intsFrom_conv)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply simp
done
lemma "primes !! 2 = 5"
unfolding primes.simps
apply (simp only: enumFrom_intsFrom_conv)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply simp
done
lemma "primes !! 3 = 7"
unfolding primes.simps
apply (simp only: enumFrom_intsFrom_conv)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (subst intsFrom.simps)
apply (simp del: filter_FF filter_TT)
done
text ‹Auxiliary lemmas about prime numbers›
lemma find_next_prime_nat:
fixes n :: nat
assumes "prime n"
shows "∃ n'. n' > n ∧ prime n' ∧ (∀k. n < k ⟶ k < n' ⟶ ¬ prime k)"
using ex_least_nat_le[of "λ k . k > n ∧ prime k"]
by (metis bigger_prime not_prime_0)
text ‹Simplification for andalso:›
lemma andAlso_Def[simp]: "((Def x) andalso (Def y)) = Def (x ∧ y)"
by (metis Def_bool2 Def_bool4 andalso_thms(1) andalso_thms(2))
text ‹This defines the bisimulation and proves it to be a list bisimulation.›
definition prim_bisim:
"prim_bisim x1 x2 = (∃ n . prime n ∧
x1 = sieve⋅(filter⋅(Λ (MkI⋅i). Def ((∀d. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅n..]) ∧
x2 = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅n..])"
lemma prim_bisim_is_bisim: "list_bisim prim_bisim"
proof -
{
fix xs ys
assume "prim_bisim xs ys"
then obtain n :: nat where
"prime n" and
"n > 1" and
"xs = sieve⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅n..])" (is "_ = sieve⋅?xs") and
"ys = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅n..]"
proof -
assume a1: "⋀n. ⟦prime n; 1 < n; xs = sieve⋅ (Data_List.filter⋅ (Λ (MkI⋅i). Def (∀d>1. d < int n ⟶ ¬ d dvd i))⋅ [MkI⋅(int n)..]); ys = Data_List.filter⋅ (Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅ [MkI⋅(int n)..]⟧ ⟹ thesis"
obtain ii where
f2: "∀is isa. (¬ prim_bisim is isa ∨ prime (ii is isa) ∧ is = sieve⋅ (Data_List.filter⋅ (Λ (MkI⋅i). Def (∀ia>1. ia < ii is isa ⟶ ¬ ia dvd i))⋅ [MkI⋅(ii is isa)..]) ∧ isa = Data_List.filter⋅ (Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅ [MkI⋅(ii is isa)..]) ∧ ((∀i. ¬ prime i ∨ is ≠ sieve⋅ (Data_List.filter⋅ (Λ (MkI⋅ia). Def (∀ib>1. ib < i ⟶ ¬ ib dvd ia))⋅ [MkI⋅i..]) ∨ isa ≠ Data_List.filter⋅ (Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅ [MkI⋅i..]) ∨ prim_bisim is isa)"
unfolding prim_bisim by (atomize_elim, (subst choice_iff[symmetric])+, blast)
then have f3: "prime (ii xs ys)"
using ‹prim_bisim xs ys› by blast
then obtain nn :: "int ⇒ nat" where
f4: "int (nn (ii xs ys)) = ii xs ys"
by (metis (no_types) prime_gt_0_int zero_less_imp_eq_int)
then have "prime (nn (ii xs ys))"
using f3 by (metis (no_types) prime_nat_int_transfer)
then show ?thesis
using f4 f2 a1 ‹prim_bisim xs ys› prime_gt_1_nat by presburger
qed
obtain n' where "n' > n" and "prime n'" and not_prime: "∀k. n < k ⟶ k < n' ⟶ ¬ prime k"
using find_next_prime_nat[OF ‹prime n›] by auto
{
fix k :: int
assume "n < k" and "k < n'"
have "k > 1" using ‹n < k› ‹n > 1› by auto
then obtain k' :: nat where "k = int k'" by (cases k) auto
then obtain p where "prime p" and "p dvd k"
using ‹k > 1› ‹k = int k'›
by (metis (full_types) less_numeral_extra(4) of_nat_1 of_nat_dvd_iff prime_factor_nat prime_nat_int_transfer)
then have "p < n'" using ‹k < n'› ‹k > 1›
using zdvd_imp_le [of p k] by simp
then have "p ≤ n" using ‹prime p› not_prime
using not_le prime_gt_0_int zero_less_imp_eq_int
by (metis of_nat_less_iff prime_nat_int_transfer)
then have "∃d::int>1. d ≤ n ∧ d dvd k" using ‹p dvd k› ‹prime p› of_nat_le_iff prime_gt_1_nat
prime_gt_1_int by auto
}
then have between_have_divisors: "⋀k::int. n < k ⟹ k < n' ⟹ ∃d::int>1. d ≤ n ∧ d dvd k".
{
fix i
{
assume small: "∀d::int>1. d ≤ n ⟶ ¬ d dvd i"
fix d
assume "1 < d" and "d dvd i" and "d < n'"
with small have "d > n" by auto
obtain d'::int where "d' > 1" and "d' ≤ n" and "d' dvd d" using between_have_divisors[OF ‹n < d› ‹d < n'›] by auto
with ‹d dvd i› small have False by (metis (full_types) dvd_trans)
}
then have "(∀d::int. d > 1 ⟶ d ≤ n ⟶ ¬ (d dvd i)) = (∀d::int. d > 1 ⟶ d < n' ⟶ ¬ (d dvd i))"
using ‹n' > n› by auto
}
then have between_not_relvant: "⋀ i. (∀d::int. d > 1 ⟶ d ≤ n ⟶ ¬ (d dvd i)) = (∀d::int. d > 1 ⟶ d < n' ⟶ ¬ (d dvd i))" .
from ‹prime n›
have "∀d::int >1. d < n ⟶ ¬ d dvd n"
unfolding prime_int_altdef using int_one_le_iff_zero_less le_less
by (simp add: prime_int_not_dvd)
then
obtain xs' where "?xs = (MkI⋅n) : xs'" and "xs' = filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..]"
by (subst (asm) intsFrom.simps[unfolded enumFrom_intsFrom_conv[symmetric]], simp add: one_Integer_def TT_def[symmetric] add.commute)
{
have "filter⋅(Λ x. neg⋅(eq⋅(mod⋅x⋅(MkI⋅n))⋅0))⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..])
= filter⋅(Λ (MkI⋅x). Def (¬ n dvd x))⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..])"
apply (rule filter_cong[rule_format])
apply (rename_tac x, case_tac x)
apply (simp)
apply (auto simp add: zero_Integer_def)
apply (rule FF_def)
apply (simp add: TT_def)
by (metis dvd_eq_mod_eq_0)
also have "... = filter⋅(Λ (MkI⋅i). Def ((¬ n dvd i) ∧ (∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..]"
by (auto intro!: filter_cong[rule_format])
also have "... = filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d ≤ n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..]"
apply (rule filter_cong[rule_format])
apply (rename_tac x, case_tac x)
using ‹n > 1›
apply auto
done
also have "... = filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d ≤ n ⟶ ¬ (d dvd i))))⋅[MkI⋅(int n+1)..]"
by (metis (no_types, lifting) of_nat_1 of_nat_add)
also have "... = filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d ≤ n ⟶ ¬ (d dvd i))))⋅[MkI⋅n'..]"
apply (rule filter_fast_forward[of n n']) using ‹n' > n›
apply (auto simp add: between_have_divisors)
done
also have "... = filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n' ⟶ ¬ (d dvd i))))⋅[MkI⋅n'..]"
by (auto intro: filter_cong[rule_format] simp add: between_not_relvant)
also note calculation
}
note tmp = this
{
have "xs = sieve⋅?xs" by fact
also have "... = sieve⋅((MkI⋅n) : xs')" using ‹?xs = _› by simp
also have "... = sieve⋅((MkI⋅n) : filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅(n+1)..])" using ‹xs' = _› by simp
also have "... = (MkI⋅n) : sieve⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n' ⟶ ¬ (d dvd i))))⋅[MkI⋅n'..])" using tmp by simp
also note calculation
}
moreover
{
have "ys = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅n..]" by fact
also have "... = (MkI⋅n) : filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅(int n+1)..]"
using ‹prime n›
by (subst intsFrom.simps[unfolded enumFrom_intsFrom_conv[symmetric]])(simp add: one_Integer_def TT_def[symmetric] add.commute)
also have "... = (MkI⋅n) : filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅n'..]"
apply (subst filter_fast_forward[of n n']) using ‹n' > n› and not_prime
apply auto
apply (metis (full_types) ‹⋀k. ⟦int n < k; k < int n'⟧ ⟹ ∃d>1. d ≤ int n ∧ d dvd k› le_less not_le prime_gt_0_int prime_int_not_dvd zdvd_imp_le)
done
also note calculation
}
moreover
note ‹prime n'›
ultimately
have "∃ p xs' ys'. xs = p : xs' ∧ ys = p : ys' ∧ prim_bisim xs' ys'" unfolding prim_bisim
using prime_nat_int_transfer by blast
}
then show ?thesis unfolding list.bisim_def by metis
qed
text ‹Now we apply coinduction:›
lemma sieve_produces_primes:
fixes n :: nat
assumes "prime n"
shows "sieve⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < n ⟶ ¬ (d dvd i))))⋅[MkI⋅n..])
= filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅n..]"
using assms
apply -
apply (rule list.coinduct[OF prim_bisim_is_bisim], auto simp add: prim_bisim)
using prime_nat_int_transfer by blast
text ‹And finally show the correctness of primes.›
theorem primes:
shows "primes = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅2..]"
proof -
have "primes = sieve⋅[2 ..]" by (rule primes.simps)
also have "... = sieve⋅(filter⋅(Λ (MkI⋅i). Def ((∀d::int. d > 1 ⟶ d < (int 2) ⟶ ¬ (d dvd i))))⋅[MkI⋅(int 2)..])"
unfolding numeral_Integer_eq
by (subst filter_TT, auto)
also have "... = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅(int 2)..]"
by (rule sieve_produces_primes[OF two_is_prime_nat])
also have "... = filter⋅(Λ (MkI⋅i). Def (prime (nat ¦i¦)))⋅[MkI⋅2..]"
by simp
finally show ?thesis .
qed
end