TAOCP 4.6.3 Exercise 26
We wish to compute the $n$th Fibonacci number $F_n$ modulo $m$, where $n$ and $m$ are large integers.
Section 4.6.3: Evaluation of Powers
Exercise 26. ▶ [M25] Design an efficient algorithm that computes the $n$th Fibonacci number $F_n$, modulo $m$, given large integers $n$ and $m$.
Verified: yes
Solve time: 1m31s
Solution
We wish to compute the $n$th Fibonacci number $F_n$ modulo $m$, where $n$ and $m$ are large integers. A naive approach, iterating the recurrence $F_0 = 0$, $F_1 = 1$, $F_k = F_{k-1} + F_{k-2}$, requires $O(n)$ additions, which is inefficient for large $n$. We derive an algorithm of time $O(\log n)$ using exponentiation by squaring on a $2 \times 2$ matrix.
Consider the Fibonacci Q-matrix
$$ Q = \begin{pmatrix}1 & 1 \ 1 & 0\end{pmatrix}. $$
It is well known that
$$ Q^n = \begin{pmatrix}F_{n+1} & F_n \ F_n & F_{n-1}\end{pmatrix}. $$
Hence, $F_n = (Q^n)_{1,2}$, the entry in the first row, second column. Computation modulo $m$ is preserved by matrix multiplication modulo $m$; that is, if all multiplications and additions are performed modulo $m$, the resulting $F_n \bmod m$ is exact.
We define the following right-to-left binary method for exponentiation of $2 \times 2$ matrices modulo $m$:
Let $N \leftarrow n$, $Y \leftarrow I$, $Z \leftarrow Q$, where $I$ is the $2 \times 2$ identity matrix. At each step, let $t \leftarrow N \bmod 2$ and $N \leftarrow \lfloor N/2 \rfloor$. If $t = 1$, set $Y \leftarrow Y \cdot Z \bmod m$. Then set $Z \leftarrow Z \cdot Z \bmod m$. Repeat until $N = 0$. Finally, return $F_n \equiv (Y)_{1,2} \bmod m$.
To justify correctness, note that the matrix recurrence preserves the identity $Y Z^N \equiv Q^n \pmod m$ at each step, analogous to Algorithm A in Section 4.6.3. When $N = 0$, $Y = Q^n \bmod m$, and the desired entry $(Y)_{1,2}$ equals $F_n \bmod m$. Each step involves at most one multiplication of $2 \times 2$ matrices modulo $m$, and the total number of steps is $\lfloor \lg n \rfloor + \nu(n)$, where $\nu(n)$ is the number of ones in the binary representation of $n$.
Explicitly, the matrix multiplication modulo $m$ is
$$ \begin{pmatrix}a & b \ c & d\end{pmatrix} \cdot \begin{pmatrix}e & f \ g & h\end{pmatrix} \equiv \begin{pmatrix} (ae + bg) \bmod m & (af + bh) \bmod m \ (ce + dg) \bmod m & (cf + dh) \bmod m \end{pmatrix} \bmod m. $$
This completes the algorithm. The asymptotic cost is $O(\log n)$ matrix multiplications modulo $m$, each involving a constant number of arithmetic operations. The method is efficient even for very large $n$ and $m$, since the size of the numbers involved remains bounded by $m$ at every step.
The algorithm can be summarized in step form, analogous to Algorithm A:
F1. [Initialize.] Set $N \leftarrow n$, $Y \leftarrow I$, $Z \leftarrow Q$.
F2. [Halve $N$.] Set $t \leftarrow N \bmod 2$ and $N \leftarrow \lfloor N/2 \rfloor$. If $t = 1$, go to F3. Otherwise, skip to F4.
F3. [Multiply $Y$ by $Z$.] Set $Y \leftarrow Y \cdot Z \bmod m$.
F4. [$N = 0$?] If $N = 0$, terminate, returning $(Y)_{1,2}$ as $F_n \bmod m$.
F5. [Square $Z$.] Set $Z \leftarrow Z \cdot Z \bmod m$ and return to F2. $\blacksquare$
This procedure achieves the desired result efficiently and correctly for all integers $n \ge 0$ and $m \ge 1$.