18.13 Modular Combinatorics

Many counting problems ask for results modulo a large number.

18.13 Modular Combinatorics

Many counting problems ask for results modulo a large number. The answer may be enormous, but the required output is usually:

answer mod M

This is common in combinatorics, dynamic programming, graph counting, probability calculations, and programming-contest problems. Modular combinatorics provides the tools for computing factorials, permutations, combinations, and related counting formulas without ever constructing the full integer answer.

The most common setting uses a large prime modulus such as:

1_000_000_007

or:

998_244_353

Prime moduli are convenient because every nonzero residue has a modular inverse.

Problem

Given integers:

n
k
M

compute expressions such as:

n!
P(n,k)
C(n,k)

modulo M.

The central example is the binomial coefficient:

C(n,k) = n! / (k!(n-k)!)

Under modular arithmetic, division must be replaced by multiplication with modular inverses.

Factorials Modulo M

A factorial is the product:

n! = 1 × 2 × 3 × ... × n

The modular version reduces after every multiplication.

def factorial_mod(n, mod):
    result = 1

    for x in range(2, n + 1):
        result = (result * x) % mod

    return result

Example:

print(factorial_mod(5, 13))

Output:

3

because:

5! = 120
120 mod 13 = 3

Precomputing Factorials

When many queries use the same modulus, precompute factorials.

def factorial_table(n, mod):
    fact = [1] * (n + 1)

    for i in range(1, n + 1):
        fact[i] = fact[i - 1] * i % mod

    return fact

Example:

fact = factorial_table(10, 1_000_000_007)

print(fact[5])
print(fact[10])

This provides O(1) factorial lookup after O(n) preprocessing.

Modular Inverse Factorials

For combinations, we need inverse factorials:

(0!)⁻¹, (1!)⁻¹, ..., (n!)⁻¹

When mod is prime, use Fermat's Little Theorem:

a⁻¹ ≡ a^(mod-2) (mod mod)

Efficient construction:

def factorial_and_inverse_tables(n, mod):
    fact = [1] * (n + 1)
    inv_fact = [1] * (n + 1)

    for i in range(1, n + 1):
        fact[i] = fact[i - 1] * i % mod

    inv_fact[n] = pow(fact[n], mod - 2, mod)

    for i in range(n, 0, -1):
        inv_fact[i - 1] = inv_fact[i] * i % mod

    return fact, inv_fact

This computes all inverse factorials in linear time after one modular exponentiation.

Binomial Coefficients

The binomial coefficient is:

C(n,k) = n! / (k!(n-k)!)

Under a prime modulus:

C(n,k) mod M
=
fact[n] × inv_fact[k] × inv_fact[n-k] mod M

Implementation:

def comb(n, k, fact, inv_fact, mod):
    if k < 0 or k > n:
        return 0

    return fact[n] * inv_fact[k] % mod * inv_fact[n - k] % mod

Example:

MOD = 1_000_000_007
fact, inv_fact = factorial_and_inverse_tables(100, MOD)

print(comb(5, 2, fact, inv_fact, MOD))

Output:

10

Why Precomputation Matters

A direct computation of C(n,k) for each query would repeatedly recompute factorials and inverses.

With precomputation:

Preprocessing: O(n)
Query:         O(1)

This is essential when handling many queries.

Example use cases:

C(100000, 50000)
C(100000, 3)
C(99999, 12345)

All become constant-time queries after tables are built.

Permutations

The number of ordered selections of k items from n is:

P(n,k) = n! / (n-k)!

Modular implementation:

def perm(n, k, fact, inv_fact, mod):
    if k < 0 or k > n:
        return 0

    return fact[n] * inv_fact[n - k] % mod

Example:

print(perm(5, 2, fact, inv_fact, MOD))

Output:

20

because there are 20 ordered ways to choose 2 items from 5.

Combinations Without Precomputation

For a single query with moderate k, compute directly:

def comb_single(n, k, mod):
    if k < 0 or k > n:
        return 0

    k = min(k, n - k)

    numerator = 1
    denominator = 1

    for i in range(1, k + 1):
        numerator = numerator * (n - k + i) % mod
        denominator = denominator * i % mod

    return numerator * pow(denominator, mod - 2, mod) % mod

Complexity:

O(k + log mod)

This is useful when n is large but only one or a few queries are needed.

Pascal Recurrence

Binomial coefficients also satisfy:

C(n,k) = C(n-1,k-1) + C(n-1,k)

This recurrence works directly modulo M.

def pascal_table(n, mod):
    C = [[0] * (n + 1) for _ in range(n + 1)]

    for i in range(n + 1):
        C[i][0] = 1
        C[i][i] = 1

        for k in range(1, i):
            C[i][k] = (C[i - 1][k - 1] + C[i - 1][k]) % mod

    return C

Complexity:

O(n²)

This is useful when n is small and many values of C(n,k) are needed, or when the modulus is not prime and modular inverses are inconvenient.

Lucas Theorem

When the modulus p is prime and n and k are very large, Lucas Theorem can compute:

C(n,k) mod p

using base-p digits.

Write:

n = n₀ + n₁p + n₂p² + ...
k = k₀ + k₁p + k₂p² + ...

Then:

C(n,k) ≡ Π C(nᵢ,kᵢ) (mod p)

Implementation:

def lucas(n, k, p):
    fact, inv_fact = factorial_and_inverse_tables(p - 1, p)

    result = 1

    while n > 0 or k > 0:
        ni = n % p
        ki = k % p

        if ki > ni:
            return 0

        result = result * comb(ni, ki, fact, inv_fact, p) % p

        n //= p
        k //= p

    return result

Lucas Theorem is especially useful when p is small and n is huge.

Stars and Bars

Many counting problems reduce to combinations.

The number of nonnegative integer solutions to:

x₁ + x₂ + ... + xₖ = n

is:

C(n+k-1, k-1)

Example:

x₁ + x₂ + x₃ = 5

Number of solutions:

C(5+3-1, 3-1) = C(7,2) = 21

Modular computation:

def stars_and_bars(n, k, fact, inv_fact, mod):
    return comb(n + k - 1, k - 1, fact, inv_fact, mod)

This pattern appears in distribution, partitioning, and dynamic programming problems.

Catalan Numbers

Catalan numbers count many recursive structures:

balanced parentheses
binary trees
triangulations
monotonic lattice paths

The nth Catalan number is:

Cat(n) = C(2n,n) / (n+1)

Modulo a prime:

def catalan(n, fact, inv_fact, mod):
    value = comb(2 * n, n, fact, inv_fact, mod)
    return value * pow(n + 1, mod - 2, mod) % mod

Example:

print(catalan(3, fact, inv_fact, MOD))

Output:

5

The first few Catalan numbers are:

1 1 2 5 14 42

Multinomial Coefficients

The multinomial coefficient counts ways to split n labeled objects into groups of sizes:

k₁, k₂, ..., kᵣ

where:

k₁ + k₂ + ... + kᵣ = n

Formula:

n! / (k₁!k₂!...kᵣ!)

Implementation:

def multinomial(parts, fact, inv_fact, mod):
    n = sum(parts)
    result = fact[n]

    for k in parts:
        if k < 0:
            return 0

        result = result * inv_fact[k] % mod

    return result

Example:

print(multinomial([2, 1, 2], fact, inv_fact, MOD))

This computes:

5! / (2!1!2!)

which equals:

30

Composite Moduli

The factorial-inverse method assumes that required denominators are invertible modulo M.

When M is prime and n < M, this usually works.

For composite moduli, denominators may not have inverses.

Example:

2⁻¹ mod 6

does not exist because:

gcd(2,6) = 2

Options include:

Pascal recurrence
prime-power decomposition
Chinese Remainder Theorem
factor tracking

For many practical problems, the modulus is deliberately chosen to be prime to avoid this complication.

Common Mistakes

Dividing Directly

This is wrong:

answer = fact[n] // (fact[k] * fact[n-k])
answer %= mod

Integer division before taking modulo does not model modular division.

Use inverse factorials.

Using Fermat with Composite Modulus

The expression:

a^(M-2) mod M

computes an inverse only when M is prime.

Forgetting Invalid k

For combinations:

C(n,k) = 0

when:

k < 0 or k > n

Handle these cases explicitly.

Precomputing Too Small

If you need C(2n,n) for Catalan numbers, factorials must be built up to 2n, not n.

Overflow Before Modulo

In fixed-width languages, multiplication may overflow before % is applied. Use a wider integer type when necessary.

Testing

Basic identities:

def test_combinations():
    MOD = 1_000_000_007
    fact, inv_fact = factorial_and_inverse_tables(100, MOD)

    assert comb(5, 0, fact, inv_fact, MOD) == 1
    assert comb(5, 5, fact, inv_fact, MOD) == 1
    assert comb(5, 2, fact, inv_fact, MOD) == 10
    assert comb(10, 3, fact, inv_fact, MOD) == 120

Symmetry:

def test_symmetry():
    MOD = 1_000_000_007
    fact, inv_fact = factorial_and_inverse_tables(100, MOD)

    for n in range(0, 100):
        for k in range(0, n + 1):
            assert comb(n, k, fact, inv_fact, MOD) == comb(n, n - k, fact, inv_fact, MOD)

Pascal identity:

def test_pascal_identity():
    MOD = 1_000_000_007
    fact, inv_fact = factorial_and_inverse_tables(100, MOD)

    for n in range(2, 100):
        for k in range(1, n):
            left = comb(n, k, fact, inv_fact, MOD)
            right = (
                comb(n - 1, k - 1, fact, inv_fact, MOD)
                + comb(n - 1, k, fact, inv_fact, MOD)
            ) % MOD

            assert left == right

These tests catch most indexing and inverse-table errors.

Takeaway

Modular combinatorics turns large counting expressions into manageable modular computations. Factorial and inverse-factorial tables allow binomial and permutation queries in constant time after linear preprocessing, provided the modulus is prime. Pascal recurrence, Lucas Theorem, Catalan numbers, stars and bars, and multinomial coefficients build on the same foundation. The central implementation rule is simple: modular division must be performed by multiplying with a valid modular inverse.