MATH
21 B
Mathematics Math21b Spring 2017
Linear Algebra and Differential Equations
Exhibit: Liber Abacus
Course Head: Oliver Knill
Office: SciCtr 432

Exhibit: Liber Abacus

There is hardly any sequence of numbers, except maybe the primes which is as famous as the Fibonacci sequence. Similarly, there is hardly any real number besides π and the Euler-e which is as famous as the golden ratio φ = (1+51/2)/2. Linear algebra relates the Fibonacci sequence with the golden ratio. If the recursion x(n+1) = x(n) + x(n-1) with the initial condition x(0)=0,x(1)=1 is solved by writing it as a discrete dynamical system v(n+1) = A v(n) written out as
    | x(n+1) |    | 1   1  |   | x(n)   |     | x(1) |   | 1 |
    | x(n)   |  = | 1   0  |   | x(n-1) |,    | x(0) | = | 0 |
The eigenvalues are φ=1.61803... and -1/φ. To solve the system, we write the initial condition as a linear combination of the eigenvectors
| 1 |      |  φ   |         | -1/φ |
| 0 |  = ( |   1  |    -    |   1  |  )/ 51/2 
then write down the closed form solution
| x(n+1) |        |  φ |              | -1/φ  |
| x(n)   | = (φn  |  1 |  -  (-1/φ)n  |    1  |  )/ 51/2
Looking at the second coordinates gives
 x(n) = [( (1+51/2)/2)n - (1-51/2)/2)n ]/51/2
In Mathematica
 phi=(1+Sqrt[5])/2; 
 F[n_]:=(phi^n - (-1/phi)^n)/Sqrt[5];
 Table[Simplify[F[n]],{n,10}]
Now, Mathematica has this already implemented. The same result can be obtained with
  Table[Fibonacci[n],{n,10}]
We could also program the recursion
  G[0]=0; G[1]=1; G[n_]:=G[n-1]+G[n-2];
  Table[G[n],{n,10}]
Now, you might think, whats the point? If the procedure is already built into the system or can be programmed with recursion of half a line, why do we need the explicit solution? Lets look how long a machine needs to compute the built in Fibonacci number for n=1000000000 or the recursion or then using the F formula above:
First[Timing[Fibonacci[10^9]]]
Now, the recursion with G fails with an error "Recursion depth of 1024 exceeded during evaluation". Well, lets increase the recursion limit:
$RecursionLimit=10^9; G[0]=0; G[1]=1; G[n_]:=G[n-1]+G[n-2]; First[Timing[G[10^9]]]
If you run it you see it kills the program. One of the ways to segfault Mathematica! What is nifty about the explicit formula is that we can just say that F[1000000000] is equal to (phi^(1000000000) - (-1/phi)^(100000000))/Sqrt[5]. This takes 50 milli seconds:
phi=(1+Sqrt[5])/2; F[n_]:=(phi^n - (-1/phi)^n)/Sqrt[5]; First[Timing[F[10^9]]]
However, doing the same command with 10^13 gives an error "Uncaught SystemException returned to top level". The memory resources to handle such large algebraic expression was too large. Folio 139v-140r of the Liber Abbaci. (Biblioteca Nazionale Centrale, Florence, Italy)

Source
Please send questions and comments to knill@math.harvard.edu
Math21b Harvard College Course ID:110989| Oliver Knill | Spring 2017 | Department of Mathematics | Faculty of Art and Sciences | Harvard University, [Canvas, for admin], Twitter