# HG changeset patch # User wenzelm # Date 1468437785 -7200 # Node ID 6e29fb72e65942174281178a4ed25684804d7c10 # Parent a4668ec480ad5f2170464c1673f29946c400e12d misc tuning and modernization; diff -r a4668ec480ad -r 6e29fb72e659 src/HOL/Library/Stirling.thy --- a/src/HOL/Library/Stirling.thy Wed Jul 13 20:48:18 2016 +0200 +++ b/src/HOL/Library/Stirling.thy Wed Jul 13 21:23:05 2016 +0200 @@ -1,5 +1,9 @@ -(* Authors: Amine Chaieb & Florian Haftmann, TU Muenchen - with contributions by Lukas Bulwahn and Manuel Eberl*) +(* Title: HOL/Library/Stirling.thy + Author: Amine Chaieb + Author: Florian Haftmann + Author: Lukas Bulwahn + Author: Manuel Eberl +*) section \Stirling numbers of first and second kind\ @@ -10,102 +14,105 @@ subsection \Stirling numbers of the second kind\ fun Stirling :: "nat \ nat \ nat" -where - "Stirling 0 0 = 1" -| "Stirling 0 (Suc k) = 0" -| "Stirling (Suc n) 0 = 0" -| "Stirling (Suc n) (Suc k) = Suc k * Stirling n (Suc k) + Stirling n k" + where + "Stirling 0 0 = 1" + | "Stirling 0 (Suc k) = 0" + | "Stirling (Suc n) 0 = 0" + | "Stirling (Suc n) (Suc k) = Suc k * Stirling n (Suc k) + Stirling n k" -lemma Stirling_1 [simp]: - "Stirling (Suc n) (Suc 0) = 1" +lemma Stirling_1 [simp]: "Stirling (Suc n) (Suc 0) = 1" by (induct n) simp_all -lemma Stirling_less [simp]: - "n < k \ Stirling n k = 0" +lemma Stirling_less [simp]: "n < k \ Stirling n k = 0" by (induct n k rule: Stirling.induct) simp_all -lemma Stirling_same [simp]: - "Stirling n n = 1" +lemma Stirling_same [simp]: "Stirling n n = 1" by (induct n) simp_all -lemma Stirling_2_2: - "Stirling (Suc (Suc n)) (Suc (Suc 0)) = 2 ^ Suc n - 1" +lemma Stirling_2_2: "Stirling (Suc (Suc n)) (Suc (Suc 0)) = 2 ^ Suc n - 1" proof (induct n) - case 0 then show ?case by simp + case 0 + then show ?case by simp next case (Suc n) have "Stirling (Suc (Suc (Suc n))) (Suc (Suc 0)) = - 2 * Stirling (Suc (Suc n)) (Suc (Suc 0)) + Stirling (Suc (Suc n)) (Suc 0)" by simp + 2 * Stirling (Suc (Suc n)) (Suc (Suc 0)) + Stirling (Suc (Suc n)) (Suc 0)" + by simp also have "\ = 2 * (2 ^ Suc n - 1) + 1" by (simp only: Suc Stirling_1) also have "\ = 2 ^ Suc (Suc n) - 1" proof - - have "(2::nat) ^ Suc n - 1 > 0" by (induct n) simp_all - then have "2 * ((2::nat) ^ Suc n - 1) > 0" by simp - then have "2 \ 2 * ((2::nat) ^ Suc n)" by simp + have "(2::nat) ^ Suc n - 1 > 0" + by (induct n) simp_all + then have "2 * ((2::nat) ^ Suc n - 1) > 0" + by simp + then have "2 \ 2 * ((2::nat) ^ Suc n)" + by simp with add_diff_assoc2 [of 2 "2 * 2 ^ Suc n" 1] - have "2 * 2 ^ Suc n - 2 + (1::nat) = 2 * 2 ^ Suc n + 1 - 2" . - then show ?thesis by (simp add: nat_distrib) + have "2 * 2 ^ Suc n - 2 + (1::nat) = 2 * 2 ^ Suc n + 1 - 2" . + then show ?thesis + by (simp add: nat_distrib) qed finally show ?case by simp qed -lemma Stirling_2: - "Stirling (Suc n) (Suc (Suc 0)) = 2 ^ n - 1" +lemma Stirling_2: "Stirling (Suc n) (Suc (Suc 0)) = 2 ^ n - 1" using Stirling_2_2 by (cases n) simp_all + subsection \Stirling numbers of the first kind\ fun stirling :: "nat \ nat \ nat" -where - "stirling 0 0 = 1" -| "stirling 0 (Suc k) = 0" -| "stirling (Suc n) 0 = 0" -| "stirling (Suc n) (Suc k) = n * stirling n (Suc k) + stirling n k" + where + "stirling 0 0 = 1" + | "stirling 0 (Suc k) = 0" + | "stirling (Suc n) 0 = 0" + | "stirling (Suc n) (Suc k) = n * stirling n (Suc k) + stirling n k" lemma stirling_0 [simp]: "n > 0 \ stirling n 0 = 0" by (cases n) simp_all -lemma stirling_less [simp]: - "n < k \ stirling n k = 0" +lemma stirling_less [simp]: "n < k \ stirling n k = 0" by (induct n k rule: stirling.induct) simp_all -lemma stirling_same [simp]: - "stirling n n = 1" +lemma stirling_same [simp]: "stirling n n = 1" by (induct n) simp_all -lemma stirling_Suc_n_1: - "stirling (Suc n) (Suc 0) = fact n" +lemma stirling_Suc_n_1: "stirling (Suc n) (Suc 0) = fact n" by (induct n) auto -lemma stirling_Suc_n_n: - shows "stirling (Suc n) n = Suc n choose 2" -by (induct n) (auto simp add: numerals(2)) +lemma stirling_Suc_n_n: "stirling (Suc n) n = Suc n choose 2" + by (induct n) (auto simp add: numerals(2)) lemma stirling_Suc_n_2: assumes "n \ Suc 0" shows "stirling (Suc n) 2 = (\k=1..n. fact n div k)" -using assms + using assms proof (induct n) - case 0 from this show ?case by simp + case 0 + then show ?case by simp next case (Suc n) show ?case proof (cases n) - case 0 from this show ?thesis by (simp add: numerals(2)) + case 0 + then show ?thesis + by (simp add: numerals(2)) next case Suc - from this have geq1: "Suc 0 \ n" by simp + then have geq1: "Suc 0 \ n" + by simp have "stirling (Suc (Suc n)) 2 = Suc n * stirling (Suc n) 2 + stirling (Suc n) (Suc 0)" by (simp only: stirling.simps(4)[of "Suc n"] numerals(2)) - also have "... = Suc n * (\k=1..n. fact n div k) + fact n" + also have "\ = Suc n * (\k=1..n. fact n div k) + fact n" using Suc.hyps[OF geq1] by (simp only: stirling_Suc_n_1 of_nat_fact of_nat_add of_nat_mult) - also have "... = Suc n * (\k=1..n. fact n div k) + Suc n * fact n div Suc n" + also have "\ = Suc n * (\k=1..n. fact n div k) + Suc n * fact n div Suc n" by (metis nat.distinct(1) nonzero_mult_divide_cancel_left) - also have "... = (\k=1..n. fact (Suc n) div k) + fact (Suc n) div Suc n" + also have "\ = (\k=1..n. fact (Suc n) div k) + fact (Suc n) div Suc n" by (simp add: setsum_right_distrib div_mult_swap dvd_fact) - also have "... = (\k=1..Suc n. fact (Suc n) div k)" by simp + also have "\ = (\k=1..Suc n. fact (Suc n) div k)" + by simp finally show ?thesis . qed qed @@ -113,52 +120,60 @@ lemma of_nat_stirling_Suc_n_2: assumes "n \ Suc 0" shows "(of_nat (stirling (Suc n) 2)::'a::field_char_0) = fact n * (\k=1..n. (1 / of_nat k))" -using assms + using assms proof (induct n) - case 0 from this show ?case by simp + case 0 + then show ?case by simp next case (Suc n) show ?case proof (cases n) - case 0 from this show ?thesis by (auto simp add: numerals(2)) + case 0 + then show ?thesis + by (auto simp add: numerals(2)) next case Suc - from this have geq1: "Suc 0 \ n" by simp + then have geq1: "Suc 0 \ n" + by simp have "(of_nat (stirling (Suc (Suc n)) 2)::'a) = - of_nat (Suc n * stirling (Suc n) 2 + stirling (Suc n) (Suc 0))" + of_nat (Suc n * stirling (Suc n) 2 + stirling (Suc n) (Suc 0))" by (simp only: stirling.simps(4)[of "Suc n"] numerals(2)) - also have "... = of_nat (Suc n) * (fact n * (\k = 1..n. 1 / of_nat k)) + fact n" + also have "\ = of_nat (Suc n) * (fact n * (\k = 1..n. 1 / of_nat k)) + fact n" using Suc.hyps[OF geq1] by (simp only: stirling_Suc_n_1 of_nat_fact of_nat_add of_nat_mult) - also have "... = fact (Suc n) * (\k = 1..n. 1 / of_nat k) + fact (Suc n) * (1 / of_nat (Suc n))" + also have "\ = fact (Suc n) * (\k = 1..n. 1 / of_nat k) + fact (Suc n) * (1 / of_nat (Suc n))" using of_nat_neq_0 by auto - also have "... = fact (Suc n) * (\k = 1..Suc n. 1 / of_nat k)" + also have "\ = fact (Suc n) * (\k = 1..Suc n. 1 / of_nat k)" by (simp add: distrib_left) finally show ?thesis . qed qed -lemma setsum_stirling: - "(\k\n. stirling n k) = fact n" +lemma setsum_stirling: "(\k\n. stirling n k) = fact n" proof (induct n) case 0 - from this show ?case by simp + then show ?case by simp next case (Suc n) have "(\k\Suc n. stirling (Suc n) k) = stirling (Suc n) 0 + (\k\n. stirling (Suc n) (Suc k))" by (simp only: setsum_atMost_Suc_shift) - also have "\ = (\k\n. stirling (Suc n) (Suc k))" by simp - also have "\ = (\k\n. n * stirling n (Suc k) + stirling n k)" by simp + also have "\ = (\k\n. stirling (Suc n) (Suc k))" + by simp + also have "\ = (\k\n. n * stirling n (Suc k) + stirling n k)" + by simp also have "\ = n * (\k\n. stirling n (Suc k)) + (\k\n. stirling n k)" by (simp add: setsum.distrib setsum_right_distrib) also have "\ = n * fact n + fact n" proof - have "n * (\k\n. stirling n (Suc k)) = n * ((\k\Suc n. stirling n k) - stirling n 0)" by (metis add_diff_cancel_left' setsum_atMost_Suc_shift) - also have "\ = n * (\k\n. stirling n k)" by (cases n) simp+ - also have "\ = n * fact n" using Suc.hyps by simp + also have "\ = n * (\k\n. stirling n k)" + by (cases n) simp_all + also have "\ = n * fact n" + using Suc.hyps by simp finally have "n * (\k\n. stirling n (Suc k)) = n * fact n" . - moreover have "(\k\n. stirling n k) = fact n" using Suc.hyps . + moreover have "(\k\n. stirling n k) = fact n" + using Suc.hyps . ultimately show ?thesis by simp qed also have "\ = fact (Suc n)" by simp @@ -166,26 +181,29 @@ qed lemma stirling_pochhammer: - "(\k\n. of_nat (stirling n k) * x ^ k) = (pochhammer x n :: 'a :: comm_semiring_1)" -proof (induction n) + "(\k\n. of_nat (stirling n k) * x ^ k) = (pochhammer x n :: 'a::comm_semiring_1)" +proof (induct n) + case 0 + then show ?case by simp +next case (Suc n) have "of_nat (n * stirling n 0) = (0 :: 'a)" by (cases n) simp_all - hence "(\k\Suc n. of_nat (stirling (Suc n) k) * x ^ k) = - (of_nat (n * stirling n 0) * x ^ 0 + - (\i\n. of_nat (n * stirling n (Suc i)) * (x ^ Suc i))) + - (\i\n. of_nat (stirling n i) * (x ^ Suc i))" + then have "(\k\Suc n. of_nat (stirling (Suc n) k) * x ^ k) = + (of_nat (n * stirling n 0) * x ^ 0 + + (\i\n. of_nat (n * stirling n (Suc i)) * (x ^ Suc i))) + + (\i\n. of_nat (stirling n i) * (x ^ Suc i))" by (subst setsum_atMost_Suc_shift) (simp add: setsum.distrib ring_distribs) also have "\ = pochhammer x (Suc n)" by (subst setsum_atMost_Suc_shift [symmetric]) - (simp add: algebra_simps setsum.distrib setsum_right_distrib pochhammer_Suc Suc [symmetric]) + (simp add: algebra_simps setsum.distrib setsum_right_distrib pochhammer_Suc Suc [symmetric]) finally show ?case . -qed simp_all +qed text \A row of the Stirling number triangle\ -definition stirling_row :: "nat \ nat list" where - "stirling_row n = [stirling n k. k \ [0.. nat list" + where "stirling_row n = [stirling n k. k \ [0.. n \ stirling_row n ! k = stirling n k" by (simp add: stirling_row_def del: upt_Suc) @@ -200,82 +218,109 @@ lemma list_ext: assumes "length xs = length ys" assumes "\i. i < length xs \ xs ! i = ys ! i" - shows "xs = ys" -using assms + shows "xs = ys" + using assms proof (induction rule: list_induct2) + case Nil + then show ?case by simp +next case (Cons x xs y ys) - from Cons.prems[of 0] have "x = y" by simp - moreover from Cons.prems[of "Suc i" for i] have "xs = ys" by (intro Cons.IH) simp + from Cons.prems[of 0] have "x = y" + by simp + moreover from Cons.prems[of "Suc i" for i] have "xs = ys" + by (intro Cons.IH) simp ultimately show ?case by simp -qed simp_all +qed subsubsection \Efficient code\ text \ - Naively using the defining equations of the Stirling numbers of the first kind to - compute them leads to exponential run time due to repeated computations. - We can use memoisation to compute them row by row without repeating computations, at - the cost of computing a few unneeded values. + Naively using the defining equations of the Stirling numbers of the first + kind to compute them leads to exponential run time due to repeated + computations. We can use memoisation to compute them row by row without + repeating computations, at the cost of computing a few unneeded values. As a bonus, this is very efficient for applications where an entire row of - Stirling numbers is needed.. + Stirling numbers is needed. \ -definition zip_with_prev :: "('a \ 'a \ 'b) \ 'a \ 'a list \ 'b list" where - "zip_with_prev f x xs = map (\(x,y). f x y) (zip (x # xs) xs)" +definition zip_with_prev :: "('a \ 'a \ 'b) \ 'a \ 'a list \ 'b list" + where "zip_with_prev f x xs = map (\(x,y). f x y) (zip (x # xs) xs)" lemma zip_with_prev_altdef: "zip_with_prev f x xs = - (if xs = [] then [] else f x (hd xs) # [f (xs!i) (xs!(i+1)). i \ [0.. [0..i. f (xs ! i) (xs ! (i + 1))) [0..a b. a + n * b) y xs @ [1]" - by (induction xs arbitrary: y) (simp_all add: zip_with_prev_def) + by (induct xs arbitrary: y) (simp_all add: zip_with_prev_def) lemma stirling_row_code [code]: "stirling_row 0 = [1]" "stirling_row (Suc n) = stirling_row_aux n 0 (stirling_row n)" -proof - +proof goal_cases + case 1 + show ?case by (simp add: stirling_row_def) +next + case 2 have "stirling_row (Suc n) = - 0 # [stirling_row n ! i + stirling_row n ! (i+1) * n. i \ [0.. [0.. Suc n" by simp - then consider "i = 0" | j where "i > 0" "i \ n" | "i = Suc n" by linarith - thus ?case + from nth have "i \ Suc n" + by simp + then consider "i = 0 \ i = Suc n" | "i > 0" "i \ n" + by linarith + then show ?case proof cases - assume i: "i > 0" "i \ n" - from this show ?thesis by (cases i) (simp_all add: nth_append nth_stirling_row) - qed (simp_all add: nth_stirling_row nth_append) - qed simp + case 1 + then show ?thesis + by (auto simp: nth_stirling_row nth_append) + next + case 2 + then show ?thesis + by (cases i) (simp_all add: nth_append nth_stirling_row) + qed + next + case length + then show ?case by simp + qed also have "0 # [stirling_row n ! i + stirling_row n ! (i+1) * n. i \ [0..a b. a + n * b) 0 (stirling_row n) @ [1]" + zip_with_prev (\a b. a + n * b) 0 (stirling_row n) @ [1]" by (cases n) (auto simp add: zip_with_prev_altdef stirling_row_def hd_map simp del: upt_Suc) also have "\ = stirling_row_aux n 0 (stirling_row n)" by (simp add: stirling_row_aux_correct) - finally show "stirling_row (Suc n) = stirling_row_aux n 0 (stirling_row n)" . -qed (simp add: stirling_row_def) + finally show ?case . +qed lemma stirling_code [code]: - "stirling n k = (if k = 0 then if n = 0 then 1 else 0 - else if k > n then 0 else if k = n then 1 - else stirling_row n ! k)" + "stirling n k = + (if k = 0 then (if n = 0 then 1 else 0) + else if k > n then 0 + else if k = n then 1 + else stirling_row n ! k)" by (simp add: nth_stirling_row) end