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