Lab 02 — Binary Exponentiation and Matrix Exponentiation for Fibonacci

Goal

Master O(log b) exponentiation in two settings: integer fast power (a^b mod p) and matrix fast power (Fibonacci F(n) mod p for n up to 10^18). By the end, you can write integer fast power in <2 minutes and matrix Fibonacci in <10 minutes from blank.

Background

Binary exponentiation is the engine behind almost every “compute a^b for huge b” subroutine in CP and cryptography. The same divide-and-conquer pattern — a^b = (a^(b/2))^2 if b even, else a · a^(b-1) — generalizes from integers to any associative operation: matrix multiplication (linear recurrences), polynomial multiplication (signal processing), function composition (rotation by an angle, repeated f(f(f(...)))). Internalize the O(log b) skeleton once and you get all these for free.

Matrix exponentiation is the canonical extension. Fibonacci, defined F(n) = F(n-1) + F(n-2), has the matrix form [F(n+1), F(n)]^T = M · [F(n), F(n-1)]^T where M = [[1, 1], [1, 0]]. Therefore [F(n+1), F(n)]^T = M^n · [F(1), F(0)]^T = M^n · [1, 0]^T. Computing M^n takes O(log n) matrix multiplications, each O(2^3) = O(8) operations. Total: O(log n) for n = 10^18, about 60 multiplications.

Interview Context

Quant interviews use exactly this pair. “Compute 2^(10^18) mod (10^9 + 7)” is a 2-minute warm-up. “Compute F(10^18) mod (10^9 + 7)” is the follow-up that filters out candidates who only know the iterative O(N) Fibonacci. A candidate who reaches for matrix exponentiation reflexively when seeing n ≤ 10^18 and a linear recurrence is signalling 1900+ CF rating, which is a strong positive signal at HFT firms.

Problem Statement

Part 1. Implement power(a, b, m) returning a^b mod m for 0 ≤ a < m, 0 ≤ b ≤ 10^18, 1 ≤ m ≤ 10^9 + 7.

Part 2. Implement fib(n, p) returning F(n) mod p for 0 ≤ n ≤ 10^18, p = 10^9 + 7, where F(0) = 0, F(1) = 1, F(n) = F(n-1) + F(n-2).

LeetCode reference: LC 50 (Pow(x, n)) — Part 1 in real-number form. LC 1137 (N-th Tribonacci) — Part 2 with three-term recurrence (analogous matrix form).

Constraints

Part 1: b up to 10^18, so naive O(b) looping is impossible.
Part 2: n up to 10^18, so naive O(n) Fibonacci is impossible.
Time limit: 2 seconds.
Memory: O(1) for Part 1, O(1) for Part 2 (matrices are 2×2).

Clarifying Questions

“Negative b for integer power?” — for modular a^b, undefined unless gcd(a, m) = 1 and you want a^(-1)^|b|. For real a^b (LC 50), return 1 / a^|b| and watch for INT_MIN overflow on -b.
“What’s 0^0?” — by convention 1.
“Fibonacci indexing — is F(1) = 1 or F(2) = 1?” — confirm; standard CP convention is F(0) = 0, F(1) = 1.
“Matrix exponentiation modulus — same p everywhere?” — yes.
“Required output format for matrices?” — only the scalar Fibonacci value, but you might be asked to return the full state vector.

Examples

power(2, 10, 1000) = 24 (2^10 = 1024).
power(3, 0, 7) = 1.
power(5, 1000000000000000000, 10^9 + 7) = some specific value (must compute).
fib(0, p) = 0, fib(1, p) = 1, fib(10, p) = 55.
fib(10^18, 10^9 + 7) = a specific value the engine must produce.

Initial Brute Force

Part 1: result = 1; for i in 1..b: result = result * a % m. O(b). Part 2: iterative two-variable Fibonacci. O(n).

def power_brute(a, b, m):
    result = 1 % m
    for _ in range(b):
        result = result * a % m
    return result

def fib_brute(n, p):
    if n == 0: return 0
    a, b = 0, 1
    for _ in range(n - 1):
        a, b = b, (a + b) % p
    return b

Brute Force Complexity

Part 1: O(b) mults, b = 10^18 is 10^17 mults/sec required — impossible. Part 2: same. Useful only as oracles on n ≤ 30.

Optimization Path

Part 1.

Naive O(b) (above).
Recursive O(log b): power(a, b) = power(a, b/2)² if b even, else a · power(a, b-1).
Iterative O(log b): process bits of b low-to-high, square a each iteration, multiply into result when the current bit of b is 1. Preferred for stack safety.

Part 2.

Naive O(n) (above).
Memoized doubling: F(2k) = F(k) · (2 · F(k+1) − F(k)), F(2k+1) = F(k)² + F(k+1)². O(log n) recursion. Beautiful but error-prone.
Matrix exponentiation: build M = [[1,1],[1,0]], compute M^n via integer-fast-power lifted to matrices, extract M^n[0][1] as F(n). O(log n) · O(8) = ~500 ops for n = 10^18. The general-purpose technique that works for any linear recurrence.

Final Expected Approach

Part 1 (iterative):

long long power(long long a, long long b, long long m) {
    long long res = 1 % m;
    a %= m;
    while (b > 0) {
        if (b & 1) res = res * a % m;
        a = a * a % m;
        b >>= 1;
    }
    return res;
}

Part 2 (matrix):

typedef vector<vector<long long>> Mat;
const long long P = 1e9 + 7;

Mat matmul(const Mat &A, const Mat &B) {
    int n = A.size();
    Mat C(n, vector<long long>(n, 0));
    for (int i = 0; i < n; ++i)
        for (int k = 0; k < n; ++k)
            if (A[i][k])
                for (int j = 0; j < n; ++j)
                    C[i][j] = (C[i][j] + A[i][k] * B[k][j]) % P;
    return C;
}

Mat matpow(Mat M, long long e) {
    int n = M.size();
    Mat res(n, vector<long long>(n, 0));
    for (int i = 0; i < n; ++i) res[i][i] = 1;  // identity
    while (e > 0) {
        if (e & 1) res = matmul(res, M);
        M = matmul(M, M);
        e >>= 1;
    }
    return res;
}

long long fib(long long n) {
    if (n == 0) return 0;
    Mat M = {{1, 1}, {1, 0}};
    Mat R = matpow(M, n);
    return R[0][1];
}

Data Structures Used

Part 1: scalars only.
Part 2: 2×2 matrices as vector<vector<long long>> or array<array<long long, 2>, 2>.

Correctness Argument

Part 1 (binary exponentiation invariant). Let b = sum b_i 2^i be the binary representation of the original exponent. After k iterations, the variable a equals a_initial^(2^k) mod m, and res equals the product of a_initial^(2^i) mod m for all i < k with b_i = 1. After all iterations, res = a_initial^b mod m. The invariant proves correctness; termination is b → 0 after floor(log₂ b) + 1 iterations.

Part 2 (matrix recurrence). Define the column vector v(k) = [F(k+1), F(k)]^T. Then M · v(k) = [[1,1],[1,0]] · [F(k+1), F(k)]^T = [F(k+1) + F(k), F(k+1)]^T = [F(k+2), F(k+1)]^T = v(k+1). By induction, v(n) = M^n · v(0) = M^n · [1, 0]^T, so F(n) = v(n)[1] = (M^n · [1, 0]^T)[1] = M^n[1][0]. Equivalently, M^n[0][1] = F(n) (by symmetry of the Fibonacci matrix). Matrix multiplication is associative, so binary exponentiation lifts directly: same invariant, just with matrices.

Complexity

Part 1: O(log b) time, O(1) space.
Part 2: O(K^3 log n) time for K × K matrices (K = 2 for Fibonacci → 8 log n ≈ 480 ops for n = 10^18), O(K^2) space.

Implementation Requirements

All multiplications mod p immediately.
Identity matrix initialization for matpow.
Handle n = 0 (return 0 directly, not M^0 · [1, 0]^T = [1, 0]^T which would give F(0) = 0 correctly anyway, but watch corner cases).
Use long long (or int64); intermediate products are up to (10^9)² ≈ 10^18, just fitting long long.

Tests

def test_power():
    assert power(2, 10, 1000) == 24
    assert power(3, 0, 7) == 1
    assert power(0, 0, 7) == 1
    assert power(0, 5, 7) == 0
    assert power(2, 63, 10**18) == 2**63 % 10**18

def test_fib():
    p = 10**9 + 7
    expected = [0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55]
    for i, e in enumerate(expected):
        assert fib(i, p) == e
    # Stress vs brute up to N = 30
    for n in range(31):
        assert fib(n, p) == fib_brute(n, p)
    # Spot-check large
    assert fib(10**18, p) == 209783453   # known value, mod 1e9+7

Follow-up Questions

“Compute F(n) for general linear recurrence f(n) = c1 f(n-1) + ... + ck f(n-k)?” — same matrix exponentiation, with a k×k companion matrix.
“Compute the number of paths of length L between two nodes in a graph?” — (adj_matrix)^L, indexed by start/end. Same O(V^3 log L) engine.
“Compute Tribonacci T(n) in O(log n)?” — 3×3 companion matrix, otherwise identical.
“Avoid O(K^3) per multiplication for huge K?” — research topics: Kitamasa’s algorithm reduces to O(K^2 log N); FFT-based polynomial multiplication reduces further. Out of scope here.

Product Extension

Cryptography (RSA encryption: compute m^e mod n with e, n 2048-bit, in milliseconds — same algorithm, with bignum). Computer graphics (rotation matrix R^n for repeated rotations). Markov chain steady-state approximation (P^n for stochastic matrix P, large n). Network reachability (adj^L for paths of length L).

Language / Runtime Follow-ups

Python: built-in pow(a, b, m) is O(log b) and ~10× faster than a hand-rolled Python loop. For matrix power: hand-roll the multiplication; numpy is overkill at K = 2 (FFI overhead exceeds compute).
Java: BigInteger.modPow exists but is 10× slower than a hand-rolled long loop. Use the hand-rolled version unless values exceed long.
Go: math/big.Int.Exp(a, b, m) is correct but slow; hand-roll for hot paths. Matrix: use 2D [][]int64 arrays, not math/big.
C++: __int128 for intermediate products if m exceeds ~3·10^9 (where (long long)² mod p overflows). For standard m = 10^9 + 7, plain long long suffices.
JS/TS: BigInt is correct but slow; for m < 2^32, use Number carefully (Math.floor after * a / m) — easy to get wrong. Matrix: same caveat.

Common Bugs

Forgetting res = 1 % m at start (returns 1 instead of 0 when m = 1).
Squaring a before checking the bit (results in extra unused multiplication; correctness OK, perf hit).
Matrix multiplication order: M^n = M · M · ... · M, but matmul is non-commutative — be deliberate about left/right.
Using int instead of long long for matrix entries; (10^9)² > 2^31.
Recursive power blowing the stack on b = 10^18 (recursion depth 60 is fine; just don’t go deeper).

Debugging Strategy

If power(a, b, m) is wrong:

Test on small cases (b ≤ 5) where you can verify by hand.
Print binary representation of b and confirm the bits are processed.
If correct on small b but wrong on large: overflow. Check that a * a % m uses long long, not int.

If fib(n, p) is wrong:

Verify M^1 = M, M^2 = [[2,1],[1,1]], M^3 = [[3,2],[2,1]]. Each entry is a Fibonacci number.
Verify fib(0..10) matches the canonical sequence.
If small n works but large doesn’t: same overflow check.

Mastery Criteria

Write power(a, b, m) from blank in <2 minutes.
Write the Fibonacci matrix-exponentiation engine from blank in <10 minutes.
Articulate the binary-exponentiation invariant in one sentence.
Articulate the matrix-Fibonacci recurrence (v(n+1) = M · v(n)) in 30 seconds.
Generalize to any linear recurrence: given the recurrence, write the companion matrix in <2 minutes.
Recognize “n up to 10^18 + linear recurrence” as the trigger for matrix exponentiation in <60 seconds.

← Lab 01 — Modular Arithmetic · Phase 7 README · Lab 03 — Sieve and Factorization →

LeetCode Interview — Extreme Coding