src/HOL/Number_Theory/Fib.thy
 author wenzelm Wed Apr 05 22:25:18 2017 +0200 (2017-04-05) changeset 65393 079a6f850c02 parent 64317 029e6247210e child 67051 e7e54a0b9197 permissions -rw-r--r--
misc tuning and modernization;
     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 = 1"

    26   by (metis One_nat_def fib1)

    27

    28 lemma fib_2 [simp]: "fib 2 = 1"

    29   using fib.simps(3) [of 0] by (simp add: numeral_2_eq_2)

    30

    31 lemma fib_plus_2: "fib (n + 2) = fib (n + 1) + fib n"

    32   by (metis Suc_eq_plus1 add_2_eq_Suc' fib.simps(3))

    33

    34 lemma fib_add: "fib (Suc (n + k)) = fib (Suc k) * fib (Suc n) + fib k * fib n"

    35   by (induct n rule: fib.induct) (auto simp add: field_simps)

    36

    37 lemma fib_neq_0_nat: "n > 0 \<Longrightarrow> fib n > 0"

    38   by (induct n rule: fib.induct) auto

    39

    40

    41 subsection \<open>More efficient code\<close>

    42

    43 text \<open>

    44   The naive approach is very inefficient since the branching recursion leads to many

    45   values of @{term fib} being computed multiple times. We can avoid this by remembering''

    46   the last two values in the sequence, yielding a tail-recursive version.

    47   This is far from optimal (it takes roughly $O(n\cdot M(n))$ time where $M(n)$ is the

    48   time required to multiply two $n$-bit integers), but much better than the naive version,

    49   which is exponential.

    50 \<close>

    51

    52 fun gen_fib :: "nat \<Rightarrow> nat \<Rightarrow> nat \<Rightarrow> nat"

    53   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 (induct 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 (induct 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) (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 "0 < n" by auto

   113     with \<open>0 < m\<close> \<open>m < n\<close> have *: "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]) (use \<open>m < n\<close> in auto)

   116     also have "\<dots> = gcd (fib m)  (fib (n - m))"

   117       by (simp add: less.hyps * \<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)"  \<comment> \<open>Law 6.111\<close>

   129   by (induct m n rule: gcd_nat_induct) (simp_all add: gcd_non_0_nat gcd.commute gcd_fib_mod)

   130

   131 theorem fib_mult_eq_sum_nat: "fib (Suc n) * fib n = (\<Sum>k \<in> {..n}. fib k * fib k)"

   132   by (induct n rule: nat.induct) (auto simp add:  field_simps)

   133

   134

   135 subsection \<open>Closed form\<close>

   136

   137 lemma fib_closed_form:

   138   fixes \<phi> \<psi> :: real

   139   defines "\<phi> \<equiv> (1 + sqrt 5) / 2"

   140     and "\<psi> \<equiv> (1 - sqrt 5) / 2"

   141   shows "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"

   142 proof (induct n rule: fib.induct)

   143   fix n :: nat

   144   assume IH1: "of_nat (fib n) = (\<phi> ^ n - \<psi> ^ n) / sqrt 5"

   145   assume IH2: "of_nat (fib (Suc n)) = (\<phi> ^ Suc n - \<psi> ^ Suc n) / sqrt 5"

   146   have "of_nat (fib (Suc (Suc n))) = of_nat (fib (Suc n)) + of_nat (fib n)" by simp

   147   also have "\<dots> = (\<phi>^n * (\<phi> + 1) - \<psi>^n * (\<psi> + 1)) / sqrt 5"

   148     by (simp add: IH1 IH2 field_simps)

   149   also have "\<phi> + 1 = \<phi>\<^sup>2" by (simp add: \<phi>_def field_simps power2_eq_square)

   150   also have "\<psi> + 1 = \<psi>\<^sup>2" by (simp add: \<psi>_def field_simps power2_eq_square)

   151   also have "\<phi>^n * \<phi>\<^sup>2 - \<psi>^n * \<psi>\<^sup>2 = \<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)"

   152     by (simp add: power2_eq_square)

   153   finally show "of_nat (fib (Suc (Suc n))) = (\<phi> ^ Suc (Suc n) - \<psi> ^ Suc (Suc n)) / sqrt 5" .

   154 qed (simp_all add: \<phi>_def \<psi>_def field_simps)

   155

   156 lemma fib_closed_form':

   157   fixes \<phi> \<psi> :: real

   158   defines "\<phi> \<equiv> (1 + sqrt 5) / 2"

   159     and "\<psi> \<equiv> (1 - sqrt 5) / 2"

   160   assumes "n > 0"

   161   shows "fib n = round (\<phi> ^ n / sqrt 5)"

   162 proof (rule sym, rule round_unique')

   163   have "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> = \<bar>\<psi>\<bar> ^ n / sqrt 5"

   164     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power_abs)

   165   also {

   166     from assms have "\<bar>\<psi>\<bar>^n \<le> \<bar>\<psi>\<bar>^1"

   167       by (intro power_decreasing) (simp_all add: algebra_simps real_le_lsqrt)

   168     also have "\<dots> < sqrt 5 / 2" by (simp add: \<psi>_def field_simps)

   169     finally have "\<bar>\<psi>\<bar>^n / sqrt 5 < 1/2" by (simp add: field_simps)

   170   }

   171   finally show "\<bar>\<phi> ^ n / sqrt 5 - of_int (int (fib n))\<bar> < 1/2" .

   172 qed

   173

   174 lemma fib_asymptotics:

   175   fixes \<phi> :: real

   176   defines "\<phi> \<equiv> (1 + sqrt 5) / 2"

   177   shows "(\<lambda>n. real (fib n) / (\<phi> ^ n / sqrt 5)) \<longlonglongrightarrow> 1"

   178 proof -

   179   define \<psi> :: real where "\<psi> \<equiv> (1 - sqrt 5) / 2"

   180   have "\<phi> > 1" by (simp add: \<phi>_def)

   181   then have *: "\<phi> \<noteq> 0" by auto

   182   have "(\<lambda>n. (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 0"

   183     by (rule LIMSEQ_power_zero) (simp_all add: \<phi>_def \<psi>_def field_simps add_pos_pos)

   184   then have "(\<lambda>n. 1 - (\<psi> / \<phi>) ^ n) \<longlonglongrightarrow> 1 - 0"

   185     by (intro tendsto_diff tendsto_const)

   186   with * show ?thesis

   187     by (simp add: divide_simps fib_closed_form [folded \<phi>_def \<psi>_def])

   188 qed

   189

   190

   191 subsection \<open>Divide-and-Conquer recurrence\<close>

   192

   193 text \<open>

   194   The following divide-and-conquer recurrence allows for a more efficient computation

   195   of Fibonacci numbers; however, it requires memoisation of values to be reasonably

   196   efficient, cutting the number of values to be computed to logarithmically many instead of

   197   linearly many. The vast majority of the computation time is then actually spent on the

   198   multiplication, since the output number is exponential in the input number.

   199 \<close>

   200

   201 lemma fib_rec_odd:

   202   fixes \<phi> \<psi> :: real

   203   defines "\<phi> \<equiv> (1 + sqrt 5) / 2"

   204     and "\<psi> \<equiv> (1 - sqrt 5) / 2"

   205   shows "fib (Suc (2 * n)) = fib n^2 + fib (Suc n)^2"

   206 proof -

   207   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"

   208     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] field_simps power2_eq_square)

   209   also

   210   let ?A = "\<phi>^(2 * n) + \<psi>^(2 * n) - 2*(\<phi> * \<psi>)^n + \<phi>^(2 * n + 2) + \<psi>^(2 * n + 2) - 2*(\<phi> * \<psi>)^(n + 1)"

   211   have "(\<phi> ^ n - \<psi> ^ n)\<^sup>2 + (\<phi> * \<phi> ^ n - \<psi> * \<psi> ^ n)\<^sup>2 = ?A"

   212     by (simp add: power2_eq_square algebra_simps power_mult power_mult_distrib)

   213   also have "\<phi> * \<psi> = -1"

   214     by (simp add: \<phi>_def \<psi>_def field_simps)

   215   then have "?A = \<phi>^(2 * n + 1) * (\<phi> + inverse \<phi>) + \<psi>^(2 * n + 1) * (\<psi> + inverse \<psi>)"

   216     by (auto simp: field_simps power2_eq_square)

   217   also have "1 + sqrt 5 > 0"

   218     by (auto intro: add_pos_pos)

   219   then have "\<phi> + inverse \<phi> = sqrt 5"

   220     by (simp add: \<phi>_def field_simps)

   221   also have "\<psi> + inverse \<psi> = -sqrt 5"

   222     by (simp add: \<psi>_def field_simps)

   223   also have "(\<phi> ^ (2 * n + 1) * sqrt 5 + \<psi> ^ (2 * n + 1) * - sqrt 5) / 5 =

   224     (\<phi> ^ (2 * n + 1) - \<psi> ^ (2 * n + 1)) * (sqrt 5 / 5)"

   225     by (simp add: field_simps)

   226   also have "sqrt 5 / 5 = inverse (sqrt 5)"

   227     by (simp add: field_simps)

   228   also have "(\<phi> ^ (2 * n + 1) - \<psi> ^ (2 * n + 1)) * \<dots> = of_nat (fib (Suc (2 * n)))"

   229     by (simp add: fib_closed_form[folded \<phi>_def \<psi>_def] divide_inverse)

   230   finally show ?thesis

   231     by (simp only: of_nat_eq_iff)

   232 qed

   233

   234 lemma fib_rec_even: "fib (2 * n) = (fib (n - 1) + fib (n + 1)) * fib n"

   235 proof (induct n)

   236   case 0

   237   then show ?case by simp

   238 next

   239   case (Suc n)

   240   let ?rfib = "\<lambda>x. real (fib x)"

   241   have "2 * (Suc n) = Suc (Suc (2 * n))" by simp

   242   also have "real (fib \<dots>) = ?rfib n^2 + ?rfib (Suc n)^2 + (?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n"

   243     by (simp add: fib_rec_odd Suc)

   244   also have "(?rfib (n - 1) + ?rfib (n + 1)) * ?rfib n = (2 * ?rfib (n + 1) - ?rfib n) * ?rfib n"

   245     by (cases n) simp_all

   246   also have "?rfib n^2 + ?rfib (Suc n)^2 + \<dots> = (?rfib (Suc n) + 2 * ?rfib n) * ?rfib (Suc n)"

   247     by (simp add: algebra_simps power2_eq_square)

   248   also have "\<dots> = real ((fib (Suc n - 1) + fib (Suc n + 1)) * fib (Suc n))" by simp

   249   finally show ?case by (simp only: of_nat_eq_iff)

   250 qed

   251

   252 lemma fib_rec_even': "fib (2 * n) = (2 * fib (n - 1) + fib n) * fib n"

   253   by (subst fib_rec_even, cases n) simp_all

   254

   255 lemma fib_rec:

   256   "fib n =

   257     (if n = 0 then 0 else if n = 1 then 1

   258      else if even n then let n' = n div 2; fn = fib n' in (2 * fib (n' - 1) + fn) * fn

   259      else let n' = n div 2 in fib n' ^ 2 + fib (Suc n') ^ 2)"

   260   by (auto elim: evenE oddE simp: fib_rec_odd fib_rec_even' Let_def)

   261

   262

   263 subsection \<open>Fibonacci and Binomial Coefficients\<close>

   264

   265 lemma sum_drop_zero: "(\<Sum>k = 0..Suc n. if 0<k then (f (k - 1)) else 0) = (\<Sum>j = 0..n. f j)"

   266   by (induct n) auto

   267

   268 lemma sum_choose_drop_zero:

   269   "(\<Sum>k = 0..Suc n. if k = 0 then 0 else (Suc n - k) choose (k - 1)) =

   270     (\<Sum>j = 0..n. (n-j) choose j)"

   271   by (rule trans [OF sum.cong sum_drop_zero]) auto

   272

   273 lemma ne_diagonal_fib: "(\<Sum>k = 0..n. (n-k) choose k) = fib (Suc n)"

   274 proof (induct n rule: fib.induct)

   275   case 1

   276   show ?case by simp

   277 next

   278   case 2

   279   show ?case by simp

   280 next

   281   case (3 n)

   282   have "(\<Sum>k = 0..Suc n. Suc (Suc n) - k choose k) =

   283      (\<Sum>k = 0..Suc n. (Suc n - k choose k) + (if k = 0 then 0 else (Suc n - k choose (k - 1))))"

   284     by (rule sum.cong) (simp_all add: choose_reduce_nat)

   285   also have "\<dots> =

   286     (\<Sum>k = 0..Suc n. Suc n - k choose k) +

   287     (\<Sum>k = 0..Suc n. if k=0 then 0 else (Suc n - k choose (k - 1)))"

   288     by (simp add: sum.distrib)

   289   also have "\<dots> = (\<Sum>k = 0..Suc n. Suc n - k choose k) + (\<Sum>j = 0..n. n - j choose j)"

   290     by (metis sum_choose_drop_zero)

   291   finally show ?case using 3

   292     by simp

   293 qed

   294

   295 end