src/HOL/Number_Theory/Fib.thy
author eberlm <eberlm@in.tum.de>
Thu Oct 20 13:53:36 2016 +0200 (2016-10-20)
changeset 64317 029e6247210e
parent 64267 b9a1486e79be
child 65393 079a6f850c02
permissions -rw-r--r--
More on Fibonacci numbers
     1 (*  Title:      HOL/Number_Theory/Fib.thy
     2     Author:     Lawrence C. Paulson
     3     Author:     Jeremy Avigad
     4     Author:     Manuel Eberl
     5 *)
     6 
     7 section \<open>The fibonacci function\<close>
     8 
     9 theory Fib
    10   imports Complex_Main
    11 begin
    12 
    13 
    14 subsection \<open>Fibonacci numbers\<close>
    15 
    16 fun fib :: "nat \<Rightarrow> nat"
    17 where
    18   fib0: "fib 0 = 0"
    19 | fib1: "fib (Suc 0) = 1"
    20 | fib2: "fib (Suc (Suc n)) = fib (Suc n) + fib n"
    21 
    22 
    23 subsection \<open>Basic Properties\<close>
    24 
    25 lemma fib_1 [simp]: "fib (1::nat) = 1"
    26   by (metis One_nat_def fib1)
    27 
    28 lemma fib_2 [simp]: "fib (2::nat) = 1"
    29   using fib.simps(3) [of 0]
    30   by (simp add: numeral_2_eq_2)
    31 
    32 lemma fib_plus_2: "fib (n + 2) = fib (n + 1) + fib n"
    33   by (metis Suc_eq_plus1 add_2_eq_Suc' fib.simps(3))
    34 
    35 lemma fib_add: "fib (Suc (n+k)) = fib (Suc k) * fib (Suc n) + fib k * fib n"
    36   by (induct n rule: fib.induct) (auto simp add: field_simps)
    37 
    38 lemma fib_neq_0_nat: "n > 0 \<Longrightarrow> fib n > 0"
    39   by (induct n rule: fib.induct) (auto simp add: )
    40 
    41 
    42 subsection \<open>More efficient code\<close>
    43 
    44 text \<open>
    45   The naive approach is very inefficient since the branching recursion leads to many
    46   values of @{term fib} being computed multiple times. We can avoid this by ``remembering''
    47   the last two values in the sequence, yielding a tail-recursive version.
    48   This is far from optimal (it takes roughly $O(n\cdot M(n))$ time where $M(n)$ is the 
    49   time required to multiply two $n$-bit integers), but much better than the naive version,
    50   which is exponential.
    51 \<close>
    52 
    53 fun gen_fib :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat" where
    54   "gen_fib a b 0 = a"
    55 | "gen_fib a b (Suc 0) = b"
    56 | "gen_fib a b (Suc (Suc n)) = gen_fib b (a + b) (Suc n)"
    57 
    58 lemma gen_fib_recurrence: "gen_fib a b (Suc (Suc n)) = gen_fib a b n + gen_fib a b (Suc n)"
    59   by (induction a b n rule: gen_fib.induct) simp_all
    60   
    61 lemma gen_fib_fib: "gen_fib (fib n) (fib (Suc n)) m = fib (n + m)"
    62   by (induction m rule: fib.induct) (simp_all del: gen_fib.simps(3) add: gen_fib_recurrence)
    63 
    64 lemma fib_conv_gen_fib: "fib n = gen_fib 0 1 n"
    65   using gen_fib_fib[of 0 n] by simp
    66 
    67 declare fib_conv_gen_fib [code]
    68 
    69 
    70 subsection \<open>A Few Elementary Results\<close>
    71 
    72 text \<open>
    73   \medskip Concrete Mathematics, page 278: Cassini's identity.  The proof is
    74   much easier using integers, not natural numbers!
    75 \<close>
    76 
    77 lemma fib_Cassini_int: "int (fib (Suc (Suc n)) * fib n) - int((fib (Suc n))\<^sup>2) = - ((-1)^n)"
    78   by (induct n rule: fib.induct)  (auto simp add: field_simps power2_eq_square power_add)
    79 
    80 lemma fib_Cassini_nat:
    81   "fib (Suc (Suc n)) * fib n =
    82      (if even n then (fib (Suc n))\<^sup>2 - 1 else (fib (Suc n))\<^sup>2 + 1)"
    83   using fib_Cassini_int [of n]  by (auto simp del: of_nat_mult of_nat_power)
    84 
    85 
    86 subsection \<open>Law 6.111 of Concrete Mathematics\<close>
    87 
    88 lemma coprime_fib_Suc_nat: "coprime (fib (n::nat)) (fib (Suc n))"
    89   apply (induct n rule: fib.induct)
    90   apply auto
    91   apply (metis gcd_add1 add.commute)
    92   done
    93 
    94 lemma gcd_fib_add: "gcd (fib m) (fib (n + m)) = gcd (fib m) (fib n)"
    95   apply (simp add: gcd.commute [of "fib m"])
    96   apply (cases m)
    97   apply (auto simp add: fib_add)
    98   apply (metis gcd.commute mult.commute coprime_fib_Suc_nat
    99     gcd_add_mult gcd_mult_cancel gcd.commute)
   100   done
   101 
   102 lemma gcd_fib_diff: "m \<le> n \<Longrightarrow> gcd (fib m) (fib (n - m)) = gcd (fib m) (fib n)"
   103   by (simp add: gcd_fib_add [symmetric, of _ "n-m"])
   104 
   105 lemma gcd_fib_mod: "0 < m \<Longrightarrow> gcd (fib m) (fib (n mod m)) = gcd (fib m) (fib n)"
   106 proof (induct n rule: less_induct)
   107   case (less n)
   108   show "gcd (fib m) (fib (n mod m)) = gcd (fib m) (fib n)"
   109   proof (cases "m < n")
   110     case True
   111     then have "m \<le> n" by auto
   112     with \<open>0 < m\<close> have pos_n: "0 < n" by auto
   113     with \<open>0 < m\<close> \<open>m < n\<close> have diff: "n - m < n" by auto
   114     have "gcd (fib m) (fib (n mod m)) = gcd (fib m) (fib ((n - m) mod m))"
   115       by (simp add: mod_if [of n]) (insert \<open>m < n\<close>, auto)
   116     also have "\<dots> = gcd (fib m)  (fib (n - m))"
   117       by (simp add: less.hyps diff \<open>0 < m\<close>)
   118     also have "\<dots> = gcd (fib m) (fib n)"
   119       by (simp add: gcd_fib_diff \<open>m \<le> n\<close>)
   120     finally show "gcd (fib m) (fib (n mod m)) = gcd (fib m) (fib n)" .
   121   next
   122     case False
   123     then show "gcd (fib m) (fib (n mod m)) = gcd (fib m) (fib n)"
   124       by (cases "m = n") auto
   125   qed
   126 qed
   127 
   128 lemma fib_gcd: "fib (gcd m n) = gcd (fib m) (fib n)"
   129     \<comment> \<open>Law 6.111\<close>
   130   by (induct m n rule: gcd_nat_induct) (simp_all add: gcd_non_0_nat gcd.commute gcd_fib_mod)
   131 
   132 theorem fib_mult_eq_sum_nat: "fib (Suc n) * fib n = (\<Sum>k \<in> {..n}. fib k * fib k)"
   133   by (induct n rule: nat.induct) (auto simp add:  field_simps)
   134 
   135 
   136 subsection \<open>Closed form\<close>
   137 
   138 lemma fib_closed_form:
   139   defines "\<phi> \<equiv> (1 + sqrt 5) / (2::real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2::real)"
   140   shows   "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"
   141 proof (induction n rule: fib.induct)
   142   fix n :: nat
   143   assume IH1: "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"
   144   assume IH2: "of_nat (fib (Suc n)) = (\<phi> ^ Suc n - \<psi> ^ Suc n) / sqrt 5"
   145   have "of_nat (fib (Suc (Suc n))) = of_nat (fib (Suc n)) + of_nat (fib n)" by simp
   146   also have "... = (\<phi>^n*(\<phi> + 1) - \<psi>^n*(\<psi> + 1)) / sqrt 5"
   147     by (simp add: IH1 IH2 field_simps)
   148   also have "\<phi> + 1 = \<phi>\<^sup>2" by (simp add: \<phi>_def field_simps power2_eq_square)
   149   also have "\<psi> + 1 = \<psi>\<^sup>2" by (simp add: \<psi>_def field_simps power2_eq_square)
   150   also have "\<phi>^n * \<phi>\<^sup>2 - \<psi>^n * \<psi>\<^sup>2 = \<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)"  
   151     by (simp add: power2_eq_square)
   152   finally show "of_nat (fib (Suc (Suc n))) = (\<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)) / sqrt 5" .
   153 qed (simp_all add: \<phi>_def \<psi>_def field_simps)
   154 
   155 lemma fib_closed_form':
   156   defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
   157   assumes "n > 0"
   158   shows   "fib n = round (\<phi> ^ n / sqrt 5)"
   159 proof (rule sym, rule round_unique')
   160   have "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> = \<bar>\<psi>\<bar> ^ n / sqrt 5"
   161     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power_abs)
   162   also {
   163     from assms have "\<bar>\<psi>\<bar>^n \<le> \<bar>\<psi>\<bar>^1"
   164       by (intro power_decreasing) (simp_all add: algebra_simps real_le_lsqrt)
   165     also have "... < sqrt 5 / 2" by (simp add: \<psi>_def field_simps)
   166     finally have "\<bar>\<psi>\<bar>^n / sqrt 5 < 1/2" by (simp add: field_simps)
   167   }
   168   finally show "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> < 1/2" .
   169 qed
   170 
   171 lemma fib_asymptotics:
   172   defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)"
   173   shows   "(\<lambda>n. real (fib n) / (\<phi> ^ n / sqrt 5)) \<longlonglongrightarrow> 1"
   174 proof -
   175   define \<psi> where "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
   176   have "\<phi> > 1" by (simp add: \<phi>_def)
   177   hence A: "\<phi> \<noteq> 0" by auto
   178   have "(\<lambda>n. (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 0"
   179     by (rule LIMSEQ_power_zero) (simp_all add: \<phi>_def \<psi>_def field_simps add_pos_pos)
   180   hence "(\<lambda>n. 1 - (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 1 - 0" by (intro tendsto_diff tendsto_const)
   181   with A show ?thesis
   182     by (simp add: divide_simps fib_closed_form [folded \<phi>_def \<psi>_def])
   183 qed
   184 
   185 
   186 subsection \<open>Divide-and-Conquer recurrence\<close>
   187 
   188 text \<open>
   189   The following divide-and-conquer recurrence allows for a more efficient computation 
   190   of Fibonacci numbers; however, it requires memoisation of values to be reasonably 
   191   efficient, cutting the number of values to be computed to logarithmically many instead of
   192   linearly many. The vast majority of the computation time is then actually spent on the 
   193   multiplication, since the output number is exponential in the input number.
   194 \<close>
   195 
   196 lemma fib_rec_odd:
   197   defines "\<phi> \<equiv> (1 + sqrt 5) / (2 :: real)" and "\<psi> \<equiv> (1 - sqrt 5) / (2 :: real)"
   198   shows   "fib (Suc (2*n)) = fib n^2 + fib (Suc n)^2"
   199 proof -
   200   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"
   201     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power2_eq_square)
   202   also have "(\<phi> ^ n - \<psi> ^ n)\<^sup>2 + (\<phi> * \<phi> ^ n - \<psi> * \<psi> ^ n)\<^sup>2 = 
   203     \<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")
   204       by (simp add: power2_eq_square algebra_simps power_mult power_mult_distrib)
   205   also have "\<phi> * \<psi> = -1" by (simp add: \<phi>_def \<psi>_def field_simps)
   206   hence "?A = \<phi>^(2*n+1) * (\<phi> + inverse \<phi>) + \<psi>^(2*n+1) * (\<psi> + inverse \<psi>)" 
   207     by (auto simp: field_simps power2_eq_square)
   208   also have "1 + sqrt 5 > 0" by (auto intro: add_pos_pos)
   209   hence "\<phi> + inverse \<phi> = sqrt 5" by (simp add: \<phi>_def field_simps)
   210   also have "\<psi> + inverse \<psi> = -sqrt 5" by (simp add: \<psi>_def field_simps)
   211   also have "(\<phi> ^ (2*n+1) * sqrt 5 + \<psi> ^ (2*n+1)* - sqrt 5) / 5 =
   212                (\<phi> ^ (2*n+1) - \<psi> ^ (2*n+1)) * (sqrt 5 / 5)" by (simp add: field_simps)
   213   also have "sqrt 5 / 5 = inverse (sqrt 5)" by (simp add: field_simps)
   214   also have "(\<phi> ^ (2*n+1) - \<psi> ^ (2*n+1)) * ... = of_nat (fib (Suc (2*n)))"
   215     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] divide_inverse)
   216   finally show ?thesis by (simp only: of_nat_eq_iff)
   217 qed
   218 
   219 lemma fib_rec_even: "fib (2*n) = (fib (n - 1) + fib (n + 1)) * fib n"
   220 proof (induction n)
   221   case (Suc n)
   222   let ?rfib = "\<lambda>x. real (fib x)"
   223   have "2 * (Suc n) = Suc (Suc (2*n))" by simp
   224   also have "real (fib ...) = ?rfib n^2 + ?rfib (Suc n)^2 + (?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n" 
   225     by (simp add: fib_rec_odd Suc)
   226   also have "(?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n = (2 * ?rfib (n + 1) - ?rfib n) * ?rfib n"
   227     by (cases n) simp_all
   228   also have "?rfib n^2 + ?rfib (Suc n)^2 + ... = (?rfib (Suc n) + 2 * ?rfib n) * ?rfib (Suc n)"
   229     by (simp add: algebra_simps power2_eq_square)
   230   also have "... = real ((fib (Suc n - 1) + fib (Suc n + 1)) * fib (Suc n))" by simp
   231   finally show ?case by (simp only: of_nat_eq_iff)
   232 qed simp
   233 
   234 lemma fib_rec_even': "fib (2*n) = (2*fib (n - 1) + fib n) * fib n"
   235   by (subst fib_rec_even, cases n) simp_all
   236 
   237 lemma fib_rec:
   238   "fib n = (if n = 0 then 0 else if n = 1 then 1 else
   239             if even n then let n' = n div 2; fn = fib n' in (2 * fib (n' - 1) + fn) * fn
   240             else let n' = n div 2 in fib n' ^ 2 + fib (Suc n') ^ 2)"
   241   by (auto elim: evenE oddE simp: fib_rec_odd fib_rec_even' Let_def)
   242 
   243 
   244 subsection \<open>Fibonacci and Binomial Coefficients\<close>
   245 
   246 lemma sum_drop_zero: "(\<Sum>k = 0..Suc n. if 0<k then (f (k - 1)) else 0) = (\<Sum>j = 0..n. f j)"
   247   by (induct n) auto
   248 
   249 lemma sum_choose_drop_zero:
   250     "(\<Sum>k = 0..Suc n. if k=0 then 0 else (Suc n - k) choose (k - 1)) = (\<Sum>j = 0..n. (n-j) choose j)"
   251   by (rule trans [OF sum.cong sum_drop_zero]) auto
   252 
   253 lemma ne_diagonal_fib: "(\<Sum>k = 0..n. (n-k) choose k) = fib (Suc n)"
   254 proof (induct n rule: fib.induct)
   255   case 1
   256   show ?case by simp
   257 next
   258   case 2
   259   show ?case by simp
   260 next
   261   case (3 n)
   262   have "(\<Sum>k = 0..Suc n. Suc (Suc n) - k choose k) =
   263         (\<Sum>k = 0..Suc n. (Suc n - k choose k) + (if k=0 then 0 else (Suc n - k choose (k - 1))))"
   264     by (rule sum.cong) (simp_all add: choose_reduce_nat)
   265   also have "\<dots> = (\<Sum>k = 0..Suc n. Suc n - k choose k) +
   266                    (\<Sum>k = 0..Suc n. if k=0 then 0 else (Suc n - k choose (k - 1)))"
   267     by (simp add: sum.distrib)
   268   also have "\<dots> = (\<Sum>k = 0..Suc n. Suc n - k choose k) +
   269                    (\<Sum>j = 0..n. n - j choose j)"
   270     by (metis sum_choose_drop_zero)
   271   finally show ?case using 3
   272     by simp
   273 qed
   274 
   275 end
   276