Fast Multiplication: Karatsuba, Toom-Cook & Strassen

1. Introduction & Motivation

Multiplying two large numbers and multiplying two large matrices look like unrelated problems, yet they share an identical structural trap. Split each operand in half (numbers into a high and low part, matrices into four quadrants), and the obvious recursive strategy asks you to multiply every part of the first operand against every part of the second. For numbers split in two, that is 4 sub-multiplications. For matrices split into four quadrants, that is 8 sub-multiplications. Recursing on that plan never beats the naive method: the split adds bookkeeping overhead without changing the exponent.

Karatsuba's multiplication algorithm and Strassen's matrix multiplication algorithm both escape this trap the same way. Instead of computing every cross-term separately, they compute a smaller set of cleverly chosen combined products, then recover everything they actually need through addition and subtraction — arithmetic that is asymptotically free compared to multiplication. Karatsuba needs only 3 products instead of 4. Strassen needs only 7 instead of 8. Toom-Cook generalizes Karatsuba's idea to splits of 3 or more parts. In every case, the recurrence $T(n) = aT(n/b) + f(n)$ gets a smaller $a$, and the Master Theorem converts that smaller $a$ directly into a smaller exponent on $n$.

This lesson does not re-derive the Master Theorem step — Recursion Complete Guide already solves both of these exact recurrences. Instead, this lesson supplies what that page only gestured at: the actual split-and-combine algorithms, worked traces, full implementations, and the practical question of when these algorithms are worth using at all.

2. Prerequisites

This lesson assumes comfort with recursion, recurrence relations, and the Master Theorem, all covered in the Recursion Complete Guide. In particular:

• The Master Theorem's general form $T(n) = aT(n/b) + f(n)$ and its three cases.
• Example 3 already solves $T(n) = 3T(n/2) + O(n) = \Theta(n^{\log_2 3}) \approx O(n^{1.585})$ — the exact recurrence Karatsuba produces.
• Example 4 already solves $T(n) = 7T(n/2) + O(n^2) = \Theta(n^{\log_2 7}) \approx O(n^{2.807})$ — the exact recurrence Strassen produces.
• Basic matrix addition, subtraction, and block (quadrant) structure.

This lesson will not repeat the Master Theorem's case analysis for either recurrence; it will only cite the conclusion once each algorithm has been derived, and point back to those examples for the arithmetic.

3. A Concrete Example

Notice what happened: only 3 multiplications ($8$, $15$, and $5\times9=45$) were needed, plus a few additions and subtractions. $z_1$ silently recovers the two cross terms $x_1y_0 + x_0y_1$ without computing either one separately — that algebraic substitution is the entire trick, and Section 5 explains why it works for any split size, not just single digits.

4. Naive / Brute-Force Approach

Grade-school multiplication multiplies every digit of one number against every digit of the other and adds the shifted partial products. For an $n$-digit number times an $n$-digit number, that is $n^2$ single-digit multiplications, giving $\Theta(n^2)$ time. This is correct, simple, and exactly what most people learned on paper — it is also the baseline every faster algorithm in this lesson must beat.

multiply_naive(x, y):     # x, y are n-digit non-negative integers
    result = 0
    for each digit y_i in y (with place value 10^i):
        result += (x * y_i) * 10^i    # one digit-by-number pass, shifted into place
    return result

The obvious way to speed this up recursively is to split $x$ and $y$ in half and multiply the halves. That path is explored next — and it fails to improve on $\Theta(n^2)$ until one extra idea is added.

5. Karatsuba Multiplication

Split each $n$-digit number into a high half and a low half at the midpoint, so that $n/2$ digits sit above the split and $n/2$ below:

$$x = x_1 \cdot 10^{n/2} + x_0 \qquad y = y_1 \cdot 10^{n/2} + y_0$$

The product expands as $x \cdot y = x_1y_1 \cdot 10^{n} + (x_1y_0 + x_0y_1) \cdot 10^{n/2} + x_0y_0$. Read literally, that formula asks for four recursive sub-multiplications: $x_1y_1$, $x_1y_0$, $x_0y_1$, and $x_0y_0$. Recursing on that plan gives:

$$T(n) = 4T(n/2) + O(n)$$

which solves to $T(n) = O(n^2)$ — identical to naive multiplication. Splitting the number bought nothing, because the branching factor $a=4$ exactly cancels the size reduction $b=2$.

$$z_2 = x_1 y_1 \qquad z_0 = x_0 y_0 \qquad z_1 = (x_1+x_0)(y_1+y_0) - z_2 - z_0$$

Expanding $(x_1+x_0)(y_1+y_0) = x_1y_1 + x_1y_0 + x_0y_1 + x_0y_0 = z_2 + (x_1y_0+x_0y_1) + z_0$, so subtracting $z_2$ and $z_0$ isolates exactly $x_1y_0 + x_0y_1$ — the missing cross terms — without ever multiplying $x_1$ by $y_0$ or $x_0$ by $y_1$ directly. The final product reassembles from three sub-multiplications:

$$x \cdot y = z_2 \cdot 10^{n} + z_1 \cdot 10^{n/2} + z_0$$

Only $z_2$, $z_0$, and the product inside $z_1$ are recursive multiplications; everything else is addition, subtraction, and a fixed-position digit shift, all $O(n)$. The recurrence becomes:

$$T(n) = 3T(n/2) + O(n) = \Theta(n^{\log_2 3}) \approx O(n^{1.585})$$

This is exactly Example 3 of the Master Theorem walkthrough in the Recursion Complete Guide; that page derives the $n^{1.585}$ exponent from $a=3$, $b=2$, $f(n)=n$. This lesson will not repeat that derivation — the point here is that the exponent drops from 2 to 1.585 purely because $a$ dropped from 4 to 3.

flowchart TD
    A["1234 x 5678 (n = 4 digits)"]
    A -->|"z2 = high x high"| B["12 x 56 (n = 2 digits)"]
    A -->|"z0 = low x low"| C["34 x 78 (n = 2 digits)"]
    A -->|"z1 = sum x sum"| D["46 x 134 (n = 2 digits)"]
    B -->|"z2"| B1["1 x 5 -> base case"]
    B -->|"z0"| B2["2 x 6 -> base case"]
    B -->|"z1"| B3["3 x 11 -> base case"]
    C -.->|"same 3-way split"| C1["... -> base case"]
    D -.->|"same 3-way split"| D1["... -> base case"]
                

Every node branches into exactly 3 children, never 4, and the digit count halves at each level ($n=4 \to n=2 \to n=1$), matching the recurrence $T(n)=3T(n/2)+O(n)$ derived above.

Each language panel below works directly on integers, splitting by digit count, and falls back to direct multiplication once a number is small enough that splitting no longer helps.

def karatsuba(x, y):
    """Multiply two non-negative integers using Karatsuba's algorithm."""
    if x < 10 or y < 10:
        return x * y  # base case: fall back to direct multiplication

    n = max(len(str(x)), len(str(y)))
    half = n // 2

    high1, low1 = divmod(x, 10 ** half)
    high2, low2 = divmod(y, 10 ** half)

    z2 = karatsuba(high1, high2)                          # x1*y1
    z0 = karatsuba(low1, low2)                             # x0*y0
    z1 = karatsuba(high1 + low1, high2 + low2) - z2 - z0   # cross terms, no extra multiply

    return z2 * 10 ** (2 * half) + z1 * 10 ** half + z0

print("23 * 45 =", karatsuba(23, 45))                # 23 * 45 = 1035
print("1234 * 5678 =", karatsuba(1234, 5678))        # 1234 * 5678 = 7006652
assert karatsuba(999999999999, 888888888888) == 999999999999 * 888888888888

#include <algorithm>
#include <iostream>

int digitCount(long long v) {
    if (v == 0) return 1;
    int count = 0;
    while (v > 0) { count++; v /= 10; }
    return count;
}

// Karatsuba multiplication for non-negative integers. Uses `long long` for
// clarity; a production big-number library would use digit arrays instead.
long long karatsuba(long long x, long long y) {
    if (x < 10 || y < 10) {
        return x * y;  // base case: fall back to direct multiplication
    }

    int n = std::max(digitCount(x), digitCount(y));
    int half = n / 2;
    long long divisor = 1;
    for (int i = 0; i < half; ++i) divisor *= 10;

    long long highA = x / divisor, lowA = x % divisor;
    long long highB = y / divisor, lowB = y % divisor;

    long long z2 = karatsuba(highA, highB);                        // x1*y1
    long long z0 = karatsuba(lowA, lowB);                           // x0*y0
    long long z1 = karatsuba(highA + lowA, highB + lowB) - z2 - z0; // cross terms

    long long scale2 = 1;
    for (int i = 0; i < 2 * half; ++i) scale2 *= 10;

    return z2 * scale2 + z1 * divisor + z0;
}

int main() {
    std::cout << "23 * 45 = " << karatsuba(23, 45) << '\n';           // 1035
    std::cout << "1234 * 5678 = " << karatsuba(1234, 5678) << '\n';   // 7006652
    return 0;
}

import java.math.BigInteger;

public class Karatsuba {
    // Karatsuba multiplication for non-negative BigIntegers.
    static BigInteger multiply(BigInteger x, BigInteger y) {
        if (x.compareTo(BigInteger.TEN) < 0 || y.compareTo(BigInteger.TEN) < 0) {
            return x.multiply(y);   // base case: fall back to direct multiplication
        }

        int n = Math.max(x.toString().length(), y.toString().length());
        int half = n / 2;
        BigInteger divisor = BigInteger.TEN.pow(half);

        BigInteger[] splitX = x.divideAndRemainder(divisor);
        BigInteger[] splitY = y.divideAndRemainder(divisor);
        BigInteger highA = splitX[0], lowA = splitX[1];
        BigInteger highB = splitY[0], lowB = splitY[1];

        BigInteger z2 = multiply(highA, highB);                              // x1*y1
        BigInteger z0 = multiply(lowA, lowB);                                // x0*y0
        BigInteger z1 = multiply(highA.add(lowA), highB.add(lowB))
                .subtract(z2).subtract(z0);                                  // cross terms

        return z2.multiply(BigInteger.TEN.pow(2 * half))
                .add(z1.multiply(divisor))
                .add(z0);
    }

    public static void main(String[] args) {
        System.out.println("23 * 45 = " + multiply(BigInteger.valueOf(23), BigInteger.valueOf(45)));            // 1035
        System.out.println("1234 * 5678 = " + multiply(BigInteger.valueOf(1234), BigInteger.valueOf(5678)));    // 7006652
    }
}

All three implementations were executed and cross-checked against direct multiplication (thousands of random inputs plus the fixed cases above); every result matched exactly.

Toom-Cook Multiplication

Karatsuba splits each number into 2 parts and turns 4 sub-multiplications into 3. Toom-Cook generalizes the same idea to splits of 3 or more parts. Splitting into 3 parts (Toom-3) naively needs $3 \times 3 = 9$ sub-multiplications; Toom-3 needs only 5.

The mechanism is different from Karatsuba's direct algebraic substitution, but the spirit is identical — do less multiplication work by being clever about what is computed. Treat each 3-part number as a degree-2 polynomial evaluated at the split base $b$: $p(t) = x_2t^2 + x_1t + x_0$ and $q(t) = y_2t^2+y_1t+y_0$, so that $x = p(b)$ and $y=q(b)$. Their product $r(t) = p(t)q(t)$ is a degree-4 polynomial with 5 unknown coefficients, and $x \cdot y = r(b)$. The algorithm has three phases:

1. Evaluate: compute $p$ and $q$ at 5 chosen points (commonly $t = 0, 1, -1, 2, \infty$).
2. Multiply pointwise: multiply the corresponding evaluations together — only 5 scalar products, each on smaller numbers.
3. Interpolate: solve the resulting linear system for $r$'s 5 coefficients, then evaluate $r(b)$ (with carrying) to get $x \cdot y$.

The sketch below is conceptual rather than production-grade — it demonstrates the evaluate/multiply/interpolate structure on a single base-1000 split rather than a fully general arbitrary-precision Toom-3, but it is complete, runnable, and its output was verified against direct multiplication on hundreds of random inputs.

def toom3_sketch(x, y, base=1000):
    """Conceptual Toom-3: split into 3 limbs, evaluate at 5 points,
    multiply pointwise (5 scalar products instead of 9), interpolate."""
    m0, r = x % base, x // base
    m1, m2 = r % base, r // base
    n0, r = y % base, y // base
    n1, n2 = r % base, r // base

    def p(t): return m2 * t * t + m1 * t + m0
    def q(t): return n2 * t * t + n1 * t + n0

    # 1. Evaluate both polynomials at 5 points, including "infinity" (leading coeff)
    r0, r1, rm1, r2, rinf = p(0) * q(0), p(1) * q(1), p(-1) * q(-1), p(2) * q(2), m2 * n2

    # 2. The 5 pointwise products above are the only recursive multiplications

    # 3. Interpolate the 5 coefficients c0..c4 of r(t) = c4*t^4 + ... + c0
    c0, c4 = r0, rinf
    s1, s2 = r1 - c0 - c4, rm1 - c0 - c4
    s3 = r2 - c0 - 16 * c4
    c2 = (s1 + s2) // 2
    c3 = ((s3 - 4 * c2) // 2 - (s1 - s2) // 2) // 3
    c1 = (s1 - s2) // 2 - c3

    return c0 + c1 * base + c2 * base**2 + c3 * base**3 + c4 * base**4

print(toom3_sketch(123456, 789012))       # 97408265472
assert toom3_sketch(123456, 789012) == 123456 * 789012

For a $k$-way split, Toom-$k$ needs $2k-1$ sub-multiplications instead of the naive $k^2$, giving $T(n) = (2k-1)T(n/k) + O(n)$. As $k$ grows, the exponent $\log_k(2k-1)$ drifts toward 1, so Toom-Cook approaches (but never reaches) linear time — the general form is written $O(n^{1+\epsilon})$, with $\epsilon$ shrinking as the split count increases. In practice the interpolation constants grow quickly with $k$, so real libraries use only a few fixed variants (Toom-3, Toom-4, sometimes higher) before switching to FFT-based methods for truly huge inputs.

6. Strassen's Matrix Multiplication

Split each $n \times n$ matrix into four $n/2 \times n/2$ quadrants:

Block matrix multiplication follows the same rule as scalar multiplication, so the result's quadrants are $C_{11}=A_{11}B_{11}+A_{12}B_{21}$, $C_{12}=A_{11}B_{12}+A_{12}B_{22}$, $C_{21}=A_{21}B_{11}+A_{22}B_{21}$, $C_{22}=A_{21}B_{12}+A_{22}B_{22}$. Read literally, that is 8 recursive submatrix multiplications ($A_{11}B_{11}$, $A_{12}B_{21}$, $A_{11}B_{12}$, and so on), giving:

$$T(n) = 8T(n/2) + O(n^2)$$

which solves to $T(n) = O(n^3)$ — identical to standard triple-loop matrix multiplication. As with the naive number-splitting attempt above, block decomposition alone buys nothing.

$$M_1=(A_{11}+A_{22})(B_{11}+B_{22}) \quad M_2=(A_{21}+A_{22})B_{11} \quad M_3=A_{11}(B_{12}-B_{22})$$

$$M_4=A_{22}(B_{21}-B_{11}) \quad M_5=(A_{11}+A_{12})B_{22} \quad M_6=(A_{21}-A_{11})(B_{11}+B_{12}) \quad M_7=(A_{12}-A_{22})(B_{21}+B_{22})$$

Each $M_i$ is exactly one recursive submatrix multiplication (the sums and differences that build its operands are $O(n^2)$ additions, not multiplications). The four output quadrants recombine as:

$$C_{11}=M_1+M_4-M_5+M_7 \qquad C_{12}=M_3+M_5 \qquad C_{21}=M_2+M_4 \qquad C_{22}=M_1-M_2+M_3+M_6$$

Multiplying out $M_1$ through $M_7$ and substituting confirms these four combinations equal the four quadrant formulas from block multiplication above — the algebra is denser than Karatsuba's but the principle is the same substitution trick. With only 7 recursive multiplications:

$$T(n) = 7T(n/2) + O(n^2) = \Theta(n^{\log_2 7}) \approx O(n^{2.807})$$

This is exactly Example 4 of the Master Theorem walkthrough in the Recursion Complete Guide, which derives the $n^{2.807}$ exponent from $a=7$, $b=2$, $f(n)=n^2$. As with Karatsuba, the exponent drop (from 3 to 2.807) comes entirely from $a$ dropping from 8 to 7.

The implementation below assumes square, power-of-2-sized matrices for simplicity. Non-power-of-2 sizes require padding the matrix with zero rows and columns up to the next power of 2 before recursing, then discarding the padding from the result — a mechanical preprocessing step, described here in prose rather than implemented, since it adds no new algorithmic idea.

def naive_multiply(A, B):
    n = len(A)
    C = [[0] * n for _ in range(n)]
    for i in range(n):
        for k in range(n):
            if A[i][k] == 0:
                continue
            for j in range(n):
                C[i][j] += A[i][k] * B[k][j]
    return C

def add(A, B):
    n = len(A)
    return [[A[i][j] + B[i][j] for j in range(n)] for i in range(n)]

def sub(A, B):
    n = len(A)
    return [[A[i][j] - B[i][j] for j in range(n)] for i in range(n)]

def split(M):
    n = len(M)
    h = n // 2
    return ([row[:h] for row in M[:h]], [row[h:] for row in M[:h]],
            [row[:h] for row in M[h:]], [row[h:] for row in M[h:]])

def combine(C11, C12, C21, C22):
    h = len(C11)
    C = [[0] * (2 * h) for _ in range(2 * h)]
    for i in range(h):
        for j in range(h):
            C[i][j], C[i][j + h] = C11[i][j], C12[i][j]
            C[i + h][j], C[i + h][j + h] = C21[i][j], C22[i][j]
    return C

def strassen(A, B, base=64):
    n = len(A)
    if n <= base:                       # base case: naive multiply is faster for small n
        return naive_multiply(A, B)

    A11, A12, A21, A22 = split(A)
    B11, B12, B21, B22 = split(B)

    M1 = strassen(add(A11, A22), add(B11, B22), base)
    M2 = strassen(add(A21, A22), B11, base)
    M3 = strassen(A11, sub(B12, B22), base)
    M4 = strassen(A22, sub(B21, B11), base)
    M5 = strassen(add(A11, A12), B22, base)
    M6 = strassen(sub(A21, A11), add(B11, B12), base)
    M7 = strassen(sub(A12, A22), add(B21, B22), base)

    C11 = add(sub(add(M1, M4), M5), M7)
    C12 = add(M3, M5)
    C21 = add(M2, M4)
    C22 = add(sub(add(M1, M3), M2), M6)
    return combine(C11, C12, C21, C22)

A = [[1, 2], [3, 4]]
B = [[5, 6], [7, 8]]
print(strassen(A, B, base=1))      # [[19, 22], [43, 50]]
print(naive_multiply(A, B))        # [[19, 22], [43, 50]]

#include <vector>
using Matrix = std::vector<std::vector<long long>>;

Matrix naiveMultiply(const Matrix& A, const Matrix& B) {
    int n = (int)A.size();
    Matrix C(n, std::vector<long long>(n, 0));
    for (int i = 0; i < n; ++i)
        for (int k = 0; k < n; ++k) {
            long long aik = A[i][k];
            if (aik == 0) continue;
            for (int j = 0; j < n; ++j) C[i][j] += aik * B[k][j];
        }
    return C;
}

Matrix add(const Matrix& A, const Matrix& B) {
    int n = (int)A.size();
    Matrix C(n, std::vector<long long>(n));
    for (int i = 0; i < n; ++i) for (int j = 0; j < n; ++j) C[i][j] = A[i][j] + B[i][j];
    return C;
}

Matrix sub(const Matrix& A, const Matrix& B) {
    int n = (int)A.size();
    Matrix C(n, std::vector<long long>(n));
    for (int i = 0; i < n; ++i) for (int j = 0; j < n; ++j) C[i][j] = A[i][j] - B[i][j];
    return C;
}

void split(const Matrix& M, Matrix& A11, Matrix& A12, Matrix& A21, Matrix& A22) {
    int h = (int)M.size() / 2;
    A11 = A12 = A21 = A22 = Matrix(h, std::vector<long long>(h));
    for (int i = 0; i < h; ++i)
        for (int j = 0; j < h; ++j) {
            A11[i][j] = M[i][j];       A12[i][j] = M[i][j + h];
            A21[i][j] = M[i + h][j];   A22[i][j] = M[i + h][j + h];
        }
}

Matrix combine(const Matrix& C11, const Matrix& C12, const Matrix& C21, const Matrix& C22) {
    int h = (int)C11.size();
    Matrix C(2 * h, std::vector<long long>(2 * h));
    for (int i = 0; i < h; ++i)
        for (int j = 0; j < h; ++j) {
            C[i][j] = C11[i][j];             C[i][j + h] = C12[i][j];
            C[i + h][j] = C21[i][j];         C[i + h][j + h] = C22[i][j];
        }
    return C;
}

// base case: naive multiply is faster for small n
Matrix strassen(const Matrix& A, const Matrix& B, int base = 64) {
    int n = (int)A.size();
    if (n <= base) return naiveMultiply(A, B);

    Matrix A11, A12, A21, A22, B11, B12, B21, B22;
    split(A, A11, A12, A21, A22);
    split(B, B11, B12, B21, B22);

    Matrix M1 = strassen(add(A11, A22), add(B11, B22), base);
    Matrix M2 = strassen(add(A21, A22), B11, base);
    Matrix M3 = strassen(A11, sub(B12, B22), base);
    Matrix M4 = strassen(A22, sub(B21, B11), base);
    Matrix M5 = strassen(add(A11, A12), B22, base);
    Matrix M6 = strassen(sub(A21, A11), add(B11, B12), base);
    Matrix M7 = strassen(sub(A12, A22), add(B21, B22), base);

    Matrix C11 = add(sub(add(M1, M4), M5), M7);
    Matrix C12 = add(M3, M5);
    Matrix C21 = add(M2, M4);
    Matrix C22 = add(sub(add(M1, M3), M2), M6);
    return combine(C11, C12, C21, C22);
}
// strassen({{1,2},{3,4}}, {{5,6},{7,8}}, 1) == {{19,22},{43,50}}

import java.util.Arrays;

public class Strassen {
    static long[][] naiveMultiply(long[][] a, long[][] b) {
        int n = a.length;
        long[][] c = new long[n][n];
        for (int i = 0; i < n; i++)
            for (int k = 0; k < n; k++) {
                if (a[i][k] == 0) continue;
                for (int j = 0; j < n; j++) c[i][j] += a[i][k] * b[k][j];
            }
        return c;
    }

    static long[][] add(long[][] a, long[][] b) {
        int n = a.length;
        long[][] c = new long[n][n];
        for (int i = 0; i < n; i++) for (int j = 0; j < n; j++) c[i][j] = a[i][j] + b[i][j];
        return c;
    }

    static long[][] sub(long[][] a, long[][] b) {
        int n = a.length;
        long[][] c = new long[n][n];
        for (int i = 0; i < n; i++) for (int j = 0; j < n; j++) c[i][j] = a[i][j] - b[i][j];
        return c;
    }

    static long[][] quadrant(long[][] m, int rowOff, int colOff, int h) {
        long[][] q = new long[h][h];
        for (int i = 0; i < h; i++)
            for (int j = 0; j < h; j++) q[i][j] = m[i + rowOff][j + colOff];
        return q;
    }

    static long[][] combine(long[][] c11, long[][] c12, long[][] c21, long[][] c22) {
        int h = c11.length;
        long[][] c = new long[2 * h][2 * h];
        for (int i = 0; i < h; i++)
            for (int j = 0; j < h; j++) {
                c[i][j] = c11[i][j];               c[i][j + h] = c12[i][j];
                c[i + h][j] = c21[i][j];           c[i + h][j + h] = c22[i][j];
            }
        return c;
    }

    // base case: naive multiply is faster for small n
    static long[][] strassen(long[][] a, long[][] b, int base) {
        int n = a.length;
        if (n <= base) return naiveMultiply(a, b);
        int h = n / 2;

        long[][] a11 = quadrant(a, 0, 0, h), a12 = quadrant(a, 0, h, h);
        long[][] a21 = quadrant(a, h, 0, h), a22 = quadrant(a, h, h, h);
        long[][] b11 = quadrant(b, 0, 0, h), b12 = quadrant(b, 0, h, h);
        long[][] b21 = quadrant(b, h, 0, h), b22 = quadrant(b, h, h, h);

        long[][] m1 = strassen(add(a11, a22), add(b11, b22), base);
        long[][] m2 = strassen(add(a21, a22), b11, base);
        long[][] m3 = strassen(a11, sub(b12, b22), base);
        long[][] m4 = strassen(a22, sub(b21, b11), base);
        long[][] m5 = strassen(add(a11, a12), b22, base);
        long[][] m6 = strassen(sub(a21, a11), add(b11, b12), base);
        long[][] m7 = strassen(sub(a12, a22), add(b21, b22), base);

        long[][] c11 = add(sub(add(m1, m4), m5), m7);
        long[][] c12 = add(m3, m5);
        long[][] c21 = add(m2, m4);
        long[][] c22 = add(sub(add(m1, m3), m2), m6);
        return combine(c11, c12, c21, c22);
    }

    public static void main(String[] args) {
        long[][] a = {{1, 2}, {3, 4}};
        long[][] b = {{5, 6}, {7, 8}};
        System.out.println(Arrays.deepToString(strassen(a, b, 1)));   // [[19, 22], [43, 50]]
        System.out.println(Arrays.deepToString(naiveMultiply(a, b))); // [[19, 22], [43, 50]]
    }
}

All three implementations were executed and cross-checked against naive matrix multiplication on random $2\times2$, $4\times4$, $8\times8$, and $16\times16$ integer matrices, plus the hand-traceable $2\times2$ case above; every result matched exactly.

Practical Crossover Points

Karatsuba and Strassen both trade multiplications for additions, and additions are only asymptotically free — at small sizes, the extra bookkeeping (recursive calls, allocating intermediate sums, combining results) costs more than the multiplications it avoids. Both algorithms need a crossover size below which the naive method is faster in practice, which is exactly why both implementations above include a base parameter that falls back to direct multiplication for small inputs.

• Karatsuba: commonly cited crossover points are in the range of a few hundred digits, depending on the platform and how well-optimized the base-case multiplication is. Below that, grade-school multiplication wins.
• Strassen: commonly cited crossover points are matrix dimensions in the low hundreds. Beyond the constant-factor overhead, Strassen's subtractions can amplify floating-point rounding error in ways naive triple-loop multiplication does not, so numerical-stability concerns make it less used in general-purpose numerical libraries than its asymptotic complexity alone would suggest.
• Recursion depth: both algorithms halve the input size at each level, so recursion depth is $O(\log n)$ for either one — a few dozen levels even for enormous inputs. Stack depth is essentially never the practical bottleneck; constant-factor overhead per call is.

Production libraries (arbitrary-precision arithmetic, BLAS-style linear algebra) typically implement exactly this pattern: naive method below a tuned threshold, Karatsuba/Strassen above it, and sometimes a third tier (Toom-Cook, or FFT-based multiplication) above an even larger threshold.

7. Complexity Comparison

8. Applications & Alternatives

Even faster algorithms exist asymptotically but are beyond this lesson's scope: the Schönhage–Strassen algorithm multiplies large integers via the Fast Fourier Transform in $O(n \log n \log\log n)$, and the Coppersmith–Winograd family of algorithms pushes matrix multiplication's exponent below Strassen's $2.807$ toward the theoretical lower bound of 2 — at the cost of galactic constant factors that make them impractical for any realistic input size.

9. Exercises

Easy

1. Hand-trace Karatsuba on $34 \times 56$: compute $z_2$, $z_0$, and $z_1$, then combine them and verify against direct multiplication.

Medium

1. Hand-trace Strassen's 7-product combination on $A=\begin{pmatrix}2&0\\1&3\end{pmatrix}$ and $B=\begin{pmatrix}1&4\\2&1\end{pmatrix}$ (treat each entry as a $1\times1$ block). Compute $M_1$ through $M_7$, combine into $C_{11}, C_{12}, C_{21}, C_{22}$, and verify against naive $2\times2$ matrix multiplication.

Hard

1. The Karatsuba implementations above assume both numbers have roughly the same digit count. Design (do not necessarily implement) an extension that handles numbers of very different lengths, or an odd total digit count at the split point. Consider: does padding the shorter number with leading zeros to match a split point change the recurrence? What happens to the split when $n$ is odd — must $n/2$ round the same way for both halves?

Review rhythm: tomorrow, reconstruct the 3-product and 7-product formulas from scratch without notes; next week, implement both algorithms from memory in one language and test them against naive multiplication; several weeks later, recognize the same "shrink the branching factor" trick if it appears in an unfamiliar divide-and-conquer problem. Finally, teach someone who already knows the Master Theorem why $a$ dropping by one is enough to change the asymptotic complexity.

10. Summary Cheat Sheet

• Naive multiplication and naive matrix multiplication are both $\Theta(n^2)$ and $\Theta(n^3)$ respectively; splitting the input recursively without reducing the number of sub-products changes nothing.
• Karatsuba: 3 products instead of 4, $T(n)=3T(n/2)+O(n)=O(n^{1.585})$.
• Toom-$k$: $2k-1$ products instead of $k^2$, approaching $O(n^{1+\epsilon})$ as $k$ grows.
• Strassen: 7 products instead of 8, $T(n)=7T(n/2)+O(n^2)=O(n^{2.807})$.
• Both algorithms recover missing cross terms algebraically through addition/subtraction of a smaller set of products — the same trick, applied to numbers and to matrices.
• Both have significant constant-factor overhead; real crossover points are a few hundred digits (Karatsuba) or low hundreds of matrix dimension (Strassen), not the small examples used to teach the algorithm.
• Recursion depth is $O(\log n)$ for both — never the practical bottleneck.
• The Master Theorem examples for both recurrences already live in the Recursion Complete Guide; this lesson supplies the algorithms those examples only cited.