src/HOL/Number_Theory/Fib.thy
 changeset 64317 029e6247210e parent 64267 b9a1486e79be child 65393 079a6f850c02
     1.1 --- a/src/HOL/Number_Theory/Fib.thy	Thu Oct 20 11:05:16 2016 +0200
1.2 +++ b/src/HOL/Number_Theory/Fib.thy	Thu Oct 20 13:53:36 2016 +0200
1.3 @@ -1,12 +1,13 @@
1.4  (*  Title:      HOL/Number_Theory/Fib.thy
1.5      Author:     Lawrence C. Paulson
1.7 +    Author:     Manuel Eberl
1.8  *)
1.9
1.10  section \<open>The fibonacci function\<close>
1.11
1.12  theory Fib
1.13 -imports Main GCD Binomial
1.14 +  imports Complex_Main
1.15  begin
1.16
1.17
1.18 @@ -38,6 +39,34 @@
1.19    by (induct n rule: fib.induct) (auto simp add: )
1.20
1.21
1.22 +subsection \<open>More efficient code\<close>
1.23 +
1.24 +text \<open>
1.25 +  The naive approach is very inefficient since the branching recursion leads to many
1.26 +  values of @{term fib} being computed multiple times. We can avoid this by remembering''
1.27 +  the last two values in the sequence, yielding a tail-recursive version.
1.28 +  This is far from optimal (it takes roughly $O(n\cdot M(n))$ time where $M(n)$ is the
1.29 +  time required to multiply two $n$-bit integers), but much better than the naive version,
1.30 +  which is exponential.
1.31 +\<close>
1.32 +
1.33 +fun gen_fib :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat" where
1.34 +  "gen_fib a b 0 = a"
1.35 +| "gen_fib a b (Suc 0) = b"
1.36 +| "gen_fib a b (Suc (Suc n)) = gen_fib b (a + b) (Suc n)"
1.37 +
1.38 +lemma gen_fib_recurrence: "gen_fib a b (Suc (Suc n)) = gen_fib a b n + gen_fib a b (Suc n)"
1.39 +  by (induction a b n rule: gen_fib.induct) simp_all
1.40 +
1.41 +lemma gen_fib_fib: "gen_fib (fib n) (fib (Suc n)) m = fib (n + m)"
1.42 +  by (induction m rule: fib.induct) (simp_all del: gen_fib.simps(3) add: gen_fib_recurrence)
1.43 +
1.44 +lemma fib_conv_gen_fib: "fib n = gen_fib 0 1 n"
1.45 +  using gen_fib_fib[of 0 n] by simp
1.46 +
1.47 +declare fib_conv_gen_fib [code]
1.48 +
1.49 +
1.50  subsection \<open>A Few Elementary Results\<close>
1.51
1.52  text \<open>
1.53 @@ -104,6 +133,114 @@
1.54    by (induct n rule: nat.induct) (auto simp add:  field_simps)
1.55
1.56
1.57 +subsection \<open>Closed form\<close>
1.58 +
1.59 +lemma fib_closed_form:
1.60 +  defines "\<phi> \<equiv> (1 + sqrt 5) / (2::real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2::real)"
1.61 +  shows   "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"
1.62 +proof (induction n rule: fib.induct)
1.63 +  fix n :: nat
1.64 +  assume IH1: "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"
1.65 +  assume IH2: "of_nat (fib (Suc n)) = (\<phi> ^ Suc n - \<psi> ^ Suc n) / sqrt 5"
1.66 +  have "of_nat (fib (Suc (Suc n))) = of_nat (fib (Suc n)) + of_nat (fib n)" by simp
1.67 +  also have "... = (\<phi>^n*(\<phi> + 1) - \<psi>^n*(\<psi> + 1)) / sqrt 5"
1.68 +    by (simp add: IH1 IH2 field_simps)
1.69 +  also have "\<phi> + 1 = \<phi>\<^sup>2" by (simp add: \<phi>_def field_simps power2_eq_square)
1.70 +  also have "\<psi> + 1 = \<psi>\<^sup>2" by (simp add: \<psi>_def field_simps power2_eq_square)
1.71 +  also have "\<phi>^n * \<phi>\<^sup>2 - \<psi>^n * \<psi>\<^sup>2 = \<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)"
1.72 +    by (simp add: power2_eq_square)
1.73 +  finally show "of_nat (fib (Suc (Suc n))) = (\<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)) / sqrt 5" .
1.74 +qed (simp_all add: \<phi>_def \<psi>_def field_simps)
1.75 +
1.76 +lemma fib_closed_form':
1.77 +  defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
1.78 +  assumes "n > 0"
1.79 +  shows   "fib n = round (\<phi> ^ n / sqrt 5)"
1.80 +proof (rule sym, rule round_unique')
1.81 +  have "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> = \<bar>\<psi>\<bar> ^ n / sqrt 5"
1.82 +    by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power_abs)
1.83 +  also {
1.84 +    from assms have "\<bar>\<psi>\<bar>^n \<le> \<bar>\<psi>\<bar>^1"
1.85 +      by (intro power_decreasing) (simp_all add: algebra_simps real_le_lsqrt)
1.86 +    also have "... < sqrt 5 / 2" by (simp add: \<psi>_def field_simps)
1.87 +    finally have "\<bar>\<psi>\<bar>^n / sqrt 5 < 1/2" by (simp add: field_simps)
1.88 +  }
1.89 +  finally show "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> < 1/2" .
1.90 +qed
1.91 +
1.92 +lemma fib_asymptotics:
1.93 +  defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)"
1.94 +  shows   "(\<lambda>n. real (fib n) / (\<phi> ^ n / sqrt 5)) \<longlonglongrightarrow> 1"
1.95 +proof -
1.96 +  define \<psi> where "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
1.97 +  have "\<phi> > 1" by (simp add: \<phi>_def)
1.98 +  hence A: "\<phi> \<noteq> 0" by auto
1.99 +  have "(\<lambda>n. (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 0"
1.101 +  hence "(\<lambda>n. 1 - (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 1 - 0" by (intro tendsto_diff tendsto_const)
1.102 +  with A show ?thesis
1.103 +    by (simp add: divide_simps fib_closed_form [folded \<phi>_def \<psi>_def])
1.104 +qed
1.105 +
1.106 +
1.107 +subsection \<open>Divide-and-Conquer recurrence\<close>
1.108 +
1.109 +text \<open>
1.110 +  The following divide-and-conquer recurrence allows for a more efficient computation
1.111 +  of Fibonacci numbers; however, it requires memoisation of values to be reasonably
1.112 +  efficient, cutting the number of values to be computed to logarithmically many instead of
1.113 +  linearly many. The vast majority of the computation time is then actually spent on the
1.114 +  multiplication, since the output number is exponential in the input number.
1.115 +\<close>
1.116 +
1.117 +lemma fib_rec_odd:
1.118 +  defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
1.119 +  shows   "fib (Suc (2*n)) = fib n^2 + fib (Suc n)^2"
1.120 +proof -
1.121 +  have "of_nat (fib n^2 + fib (Suc n)^2) = ((\<phi> ^ n - \<psi> ^ n)\<^sup>2 + (\<phi> * \<phi> ^ n - \<psi> * \<psi> ^ n)\<^sup>2)/5"
1.122 +    by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power2_eq_square)
1.123 +  also have "(\<phi> ^ n - \<psi> ^ n)\<^sup>2 + (\<phi> * \<phi> ^ n - \<psi> * \<psi> ^ n)\<^sup>2 =
1.124 +    \<phi>^(2*n) + \<psi>^(2*n) - 2*(\<phi>*\<psi>)^n + \<phi>^(2*n+2) + \<psi>^(2*n+2) - 2*(\<phi>*\<psi>)^(n+1)" (is "_ = ?A")
1.125 +      by (simp add: power2_eq_square algebra_simps power_mult power_mult_distrib)
1.126 +  also have "\<phi> * \<psi> = -1" by (simp add: \<phi>_def \<psi>_def field_simps)
1.127 +  hence "?A = \<phi>^(2*n+1) * (\<phi> + inverse \<phi>) + \<psi>^(2*n+1) * (\<psi> + inverse \<psi>)"
1.128 +    by (auto simp: field_simps power2_eq_square)
1.129 +  also have "1 + sqrt 5 > 0" by (auto intro: add_pos_pos)
1.130 +  hence "\<phi> + inverse \<phi> = sqrt 5" by (simp add: \<phi>_def field_simps)
1.131 +  also have "\<psi> + inverse \<psi> = -sqrt 5" by (simp add: \<psi>_def field_simps)
1.132 +  also have "(\<phi> ^ (2*n+1) * sqrt 5 + \<psi> ^ (2*n+1)* - sqrt 5) / 5 =
1.133 +               (\<phi> ^ (2*n+1) - \<psi> ^ (2*n+1)) * (sqrt 5 / 5)" by (simp add: field_simps)
1.134 +  also have "sqrt 5 / 5 = inverse (sqrt 5)" by (simp add: field_simps)
1.135 +  also have "(\<phi> ^ (2*n+1) - \<psi> ^ (2*n+1)) * ... = of_nat (fib (Suc (2*n)))"
1.136 +    by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] divide_inverse)
1.137 +  finally show ?thesis by (simp only: of_nat_eq_iff)
1.138 +qed
1.139 +
1.140 +lemma fib_rec_even: "fib (2*n) = (fib (n - 1) + fib (n + 1)) * fib n"
1.141 +proof (induction n)
1.142 +  case (Suc n)
1.143 +  let ?rfib = "\<lambda>x. real (fib x)"
1.144 +  have "2 * (Suc n) = Suc (Suc (2*n))" by simp
1.145 +  also have "real (fib ...) = ?rfib n^2 + ?rfib (Suc n)^2 + (?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n"
1.146 +    by (simp add: fib_rec_odd Suc)
1.147 +  also have "(?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n = (2 * ?rfib (n + 1) - ?rfib n) * ?rfib n"
1.148 +    by (cases n) simp_all
1.149 +  also have "?rfib n^2 + ?rfib (Suc n)^2 + ... = (?rfib (Suc n) + 2 * ?rfib n) * ?rfib (Suc n)"
1.150 +    by (simp add: algebra_simps power2_eq_square)
1.151 +  also have "... = real ((fib (Suc n - 1) + fib (Suc n + 1)) * fib (Suc n))" by simp
1.152 +  finally show ?case by (simp only: of_nat_eq_iff)
1.153 +qed simp
1.154 +
1.155 +lemma fib_rec_even': "fib (2*n) = (2*fib (n - 1) + fib n) * fib n"
1.156 +  by (subst fib_rec_even, cases n) simp_all
1.157 +
1.158 +lemma fib_rec:
1.159 +  "fib n = (if n = 0 then 0 else if n = 1 then 1 else
1.160 +            if even n then let n' = n div 2; fn = fib n' in (2 * fib (n' - 1) + fn) * fn
1.161 +            else let n' = n div 2 in fib n' ^ 2 + fib (Suc n') ^ 2)"
1.162 +  by (auto elim: evenE oddE simp: fib_rec_odd fib_rec_even' Let_def)
1.163 +
1.164 +
1.165  subsection \<open>Fibonacci and Binomial Coefficients\<close>
1.166
1.167  lemma sum_drop_zero: "(\<Sum>k = 0..Suc n. if 0<k then (f (k - 1)) else 0) = (\<Sum>j = 0..n. f j)"