Calculating Fibonacci Numbers, Quickly and Exactly

The well-known Fibonacci series \(F_n\) can be defined as follows:

\(F_n =
\begin{cases}
0 & n = 0 \\
1 & n = 1 \\
F_{n-2} + F_{n-1} & n \ge 2\\
\end{cases}\)

Let’s use a few facts about matrices to find a quick way to calculate terms in this famous series.

Lemma

Let \(A = \begin{bmatrix}
0 & 1 \\
1 & 1
\end{bmatrix}\). Then \(A^n = \begin{bmatrix} F_{n-1} & F_n \\ F_n & F_{n+1} \end{bmatrix}\).

Proof

The \(n = 1\) case follows immediately from the definitions of \(A \text{ and } F_n\).

Suppose the statement is true for n. Then
$$\begin{align}
A^{n+1} & = A^n A \\
& = \begin{bmatrix} F_{n-1} & F_n \\ F_n & F_{n+1} \end{bmatrix} \begin{bmatrix} 0 & 1 \\ 1 & 1 \end{bmatrix} \\
& = \begin{bmatrix} F_n & F_{n-1} + F_n \\ F_{n+1} & F_n + F_{n+1} \end{bmatrix} \\
& = \begin{bmatrix} F_n & F_{n+1} \\ F_{n+1} & F_{n+2} \end{bmatrix}
\end{align}$$

And our result follows by induction. QED

Now, notice that

\(\begin{align}
\begin{bmatrix} F_{2n-1} & F_{2n} \\ F_{2n} & F_{2n+1} \end{bmatrix} & = A^{2n} \\
& = A^n A^n \\
& = \begin{bmatrix} F_{n-1} & F_n \\ F_n & F_{n+1} \end{bmatrix} \begin{bmatrix} F_{n-1} & F_n \\ F_n & F_{n+1} \end{bmatrix} \\
& = \begin{bmatrix} F_{n-1}^2 + F_n ^ 2 & F_{n-1} F_n + F_n F_{n+1} \\ F_{n-1} F_n + F_n F_{n+1} & F_n^2 + F_{n+1} ^ 2 \end{bmatrix}
\end{align}\)

Using this identity, we can write \(F_{2n}\) and \(F_{2n+1}\) in terms of \(F_n\) and \(F_{n+1}\).

\(\begin{align}
F_{2n} & = F_{n-1} F_n + F_n F_{n+1} \\
& = (F_{n+1} – F_n) F_n + F_n F_{n+1} \\
& = F_n F_{n+1} – F_n^2 + F_n F_{n+1} \\
& = 2 F_{n+1} F_n – F_n ^ 2
\end{align} \)

\(F_{2n+1} = F_{n}^2 + F_{n+1}^2\)

We now have a way of calculating \(F_{2n}\) and \(F_{2n+1}\) by calculating only a few of the smaller terms in the sequence. This yields some wildly efficient Python code:

def fib2(N):
    if N == 0: return (0, 1)
    half_N, is_N_odd = divmod(N, 2)
    f_n, f_n_plus_1 = fib2(half_N)
    f_n_squared = f_n * f_n
    f_n_plus_1_squared = f_n_plus_1 * f_n_plus_1
    f_2n = 2 * f_n * f_n_plus_1 - f_n_squared
    f_2n_plus_1 = f_n_squared + f_n_plus_1_squared
    if is_N_odd:
        return (f_2n_plus_1, f_2n + f_2n_plus_1)
    return (f_2n, f_2n_plus_1)

def fib(N):
    return fib2(N)[0]

And it’s fast! It’s \(O(\log N)\) in fact:

$ pypy -m timeit -s 'import expfib' 'expfib.fib(500000)'
10 loops, best of 3: 22.4 msec per loop

Compare to the naïve, iterative version, which is \(O(N)\):

def fib2_linear(N):
    f_n, f_n_plus_1 = (0, 1)
    while N > 0:
        N -= 1
        f_n, f_n_plus_1 = f_n_plus_1, f_n + f_n_plus_1
    return (f_n, f_n_plus_1)

def fib_linear(N):
    return fib2_linear(N)[0]

$ pypy -m timeit -s 'import myfib' 'myfib.fib_linear(500000)'
10 loops, best of 3: 2.76 sec per loop

This post was inspired by a post by Lee Phillips and as an excuse to play around with MathJax. The best guide I found for getting started is on Stack Exchange.

About Richard

Richard Kiss was born and raised in Canada. In 1991, he moved to California, ultimately receiving an (or is it "a") M.A. in mathematics. This goes a long way to explaining his obsession with numbers, structure and literalism. His personality is tempered with a heavy dose of wit and opinion. He enjoys donuts and writing in the third person. He also performs stand-up comedy.
This entry was posted in Computers. Bookmark the permalink.

3 Responses to Calculating Fibonacci Numbers, Quickly and Exactly

  1. Great code. I love it!

  2. George says:

    It took me a little effort to convince myself that your final algorithm is logarithmic and not linear in time. A naïve implementation would be linear, but in a more interesting way than usual: the call tree would have a branch factor of two but would now only be log2(n) deep.

    But I see that you take advantage of the two right-hand-side terms always being adjacent to reduce the branch factor to one. Very nice.

Leave a Reply

Your email address will not be published. Required fields are marked *

You may use these HTML tags and attributes: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <strike> <strong>