src/HOL/Number_Theory/Fib.thy
 author eberlm 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