Fibonacci Numbers

We are going to study a sequence of numbers known as the Fibonacci numbers. These numbers arise in many situations, one of which we now describe.
Suppose on a certain day we have one pair of baby rabbits, and we know that in one month they will have matured to adulthood, and then one month later will produce a pair of baby rabbits, then another pair of babies two months later, another pair three months later, and so on. Furthermore, each new pair of baby rabbits will reach adulthood in one month and then produce a pair of babies one month later and each month thereafter. How many adult pairs will there be after 5 months, after 30 months, after any specific number of months? The numbers that answer these questions are called the Fibonacci numbers. If we denote the number of adult pairs after n months by , then we want to find or calculate . The first few Fibonacci numbers are easy to calculate by following the generations of our rabbit colony month-by-month; we find that = 0, = 1, = 1, = 2, and = 3. We would like, however, to have a simple, systematic way to find for any n; we will, in fact, discuss several ways to calculate .

The Built-in Function Fibonacci

Let's see if Mathematica has a built-in function for computing Fibonacci numbers.

?Fib*

So we see that Fibonacci[n] gives We illustrate the use of this command.

Fibonacci[0]
Fibonacci[1]
Fibonacci[2]
Fibonacci[3]

These results confirm the values previously calculated. And the following results answer the specific questions asked above.

Fibonacci[5]
Fibonacci[30]

So there are 5 pairs of adult rabbits after 5 months, and 832,040 adult pairs after 30 months. We can also make a table of Fibonacci numbers with the commands Table and TableForm. Here is a table of the first 16 Fibonacci numbers.

fibtable = Table[{n, Fibonacci[n]}, {n, 0, 15}]

TableForm[fibtable, TableDirections -> {Row, Column}, TableHeadings -> {None, {"n", "F[n]"}}, TableSpacing ->{1,1}]

The Fundamental Recurrence Relation

We have used the Fibonacci command as a "black box", with no idea how the calculations are done. We now turn to a discussion of the computation of Fibonacci numbers. In addition to explaining methods of computation, this discussion will introduce several useful Mathematica commands. In this section we show that the 's satisfy a recurrence relation. From the earlier month-by-month calculation we see that the number of adult pairs after n months equals the number of adult pairs after n - 1 months plus the number of baby pairs after n - 1 months. Thus, since the number of baby pairs after n - 1 months equals the number of adult pairs after n - 2 months, we see that the number of adult pairs after n months equals the number of adult pairs after n - 1 months plus the number of adult pairs after n - 2 months.
If we let denote the number of adult pairs after n months, we can state this principle as
(1a) = + , for all n >= 2.
This formula, which is called a recurrence relation, together with the values
(1b) = 0, = 1,
which are called initial values for the sequence , provides a procedure for calculating for any n: = 3, = 5, and so on. We can thus calculate any Fibonacci number, and hence the number of adult pairs of rabbits after any number of months.
It would be useful, however, to have a formula for that depends only on n, and not on the two previous Fibonacci numbers. Such a formula can be obtained by solving the recurrence relation (1a) with initial values (1b), and this can be done with the RSolve command, which is part of the DiscreteMath package. We first load the command and then use it to solve (1a, b).

<<DiscreteMath`RSolve`

sol = RSolve[{F[n] == F[n-1] + F[n-2], F[0] == 0, F[1] == 1}, F[n], n]

fibtemp[n_] = F[n] /. First[sol]

fibtemp[2]

Though fibtemp[n] is the correct formula, it does not present its value in simplified form. We thus modify fibtemp by incorporating the Simplify command using delayed evaluation:

We illustrate the use of fib by computing a few values.

fib[0]
fib[3]
fib[5]

We see that fib produces the same values as the Fibonacci command.

Calculation Via Powers of a Matrix

The recurrence relation (1a, b) can be conveniently expressed in vector-matrix form. Introduce the matrix A = and the vector = , and consider the vector-matrix recurrence relation
(2a) = = A = , for all n >= 1.
Carrying out the matrix multiplication, we see that (2a) reduces to = + and = . The first of these equations is the recurrence relation (1a) and the second is automatically true. So (2a) is just a restatement of (1a). The vector form of the initial values (1b) is
(2b) = = .
Using (2a) repeatedly, together with (2b), we see that = A = A = A(A) = , and in general
(3) =, for all n >= 1.
That is, = is the solution of the vector-matrix recurrence relation (2a, b). So, we can calculate by using (3) to calculate and then in the second component of . Our problem has been reduced to calculating powers of A, and multiplying . Let's carry out this procedure.

A = {{1,1},{1,0}}  

To write A in usual matrix notation, use the MatrixForm command.

MatrixForm[A]

Powers of a matrix are calculated with the MatrixPower command, and not with A^n. We check the online help for MatrixPower, redefine fib using powers of A, and then compute a few Fibonacci numbers.

?MatrixPower

Clear[fib]

The use of Last extracts the second component of the vector . Note that matrix-vector multiplication is indicated by . rather than by *.

fib[1]
fib[5]
fib[30]

Using Eigenpairs to Calculate

Powers of a matrix can be readily calculated if we know the eigenpairs of the matrix. Recall that a pair (&lgr;, x) consisting of a number &lgr; and a vector x is an eigenpair of A if Ax = &lgr;x and x ~= 0. The number &lgr; is called an eigenvalue of A and the vector x is called a corresponding eigenvector. The eigenpairs of A are found with the Eigensystem command.

vals = Eigensystem[A]

The output is a nested list containing first a list of the eigenvalues, and then a list of the corresponding eigenvectors. The eigenvalues are

&lgr; = First[First[vals]]

and

&mgr; = Last[First[vals]]

The corresponding eigenvectors are

x = First[Last[vals]]

and

y = Last[Last[vals]]

The importance of eigenvectors in our situation is that powers of A act on them in a very simple manner. From the definition of eigenpairs we know that Ax = &lgr;x, and from this we learn that x =A(&lgr;x) = &lgr;(Ax) = x. In general, we have x = x. Likewise we have y = y. So a power of A acts on an eigenvector by multiplying the eigenvector by the corresponding eigenvalue raised to the power. But we are interested in = . This can also be calculated, by first expressing the vector as a linear combination of x and y (the eigenvectors are a basis for two dimensional vectors), and applying to this linear combination. To express s a linear combination of x and y we find coefficients a and b so that
= ax + by.
It is easily seen that c = is the solution of the linear system mc = , where

m = {{&lgr;, &mgr;}, {1, 1}} 

So

As a check, let's multiply m by c.

m.c 

Simplify[%]

So = (-1/)x + (1/)y, and thus
= = (-1/) x + (1/)y.
We have reduced the calculation of to the calculation of and . Since is the second component of we see that
(4) = (-1/) + (1/)
= (-1/) + (1/).
We have thus derived the formula produced by RSolve.
We have one further formula. Note that

N[&lgr;]
N[&mgr;]

Since ö&lgr;ö < 1, we see that is small for n large, and hence that the first term in formula (4) is small. Thus is approximately equal to the second term in (4). But, since is an integer, it must be the nearest integer to the second term. Let's compute these terms for n = 3.

N[(-1/Sqrt[5])&lgr;^3]   

N[(1/Sqrt[5])&mgr;^3]

We see that is 2, the closest integer to the second term in (4). The closest integer can be obtained with the Round command.

Clear[fib]

fib[n_] := Round[(1/Sqrt[5])((1+Sqrt[5])/2)^n]

We illustrate the use of this formula by computing a couple of Fibonacci numbers.

fib[30]
fib[100]

For further information on Fibonacci numbers, see V. E. Hoggatt, Jr., Fibonacci and Lucas Numbers, Houghton Mifflin Company, 1969.