Applying Fibonacci Fast Doubling Identities

Question

So I sort of understand of how these identities came about from reading this article.

• $F_{2n+1} = F_{n}^2 + F_{n+1}^2$
• $F_{2n} = 2F_{n+1}F_{n}-F_{n}^2 $

But I don't understand how to apply them. I can certainly derive these identities myself and calculate them, but once I have them how do I use them to calculate the nth Fibonacci number? Why is there a branch statement for even and odd? I tried to look up more information, but there wasn't all that much on it with the exception of some academic papers which are written so esoterically and are way beyond me.

Answer 1

Let $\mathbb{Z}[\phi]$ be the ring of integers extended by $\phi = \frac{1 + \sqrt{5}}{2}$. Note that $\phi^2 = 1 + \phi$.

In this ring, $\phi^n = \mathrm{F}_n\phi + \mathrm{F}_{n-1}$, where $\mathrm{F}_n$ is the $n$th Fibonacci number (assuming $\mathrm{F}_0, \mathrm{F}_1, \mathrm{F}_2 = 0, 1, 1$). Proof by induction: $\phi = \mathrm{F}_1\phi + \mathrm{F}_0$, and $\phi^{n+1} = \phi\phi^n = \phi(\mathrm{F}_n\phi + \mathrm{F}_{n-1}) = \mathrm{F}_n\phi^2 + \mathrm{F}_{n-1}\phi = (\mathrm{F}_{n-1} + \mathrm{F}_n)\phi + \mathrm{F}_n = \mathrm{F}_{n+1}\phi + \mathrm{F}_n$.

Comparing the $\phi$ and constant terms at the two ends of this calculation give the two doubling formulae: $\mathrm{F}_{2n}\phi + \mathrm{F}_{2n-1} = \phi^{2n} = (\phi^n)^2 = (\mathrm{F}_n\phi + \mathrm{F}_{n-1})^2 = \mathrm{F}_n^2\phi^2 + 2\mathrm{F}_n\mathrm{F}_{n-1}\phi + \mathrm{F}_{n-1}^2 = (\mathrm{F}_n^2 + 2\mathrm{F}_n\mathrm{F}_{n-1})\phi + \mathrm{F}_n^2 + \mathrm{F}_{n-1}^2$.

But one can finesse the doubling formulae (using them indirectly), by computing $\phi^n$ in $\mathbb{Z}[\phi]$ in code.

First, define a generic exponentiation-by-doubling power function pow(x, n) where n is a positive integer, defined in terms of a multiplication operator mult:

def power(x, n):
    if n == 1: return x
    if n % 2 == 0:
        p = power(x, n//2)
        return mult(p, p)
    return mult(x, power(x, n-1))

The argument n halves at most every two recursive calls, so this performs $O(\mathrm{log}\ n)$ multiplications.

Now one needs to define mult on our ring $\mathbb{Z}[\phi]$. $(a\phi+b)(c\phi+d) = ac\phi^2 + (ad+bc)\phi + bd = (ad+bc+ac)\phi + (bd+ac)$. So, representing $a\phi + b$ as the pair (a, b).

def mult((a, b), (c, d)):
    return (a*d + b*c + a*c, b*d + a*c)

Now we can compute Fibonacci numbers by calculating pow((1, 0), n). The result will be the pair $(\mathrm{F}_n, \mathrm{F}_{n-1})$, so we just output the first item of the tuple.

for n in xrange(1, 21):
    print power((1, 0), n)[0],

The output is the first 20 Fibonacci numbers:

1 1 2 3 5 8 13 21 34 55 89 144 233 377 610 987 1597 2584 4181 6765

As a side note, and somewhat surprisingly, one can partially emulate $\mathbb{Z}[\phi]$ in the integers. Let $X$ be a large number, and express $a\phi + b$ as $aX + b$. Perform calculations modulo $X^2 - X - 1$. When $a$ and $b$ are small relative to $X$, it's easy to check that multiplication works the same as in the true ring. As before $X^n = \mathrm{F}_nX + \mathrm{F}_{n-1}$ (when $n$ is small enough that $\mathrm{F}_n < X - 1$).

Then one can compute Fibonacci numbers efficiently, and with barely any code:

X = 1000000

for i in xrange(1, 21):
    print pow(X, i, X*X-X-1) // X,

(Note that the 3-argument pow(x, n, m) of Python computes $x^n$ modulo $m$). Increasing $X$ allows a larger range of Fibonacci numbers to be calculated.

Answer 2

Suppose you want to calculate $F(n)$ for a given $n$. Here are the key ideas:

• Start with $k=0$.
• Keep track of the pair $(F(k), F(k+1))$.
• Either take one step forward ($k' = k + 1$) or take a doubling step ($k' = 2k$).

Actual example - let's calculate $F(5)$ using these ideas:

• Start with the known $k = 0: (F(k), F(k+1)) = (0, 1)$.
• Take one step forward.
• $k = 1: (F(k), F(k+1)) = (1, 1)$.
• Take a doubling step.
• $k = 2: (F(k), F(k+1)) = (1, 2)$.
• Take a doubling step.
• $k = 4: (F(k), F(k+1)) = (3, 5)$.
• Take one step forward.
• $k = 5: (F(k), F(k+1)) = (5, 8)$.

So now we have found that $F(5) = 5$.

Further reading - math and code: https://www.nayuki.io/page/fast-fibonacci-algorithms

Answer 3

Here is (to my mind) a clean way of both presenting and computing with these identities, although you will have to know how to multiply matrices. (This is one of the most useful things you could possibly learn how to do anyway, so it's well worth doing if you haven't already.)

Start with the matrix $M = \left[ \begin{array}{cc} 1 & 1 \\ 1 & 0 \end{array} \right]$. It's a nice exercise to prove by induction that

$$M^k = \left[ \begin{array}{cc} F_{k+1} & F_k \\ F_k & F_{k-1} \end{array} \right]$$

where my convention is that $F_0 = 0, F_1 = 1$ and hence that $F_{-1} = 1$ (no, really). Then both doubling formulas together are equivalent to the single doubling formula

$$M^{2k} = M^k M^k.$$

To compute $F_{k+1}$ for a very large $k$, the way you use this version of the doubling formula is to compute $M^k$ using binary exponentiation, then read off the top left entry.

Example. Let's compute $F_9$ this way. This is the top left entry of the matrix $M^8$, which we will compute by repeatedly squaring:

$$M^2 = \left[ \begin{array}{cc} 2 & 1 \\ 1 & 1 \end{array} \right]$$ $$M^4 = \left[ \begin{array}{cc} 2 & 1 \\ 1 & 1 \end{array} \right]^2 = \left[ \begin{array}{cc} 5 & 3 \\ 3 & 2 \end{array} \right]$$ $$M^8 = \left[ \begin{array}{cc} 5 & 3 \\ 3 & 2 \end{array} \right]^2 = \left[ \begin{array}{cc} 34 & 21 \\ 21 & 13 \end{array} \right].$$

Hence $F_9 = 34$.

Exercise. Show that this answer is equivalent to Nayuki Minase's answer by first computing that $M^2 = M + 1$ and then showing by induction that this implies

$$M^k = F_k M + F_{k-1}.$$

Hence keeping track of the entries of $M^k$ is equivalent to keeping track of the pair $(F_k, F_{k-1})$.