Computing huge Fibonacci numbers

In this post, we will show one rather efficient technique for computing large Fibonacci numbers. The resulting Java implementation computes $F_{1 000 000}$ in less than two seconds on a 2.5 GHz CPU. The key idea is that

$\displaystyle A^n = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^n = \begin{pmatrix} F_{n + 1} & F_n \\ F_n & F_{n - 1} \end{pmatrix},$

where

$\displaystyle F_n = \begin{cases} 0 & \text{if } n = 0, \\ 1 & \text{if } n = 1, \\ F_{n - 1} + F_{n - 2} & \text{if } n > 1. \end{cases}$

Let us use induction in order to prove the matrix equation.

The base case

$\displaystyle \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^1 = \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} = \begin{pmatrix} F_2 & F_1 \\ F_1 & F_0 \end{pmatrix}.$

Inductive step

Suppose that

$\displaystyle \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^n = \begin{pmatrix} F_{n + 1} & F_n \\ F_n & F_{n - 1} \end{pmatrix}.$

Now we have that

\displaystyle \begin{aligned} \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^{n + 1} &= \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix}^n \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \\ &= \begin{pmatrix} F_{n + 1} & F_n \\ F_n & F_{n - 1} \end{pmatrix} \begin{pmatrix} 1 & 1 \\ 1 & 0 \end{pmatrix} \\ &= \begin{pmatrix} F_{n + 1} + F_n & F_{n + 1} \\ F_n + F_{n - 1} & F_n \end{pmatrix} \\ &= \begin{pmatrix} F_{n + 2} & F_{n + 1} \\ F_{n + 1} & F_n \end{pmatrix}, \end{aligned}

which concludes the proof. Since we want to deal with large Fibonacci numbers, we have to opt to java.math.BigInteger. Using BigInteger not just allows dealing with integers that do not fit in primitive integer types, but also provides more efficient methods for multiplying two integers: as of now, JDK 8 implementation of BigInteger relies on two distinct multiplication algorithms

• Toom-Cook algorithm that runs in time $\Theta(n^{1.465}),$
• Karatsuba algorithm that runs in time $\mathcal{O}(n^{1.585}).$

(The naïve multiplication runs in $\Theta(n^2)$.) Another optimization is to compute the matrix exponent by divide-and-conquer technique that requires $\Theta(\log n)$ actual multiplications: basically, if $n$ is even, compute recursively $B = A^{n / 2}$ and return $B^2$; if $n$ is odd, compute recursively $C = A^{n - 1}$ and return $AC$. (The recursion is terminated at $n = 1$, where we return the input matrix unmodified.) Since $F_n = \Theta(c^n)$ for some $c > 1$, the length of the integer is $\Theta(\log F_n) = \Theta(\log c^n) = \Theta(n)$. Also, the worst case complexity of multiplying two BigIntegers of length $n$ is $\mathcal{O}(n^{1.585}).$ As we need logarithmic amount of matrix multiplications and each matrix multiplication involves a constant number of BigInteger multiplications, we have the upper bound of $\mathcal{O}(n^{1.585} \log n)$ for computing $F_n$.

Source code

import java.math.BigInteger;

public class LargeFibonacciNumbers {

public static BigInteger fibonacci(int n) {
if (n < 0) {
throw new IllegalArgumentException(
"The (" + n + ")th Fibonacci number is not defined.");
}

switch (n) {
case 0:
return BigInteger.ZERO;

case 1:
return BigInteger.ONE;
}

BigInteger[][] matrix = new BigInteger[2][2];

matrix[0][0] = BigInteger.ONE;
matrix[0][1] = BigInteger.ONE;
matrix[1][0] = BigInteger.ONE;
matrix[1][1] = BigInteger.ZERO;

return fibonacciMatrix(matrix, n)[0][1];
}

private static BigInteger[][] multiply(BigInteger[][] matrix1,
BigInteger[][] matrix2) {
BigInteger[][] ret = new BigInteger[2][2];

ret[0][0] = matrix1[0][0].multiply(matrix2[0][0])

ret[0][1] = matrix1[0][0].multiply(matrix2[0][1])

ret[1][0] = matrix1[1][0].multiply(matrix2[0][0])

ret[1][1] = matrix1[1][0].multiply(matrix2[0][1])

return ret;
}

private static BigInteger[][]
fibonacciMatrix(BigInteger[][] matrix, int n) {
if (n == 1) {
// End the recursion.
return matrix;
}

if ((n & 1) == 1) {
// 'n' is odd.
return multiply(matrix, fibonacciMatrix(matrix, n - 1));
}

// 'n' is even.
BigInteger[][] tmp = fibonacciMatrix(matrix, n >> 1);
return multiply(tmp, tmp);
}

public static void main(String[] args) {
int n = -1;

if (args.length > 0) {
try {
n = Integer.parseInt(args[0]);
} catch (NumberFormatException ex) {
System.err.println("Could not parse \"" + args[0] +
"\" as an integer.");
return;
}

try {
System.out.println(fibonacci(n));
} catch (IllegalArgumentException ex) {
System.err.println(ex.getMessage());
}
} else {
System.out.println("Usage: java -jar File.jar N");
System.out.println("Where N is the index of the Fibonacci number " +
"to compute.");
}
}
}