Lab 07 — Combinatorics mod p

Goal

Extend the fact[]/inv_fact[] engine from lab-01 (modular arithmetic) into the full counting toolkit: O(1) binomial coefficients, permutations, multinomials, Catalan numbers, stars and bars, and Lucas’ theorem for n up to 10¹⁸. Every CF round with the line “output the answer modulo 10⁹+7” assumes this is preloaded. Target: write the Comb class from blank in <4 minutes; map a counting problem to the right formula in <2 minutes.

Background Concepts

The engine. Precompute fact[i] = i! mod p forward in O(n). Compute inv_fact[n] = fact[n]^(p−2) mod p with one binary exponentiation (Fermat’s little theorem, prime p — see lab-02), then sweep backward: inv_fact[i−1] = inv_fact[i] · i mod p, because (i−1)!⁻¹ = i!⁻¹ · i. Total O(n + log p). After that:

  • nCr(n, r) = fact[n] · inv_fact[r] · inv_fact[n−r] — O(1).
  • nPr(n, r) = fact[n] · inv_fact[n−r] — ordered selections.
  • Multinomial n! / (k₁! k₂! … k_m!) = fact[n] · Π inv_fact[kᵢ] — arrangements with repeated items.

Stars and bars. Number of ways to write n as an ordered sum of k non-negative integers = C(n+k−1, k−1) (place k−1 bars among n stars). Positive integers: substitute xᵢ = yᵢ+1 → C(n−1, k−1).

Catalan numbers. Cat(n) = C(2n, n)/(n+1) = (2n)! / (n!(n+1)!) = C(2n, n) − C(2n, n+1). They count: valid parentheses sequences of n pairs (the count behind LC 22), distinct BSTs on n nodes (LC 96), Dyck paths (lattice paths that never dip below the diagonal), non-crossing chord pairings, and triangulations of an (n+2)-gon.

Lucas’ theorem. For prime p: C(n, r) ≡ Π C(nᵢ, rᵢ) mod p, where nᵢ, rᵢ are the base-p digits. Any digit with rᵢ > nᵢ zeroes the product. Lets you answer C(10¹⁸, r) mod p for small prime p with a table of size p.

Derangements. D(n) = n! · Σ (−1)ᵏ/k! = (n−1)(D(n−1) + D(n−2)); D(n) is the nearest integer to n!/e. Classic inclusion-exclusion output (phase-12 lab-10).

Composite modulus. Fermat fails. Factor m (lab-03 sieve), compute mod each prime power (Lucas generalizes via Andrew Granville / factorial-with-p-stripped tricks), recombine with CRT — machinery in phase-13 lab-07 (ext-gcd/CRT). Name it in the interview; implement it only in CP.

Interview Context

FAANG rarely asks nCr mod p cold, but LC 62 rewards the formula answer over DP, LC 96 is Catalan recognition, and LC 920 (Hard) is unsolvable in time without modular counting fluency. Quant/HFT interviews use these as direct filters — “how many monotone lattice paths avoid the diagonal?” expects Catalan by reflection, not DP. The strongest signal: reaching for inv_fact without narrating it, the way you’d use a hash map.

Problem Statement (LC 62, LC 96, LC 920)

Canonical framing: count paths / sequences / arrangements, output mod 10⁹+7. Three graded instances:

  • LC 62 (Unique Paths): paths from (0,0) to (m−1,n−1) moving only right/down = arrangements of m−1 downs and n−1 rights = C(m+n−2, m−1).
  • LC 96 (Unique BSTs): answer is Cat(n).
  • LC 920 (Number of Music Playlists): n distinct songs, playlist length L, every song played, a song may repeat only after k others — DP over distinct-songs-used whose transitions are modular products; closed form uses C(n, j) with inclusion-exclusion.

Constraints

  • Table size N ≤ 10⁶–10⁷ (build once, O(N)).
  • Queries Q ≤ 10⁶, O(1) each.
  • p = 10⁹+7 (prime) unless stated; Lucas variant: n, r ≤ 10¹⁸, p ≤ 10⁵ prime.
  • LC 62: 1 ≤ m, n ≤ 100. LC 920: 0 ≤ k < n ≤ L ≤ 100.

Clarifying Questions

  1. Is the modulus prime? (Yes, 10⁹+7 — Fermat inverses valid; if composite, ext-gcd inverses + CRT, phase-13 lab-07.)
  2. How large can n get relative to my table? (If n exceeds table size but p is small, Lucas; otherwise size the table to max n — and 2n for Catalan.)
  3. Are r < 0 or r > n possible? (Yes at formula boundaries — return 0, don’t crash.)
  4. Ordered or unordered selections? (Determines nPr vs nCr — say it out loud before coding.)
  5. Do “objects” repeat? (Repetition → multinomial or stars and bars, not plain nCr.)

Examples

C(5, 2) = 10          P(5, 2) = 20         multinomial(4; 2,1,1) = 12
stars&bars: x+y+z = 5, xᵢ ≥ 0 → C(7, 2) = 21
Cat(3) = 5   (LC 96: numTrees(3) = 5; LC 22: 5 valid strings of 3 pairs)
LC 62: m=3, n=7 → C(8, 2) = 28
Lucas: C(10¹⁸, 3) mod 7 = 0  (10¹⁸ ≡ 1 mod 7, so the n−1 factor ≡ 0);  C(10¹⁸+5, 3) mod 7 = 6

Initial Brute Force

Pascal’s triangle DP: C[n][r] = C[n−1][r−1] + C[n−1][r] mod p, O(n²) table. Or Python big-int math.factorial(n) // (…) then mod. Both are fine oracles for n ≤ 2000.

Brute Force Complexity

Pascal: O(n²) time and space — 10¹² cells at n = 10⁶. Big-int factorial: O(n) multiplications of numbers with O(n log n) digits ≈ O(n² log n) bit ops per query. Both infeasible at CP scale; keep them as stress oracles.

Optimization Path

  1. Pascal O(n²) → correctness oracle.
  2. fact[] table + Fermat inverse per query: O(n + Q log p). Passes, wasteful.
  3. fact[] + backward inv_fact[] sweep: O(n + log p) build, O(1) query. Optimal — one pow total.
  4. n beyond any table, small prime p → Lucas, O(p + log_p n) per query.

Final Expected Approach

MOD = 10**9 + 7

class Comb:
    """O(1) nCr/nPr/Catalan after O(n) precomputation, prime modulus."""
    def __init__(self, n, mod=MOD):
        self.mod = mod
        self.fact = fact = [1] * (n + 1)
        for i in range(1, n + 1):
            fact[i] = fact[i - 1] * i % mod
        self.inv_fact = inv = [1] * (n + 1)
        inv[n] = pow(fact[n], mod - 2, mod)          # the single Fermat inverse
        for i in range(n, 0, -1):
            inv[i - 1] = inv[i] * i % mod            # backward sweep

    def ncr(self, n, r):
        if r < 0 or r > n:
            return 0
        return self.fact[n] * self.inv_fact[r] % self.mod * self.inv_fact[n - r] % self.mod

    def npr(self, n, r):
        if r < 0 or r > n:
            return 0
        return self.fact[n] * self.inv_fact[n - r] % self.mod

    def catalan(self, n):
        # Cat(n) = (2n)! / (n! (n+1)!) — table must reach 2n
        return self.fact[2 * n] * self.inv_fact[n] % self.mod * self.inv_fact[n + 1] % self.mod

def lucas(n, r, p):
    """C(n, r) mod prime p, p small enough to table; n, r up to 10^18."""
    fact = [1] * p
    for i in range(1, p):
        fact[i] = fact[i - 1] * i % p
    res = 1
    while n or r:
        ni, ri = n % p, r % p
        if ri > ni:
            return 0
        res = res * fact[ni] % p
        res = res * pow(fact[ri] * fact[ni - ri] % p, p - 2, p) % p
        n //= p
        r //= p
    return res

def unique_paths_dp(m, n):                            # LC 62, O(mn) oracle
    row = [1] * n
    for _ in range(m - 1):
        for j in range(1, n):
            row[j] = (row[j] + row[j - 1]) % MOD
    return row[-1]

def unique_paths_formula(m, n, C):                    # LC 62, O(1) after tables
    return C.ncr(m + n - 2, m - 1)

if __name__ == "__main__":
    C = Comb(4000)
    for m in range(1, 41):
        for n in range(1, 41):
            assert unique_paths_dp(m, n) == unique_paths_formula(m, n, C)
    assert unique_paths_formula(3, 7, C) == 28
    assert [C.catalan(i) for i in range(6)] == [1, 1, 2, 5, 14, 42]
    assert lucas(10**18, 3, 7) == 0 and lucas(10**18 + 5, 3, 7) == 6

Data Structures Used

  • Two flat int arrays fact[0..N], inv_fact[0..N].
  • Lucas: one array of size p plus the base-p digit loop — no big tables.

Correctness Argument

Backward sweep: i! · (i!)⁻¹ ≡ 1 at i = n by Fermat (p prime, fact[n] ≢ 0 since n < p). Inductively, inv_fact[i−1] = inv_fact[i] · i = (i!)⁻¹ · i = ((i−1)!)⁻¹ because (i−1)! · i = i!. So every entry is a true inverse with a single pow.

nCr: n! / (r!(n−r)!) is an integer; modulo a prime, division by r!(n−r)! is multiplication by its inverse, valid because r, n−r < p keeps both factorials nonzero mod p.

Lucas (sketch): in F_p[x], (1+x)^p ≡ 1 + x^p (freshman’s dream — binomials C(p,k) vanish for 0<k<p). Writing n = Σ nᵢ pⁱ gives (1+x)^n ≡ Π (1+x^{pⁱ})^{nᵢ}; extracting the coefficient of x^r digit by digit yields Π C(nᵢ, rᵢ).

LC 62: every path is a length-(m+n−2) string over {D, R} with exactly m−1 D’s; strings ↔ paths is a bijection, hence C(m+n−2, m−1). The DP-vs-formula assertion over 1600 cases is the machine check.

Complexity

OperationTimeSpace
Build Comb(N)O(N + log p)O(N)
ncr / npr / catalanO(1)—
lucas(n, r, p)O(p + log_p n · log p)O(p)

Implementation Requirements

  • Exactly one pow call in the constructor — a pow per entry is the classic O(n log p) regression.
  • Guard r < 0 or r > n → 0 in every query path.
  • Catalan requires the table sized to 2n, and inv_fact[n+1] — size for it up front.
  • Lucas must reject composite p (document it; assert if cheap — Miller–Rabin in phase-13 lab-06).
  • Keep everything in one Comb class you can retype from memory.

Tests

  • ncr/npr vs math.comb/factorial oracle for all n ≤ 200, spot checks at n ≈ 4000.
  • Boundaries: C(0,0)=1, C(n,0)=C(n,n)=1, C(5,7)=0, C(5,−1)=0.
  • Catalan vs 1, 1, 2, 5, 14, 42, 132, … and vs comb(2n,n)//(n+1) for n ≤ 200.
  • Lucas vs math.comb % p exhaustively for n ≤ 120, p ∈ {2,3,5,7,13}; large-n spot check C(10¹⁸,3) mod 7 = 0 via the polynomial n(n−1)(n−2)/6 mod 7.
  • Parity law: C(n,r) odd ⇔ r AND (n−r) = 0 — checked against lucas(·,·,2) for n < 64.
  • LC 62 DP vs formula for all m, n ≤ 40. (All of the above executed for this lab.)

Follow-up Questions

  • Why build inv_fact backwards instead of n modular inverses? → n Fermat calls cost O(n log p) — 30× slower at n = 10⁷. One pow at the top plus the O(n) identity (i−1)!⁻¹ = i!⁻¹ · i gets every inverse. (Bonus: inv[i] = inv_fact[i] · fact[i−1] then gives plain inverses of 1..n in O(n) too.)
  • When is C(n, r) odd? → Lucas base 2: odd iff every bit of r is ≤ the corresponding bit of n, i.e. r AND (n−r) = 0. Kummer’s theorem generalizes: the exponent of p in C(n, r) equals the number of carries when adding r and n−r in base p.
  • Derive stars and bars. → arrange n stars and k−1 bars in a row; each arrangement ↔ one ordered solution of x₁+…+x_k = n, xᵢ ≥ 0; choose bar positions: C(n+k−1, k−1). LC disguises: distribute candies, count non-decreasing digit strings, bounded knapsack counting after inclusion-exclusion on upper bounds.
  • How do you recognize a Catalan problem? → checklist: balanced open/close structure, a “never dip below zero” prefix condition, non-crossing pairings, or recursive splitting into left/right independent subproblems (Cat(n) = Σ Cat(i)Cat(n−1−i)). Any one → try C(2n,n)/(n+1) before writing DP.
  • Count surjections from n elements onto k labels. → inclusion-exclusion over excluded labels: Σ (−1)ʲ C(k,j)(k−j)ⁿ. Same alternating pattern gives derangements. Phase-12 lab-10 drills the technique; this lab supplies the O(1) C(k,j).
  • Modulus 10⁹+6 — now what? → composite (2 × 500000003). Factor, solve mod each prime power, recombine with CRT (phase-13 lab-07); non-invertible factorials require carrying the exponent of each prime separately.

Product Extension

A poker/board-game backend computing exact win equities evaluates thousands of C(n, r) per hand — precompute tables once at service start, never at request time. Same shape in ranking systems: “how many ways can k of n experiments beat control” is a binomial tail; a fintech risk engine summing hypergeometric probabilities is ncr three times per term. The O(n)-build/O(1)-query split is exactly a cache-warming pattern.

Language/Runtime Follow-ups

  • Python: pow(a, -1, mod) (3.8+) is the idiomatic single inverse; math.comb is exact big-int — fine as oracle, too slow as engine at 10⁶ queries.
  • Java: keep everything long; a · b with both < 10⁹+7 fits in 63 bits — but three-factor products need an intermediate % mod.
  • C++: long long with % MOD after every multiply, or a Montgomery/mint wrapper; __int128 only if you must multiply two 62-bit values.
  • Go: uint64 works with the same discipline; no operator overloading → write mulmod helpers.
  • JS/TS: Number breaks past 2⁵³ — a single product of two ~10⁹ values already overflows. Use BigInt throughout or split multiplies.

Common Bugs

  1. Forward inverse loop: computing inv_fact[i] from inv_fact[i−1] — the identity runs backward; forward needs a division you don’t have.
  2. Table one short: catalan(n) touches fact[2n] and inv_fact[n+1]; a table sized to n throws IndexError (best case) or reads garbage (C++).
  3. Missing range guard: ncr(5, 7) must return 0; without the guard, inv_fact[−2] silently indexes from the end in Python — wrong answer, no crash.
  4. Negative after subtraction: C(2n,n) − C(2n,n+1) can be negative pre-mod in C++/Java — add mod before %. (Python % is already non-negative.)
  5. Per-query pow: O(Q log p) instead of O(Q); passes small tests, TLEs at Q = 10⁶ — the regression is invisible until scale.
  6. Lucas with composite p: both Fermat inverses and the theorem itself require prime p; composite moduli silently produce wrong digits products.

Debugging Strategy

First check the two identities fact[i] · inv_fact[i] ≡ 1 and ncr(n, r) == math.comb(n, r) % p for all n ≤ 50 — this isolates table bugs from formula bugs. For Lucas, print the digit pairs (nᵢ, rᵢ) per iteration: C(10, 4) mod 3 → digits (1,1),(0,1),(1,0) → the (0,1) pair zeroes it; math.comb(10,4)=210=3·70 confirms. For a wrong LC 62/920 answer, shrink to m, n ≤ 5 and diff formula vs Pascal DP row by row — the first divergent cell names the wrong term.

Mastery Criteria

  • Comb class (ncr, npr, catalan) from blank, correct first run, in <4 minutes.
  • Explained the backward inv_fact sweep and why one pow suffices, in <1 minute.
  • LC 62 via formula and LC 96 via Catalan, each in <5 minutes including the mapping argument.
  • LC 920 solved in <35 minutes with modular DP.
  • lucas written and hand-verified on C(10, 4) mod 3 in <8 minutes.
  • Stated the composite-modulus escape hatch (CRT, phase-13 lab-07) and the Kummer carry rule without notes.