Appendix E: Math Foundations for Competitive Programming

💡 About This Appendix: Competitive programming often requires mathematical tools beyond basic arithmetic. This appendix covers the essential math you'll encounter in USACO Bronze, Silver, and Gold — with contest-ready code templates for each topic.

E.1 Modular Arithmetic

Why Do We Need Modular Arithmetic?

Many problems ask you to output an answer "modulo 10⁹ + 7". This isn't arbitrary — it prevents integer overflow when answers are astronomically large.

Consider: "How many permutations of N elements?" Answer: N! For N = 20, that's 2,432,902,008,176,640,000 — larger than long long's max (~9.2 × 10¹⁸). For N = 100, it's completely unrepresentable.

Solution: Compute everything modulo a prime M (typically 10⁹ + 7).

Common MOD Values

Basic Modular Operations Template

// Solution: Modular Arithmetic Basics
#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
const ll MOD = 1e9 + 7;  // standard competitive programming MOD

// Safe addition: (a + b) % MOD
ll addMod(ll a, ll b) {
    return (a % MOD + b % MOD) % MOD;
}

// Safe subtraction: (a - b + MOD) % MOD (handle negative result)
ll subMod(ll a, ll b) {
    return ((a % MOD) - (b % MOD) + MOD) % MOD;  // +MOD prevents negative!
}

// Safe multiplication: (a * b) % MOD
// Key: a and b are at most MOD-1 ≈ 10^9, so a*b ≈ 10^18 which fits long long
ll mulMod(ll a, ll b) {
    return (a % MOD) * (b % MOD) % MOD;
}

// Example: Compute sum of first N integers modulo MOD
ll sumFirstN(ll n) {
    // Formula: n*(n+1)/2, but careful with division — need modular inverse!
    // For now: just accumulate with addMod
    ll result = 0;
    for (ll i = 1; i <= n; i++) {
        result = addMod(result, i);
    }
    return result;
}

⚠️ Critical Bug: (a - b) % MOD can be negative in C++ if a < b! Always use (a - b + MOD) % MOD.

E.1.1 Fast Exponentiation (Binary Exponentiation)

Computing a^n mod M naively takes O(N) multiplications. Fast exponentiation (exponentiation by squaring) does it in O(log N).

Key insight: a^n = a^(n/2) × a^(n/2)          if n is even
              a^n = a × a^((n-1)/2) × a^((n-1)/2)  if n is odd

Example: a^13 = a^(1101 in binary)
       = a^8 × a^4 × a^1
       = 3 multiplications instead of 12!

// Solution: Fast Modular Exponentiation — O(log n)
// Computes (base^exp) % mod
ll power(ll base, ll exp, ll mod = MOD) {
    ll result = 1;
    base %= mod;                  // reduce base first
    
    while (exp > 0) {
        if (exp & 1) {            // if current bit is 1
            result = result * base % mod;
        }
        base = base * base % mod; // square the base
        exp >>= 1;                // shift to next bit
    }
    return result;
}

// Example usage:
// power(2, 10) = 1024 % MOD = 1024
// power(2, 100, MOD) = 2^100 mod (10^9+7)

E.1.2 Modular Inverse (Fermat's Little Theorem)

The modular inverse of a modulo M is a number a⁻¹ such that a × a⁻¹ ≡ 1 (mod M).

This lets us do modular division: a / b mod M = a × b⁻¹ mod M.

Fermat's Little Theorem: If M is prime and gcd(a, M) = 1, then:

// Solution: Modular Inverse using Fermat's Little Theorem
// Only works when MOD is PRIME and gcd(a, MOD) = 1
ll modInverse(ll a, ll mod = MOD) {
    return power(a, mod - 2, mod);
}

// Division with modular arithmetic:
ll divMod(ll a, ll b) {
    return mulMod(a, modInverse(b));
}

// Example: (n! / k!) mod MOD
// = n! × (k!)^(-1) mod MOD
// = n! × modInverse(k!) mod MOD

E.1.3 Precomputing Factorials and Inverses

For problems requiring many combinations C(n, k):

// Solution: Precomputed Factorials for O(1) Combination Queries
const int MAXN = 1000005;
ll fact[MAXN], inv_fact[MAXN];

void precompute() {
    fact[0] = 1;
    for (int i = 1; i < MAXN; i++) {
        fact[i] = fact[i-1] * i % MOD;
    }
    inv_fact[MAXN-1] = modInverse(fact[MAXN-1]);
    for (int i = MAXN-2; i >= 0; i--) {
        inv_fact[i] = inv_fact[i+1] * (i+1) % MOD;
    }
}

// C(n, k) = n! / (k! * (n-k)!)
ll C(int n, int k) {
    if (k < 0 || k > n) return 0;
    return fact[n] * inv_fact[k] % MOD * inv_fact[n-k] % MOD;
}

// Usage: precompute() once, then C(n, k) in O(1)

E.2 GCD and LCM

Euclidean Algorithm

The Greatest Common Divisor (GCD) of two numbers is the largest number that divides both.

Euclidean Algorithm: Based on gcd(a, b) = gcd(b, a % b).

// Solution: GCD — O(log(min(a,b)))
int gcd(int a, int b) {
    while (b != 0) {
        a %= b;
        swap(a, b);
    }
    return a;
}
// Or recursively:
// int gcd(int a, int b) { return b == 0 ? a : gcd(b, a % b); }

// C++17: std::gcd from <numeric>
// int g = gcd(a, b);           // std::gcd, C++17 (recommended)
// int g = __gcd(a, b);         // legacy GCC built-in, still works

Trace: gcd(48, 18):

gcd(48, 18) → gcd(18, 48%18=12) → gcd(12, 18%12=6) → gcd(6, 0) = 6

LCM and the Overflow Trap

// Solution: LCM — be careful with overflow!

// WRONG: overflows for large a, b
long long lcmWrong(long long a, long long b) {
    return a * b / gcd(a, b);  // a*b can overflow even long long!
}

// CORRECT: divide first, then multiply
long long lcm(long long a, long long b) {
    return a / gcd(a, b) * b;  // divide BEFORE multiplying
}
// a / gcd(a,b) is always an integer, so no precision loss
// Then * b: max value is around 10^18 which fits in long long

⚠️ Always divide before multiplying to avoid overflow!

Extended Euclidean Algorithm

Finds integers x, y such that ax + by = gcd(a, b) — useful for modular inverse when MOD is not prime:

// Solution: Extended Euclidean Algorithm — O(log(min(a,b)))
// Returns gcd(a,b), and sets x,y such that a*x + b*y = gcd(a,b)
long long extgcd(long long a, long long b, long long &x, long long &y) {
    if (b == 0) { x = 1; y = 0; return a; }
    long long x1, y1;
    long long g = extgcd(b, a % b, x1, y1);
    x = y1;
    y = x1 - (a / b) * y1;
    return g;
}

// Modular inverse using extgcd (works even when MOD is not prime):
long long modInverseExtGcd(long long a, long long mod) {
    long long x, y;
    long long g = extgcd(a, mod, x, y);
    if (g != 1) return -1;  // no inverse exists (gcd != 1)
    return (x % mod + mod) % mod;
}

E.3 Prime Numbers and Sieves

Trial Division

// Solution: Trial Division Primality Test — O(sqrt(N))
bool isPrime(long long n) {
    if (n < 2) return false;
    if (n == 2) return true;
    if (n % 2 == 0) return false;
    for (long long i = 3; i * i <= n; i += 2) {
        if (n % i == 0) return false;
    }
    return true;
}
// Efficient because: if n has a factor > sqrt(n), it must also have one <= sqrt(n)
// Only check odd numbers after 2 (halves the iterations)

Sieve of Eratosthenes

Find all primes up to N efficiently:

// Solution: Sieve of Eratosthenes — O(N log log N) time, O(N) space
// After running, isPrime[i] = true iff i is prime
const int MAXN = 1000005;
bool isPrime[MAXN];

void sieve(int n) {
    fill(isPrime, isPrime + n + 1, true);  // assume all prime initially
    isPrime[0] = isPrime[1] = false;        // 0 and 1 are not prime
    
    for (int i = 2; (long long)i * i <= n; i++) {
        if (isPrime[i]) {
            // Mark all multiples of i as composite
            for (int j = i * i; j <= n; j += i) {
                isPrime[j] = false;
                // Start from i*i (smaller multiples already marked by smaller primes)
            }
        }
    }
}

// Count primes up to N:
void countPrimes(int n) {
    sieve(n);
    int count = 0;
    for (int i = 2; i <= n; i++) {
        if (isPrime[i]) count++;
    }
    cout << count << "\n";
}

Why start inner loop at i²? All multiples of i smaller than i² (i.e., 2i, 3i, ..., (i-1)i) were already marked by smaller primes (2, 3, ..., i-1).

Linear Sieve (Euler Sieve) — O(N)

The Euler sieve marks each composite number exactly once:

// Solution: Linear Sieve (Euler Sieve) — O(N) time
// Also computes smallest prime factor (SPF) for each number
const int MAXN = 1000005;
int spf[MAXN];      // smallest prime factor
vector<int> primes;

void linearSieve(int n) {
    fill(spf, spf + n + 1, 0);
    for (int i = 2; i <= n; i++) {
        if (spf[i] == 0) {          // i is prime
            spf[i] = i;
            primes.push_back(i);
        }
        for (int j = 0; j < (int)primes.size() && primes[j] <= spf[i] && (long long)i * primes[j] <= n; j++) {
            spf[i * primes[j]] = primes[j];  // mark composite
        }
    }
}

// Fast prime factorization using SPF:
// O(log N) per factorization
vector<int> factorize(int n) {
    vector<int> factors;
    while (n > 1) {
        factors.push_back(spf[n]);
        n /= spf[n];
    }
    return factors;
}

E.4 Binary Representations and Bit Manipulation

Fundamental Bit Operations

// Solution: Common Bit Operations Reference
int n = 42;   // binary: 101010

// ── AND (&): both bits must be 1 ──
int a = 6 & 3;     // 110 & 011 = 010 = 2

// ── OR (|): at least one bit is 1 ──
int b = 6 | 3;     // 110 | 011 = 111 = 7

// ── XOR (^): exactly one bit is 1 ──
int c = 6 ^ 3;     // 110 ^ 011 = 101 = 5

// ── NOT (~): flip all bits (two's complement) ──
int d = ~6;        // = -7 (in two's complement)

// ── Left shift (<<): multiply by 2^k ──
int e = 1 << 4;    // = 16 = 2^4

// ── Right shift (>>): divide by 2^k (arithmetic) ──
int f = 32 >> 2;   // = 8 = 32/4

Essential Bit Tricks

// Solution: Competitive Programming Bit Tricks

// ── Check if n is odd ──
bool isOdd(int n) { return n & 1; }  // last bit is 1 iff odd

// ── Check if n is a power of 2 ──
bool isPow2(int n) { return n > 0 && (n & (n-1)) == 0; }
// Why? Powers of 2: 1=001, 2=010, 4=100. n-1 flips all lower bits.
// 4 & 3 = 100 & 011 = 000. Non-powers: 6 & 5 = 110 & 101 = 100 ≠ 0.

// ── Get k-th bit (0-indexed from right) ──
bool getBit(int n, int k) { return (n >> k) & 1; }

// ── Set k-th bit to 1 ──
int setBit(int n, int k) { return n | (1 << k); }

// ── Clear k-th bit (set to 0) ──
int clearBit(int n, int k) { return n & ~(1 << k); }

// ── Toggle k-th bit ──
int toggleBit(int n, int k) { return n ^ (1 << k); }

// ── lowbit: lowest set bit (used in Fenwick tree!) ──
int lowbit(int n) { return n & (-n); }
// Example: lowbit(12) = lowbit(1100) = 0100 = 4

// ── Count number of set bits (popcount) ──
int popcount(int n) { return __builtin_popcount(n); }   // use built-in!
// For long long: __builtin_popcountll(n)

// ── Swap two numbers without temp variable ──
void swapXOR(int &a, int &b) {
    a ^= b;
    b ^= a;
    a ^= b;
}
// (usually just use std::swap — this is mainly a curiosity)

// ── Find position of lowest set bit ──
int lowestBitPos(int n) { return __builtin_ctz(n); }  // count trailing zeros
// __builtin_clz(n) = count leading zeros

Subset Enumeration

A powerful technique: enumerate all subsets of a set represented as a bitmask.

// Solution: Subset Enumeration with Bitmasks
// Enumerate all subsets of an N-element set

void enumerateAllSubsets(int n) {
    // Total subsets = 2^n
    for (int mask = 0; mask < (1 << n); mask++) {
        // 'mask' represents a subset: bit i set = element i is included
        cout << "Subset: {";
        for (int i = 0; i < n; i++) {
            if (mask & (1 << i)) {
                cout << i << " ";
            }
        }
        cout << "}\n";
    }
}

// Enumerate all NON-EMPTY subsets of a given set 'S'
void enumerateSubsetsOf(int S) {
    for (int sub = S; sub > 0; sub = (sub - 1) & S) {
        // Process subset 'sub'
        // The trick: (sub-1) & S gives the "next smaller" subset of S
        // This enumerates all 2^|S| subsets of S in O(1) amortized per step
    }
}

// Classic use: bitmask DP
// dp[mask] = minimum cost to visit the set of cities represented by mask
// dp[0] = 0 (start: no cities visited)
// dp[mask | (1 << v)] = min(dp[mask | (1 << v)], dp[mask] + cost[last][v])

E.5 Combinatorics Basics

Counting Formulas

// Solution: Combinatorics with Modular Arithmetic
// Assumes precompute() from E.1.3 has been called

// C(n, k) = n! / (k! * (n-k)!)
ll combination(int n, int k) {
    if (k < 0 || k > n) return 0;
    return fact[n] * inv_fact[k] % MOD * inv_fact[n-k] % MOD;
}

// P(n, k) = n! / (n-k)!
ll permutation(int n, int k) {
    if (k < 0 || k > n) return 0;
    return fact[n] * inv_fact[n-k] % MOD;
}

// Stars and Bars: number of ways to put n identical balls into k distinct boxes
// = C(n + k - 1, k - 1)
ll starsAndBars(int n, int k) {
    return combination(n + k - 1, k - 1);
}

Pascal's Triangle — Computing C(n, k) without Precomputation

When n is small (n ≤ 2000), Pascal's triangle is simpler:

// Solution: Pascal's Triangle DP — O(n^2) precomputation
const int MAXN = 2005;
ll C[MAXN][MAXN];

void buildPascal() {
    for (int i = 0; i < MAXN; i++) {
        C[i][0] = C[i][i] = 1;
        for (int j = 1; j < i; j++) {
            C[i][j] = (C[i-1][j-1] + C[i-1][j]) % MOD;
        }
    }
}
// Then C[n][k] is the answer for any 0 <= k <= n < MAXN
// This avoids modular inverse entirely — useful when MOD might not be prime

Pascal's Rule: C(n, k) = C(n-1, k-1) + C(n-1, k)

This comes from: "choose k items from n" = "include item n and choose k-1 from n-1" + "exclude item n and choose k from n-1".

Key Combinatorial Identities

// Useful identities in competitive programming:

// Hockey Stick Identity: sum of C(r+k, k) for k=0..n = C(n+r+1, n)
// Useful for: 2D prefix sums, polynomial evaluations

// Vandermonde's Identity: sum_k C(m,k)*C(n,r-k) = C(m+n, r)
// Useful for: counting problems with two groups

// Inclusion-Exclusion:
// |A ∪ B| = |A| + |B| - |A ∩ B|
// |A ∪ B ∪ C| = |A| + |B| + |C| - |A∩B| - |A∩C| - |B∩C| + |A∩B∩C|
// Generalizes to n sets with 2^n terms (or bitmask enumeration)

E.6 Common Mathematical Results for Complexity Analysis

Harmonic Series

This explains why the Sieve of Eratosthenes runs in O(N log log N):

• Total work = N/2 + N/3 + N/5 + N/7 + ... (for each prime p, mark N/p multiples)
• Sum over primes ≈ N × ln(ln(N))

And why Fenwick tree operations are O(log N): the lowbit operation advances by 1, 2, 4, ... bits.

Key Estimates

Operations Per Second Estimate

Using this, you can estimate if your algorithm is fast enough:

• N = 10⁵, O(N log N) → ~1.7 × 10⁶ ops → fast
• N = 10⁵, O(N²) → 10¹⁰ ops → too slow
• N = 10⁵, O(N√N) → ~3 × 10⁷ ops → borderline (usually OK with 2s limit)

E.7 Complete Math Template

Here's a single file with all the templates from this appendix:

// Solution: Complete Math Template for Competitive Programming
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;

// ═══════════════════════════════════════════════
// MODULAR ARITHMETIC
// ═══════════════════════════════════════════════
const ll MOD = 1e9 + 7;

ll power(ll base, ll exp, ll mod = MOD) {
    ll result = 1;
    base %= mod;
    while (exp > 0) {
        if (exp & 1) result = result * base % mod;
        base = base * base % mod;
        exp >>= 1;
    }
    return result;
}

ll modInverse(ll a, ll mod = MOD) {
    return power(a, mod - 2, mod);
}

// ═══════════════════════════════════════════════
// FACTORIALS (precomputed up to MAXN)
// ═══════════════════════════════════════════════
const int MAXN = 1000005;
ll fact[MAXN], inv_fact[MAXN];

void precomputeFactorials() {
    fact[0] = 1;
    for (int i = 1; i < MAXN; i++) fact[i] = fact[i-1] * i % MOD;
    inv_fact[MAXN-1] = modInverse(fact[MAXN-1]);
    for (int i = MAXN-2; i >= 0; i--) inv_fact[i] = inv_fact[i+1] * (i+1) % MOD;
}

ll C(int n, int k) {
    if (k < 0 || k > n) return 0;
    return fact[n] * inv_fact[k] % MOD * inv_fact[n-k] % MOD;
}

// ═══════════════════════════════════════════════
// GCD / LCM
// ═══════════════════════════════════════════════
ll gcd(ll a, ll b) { return b == 0 ? a : gcd(b, a % b); }
ll lcm(ll a, ll b)  { return a / gcd(a, b) * b; }

// ═══════════════════════════════════════════════
// PRIME SIEVE
// ═══════════════════════════════════════════════
const int MAXP = 1000005;
bool notPrime[MAXP];
vector<int> primes;

void sieve(int n = MAXP - 1) {
    notPrime[0] = notPrime[1] = true;
    for (int i = 2; i <= n; i++) {
        if (!notPrime[i]) {
            primes.push_back(i);
            for (long long j = (long long)i*i; j <= n; j += i)
                notPrime[j] = true;
        }
    }
}

bool isPrime(int n) { return n >= 2 && !notPrime[n]; }

// ═══════════════════════════════════════════════
// BIT TRICKS
// ═══════════════════════════════════════════════
bool isOdd(int n)       { return n & 1; }
bool isPow2(int n)      { return n > 0 && !(n & (n-1)); }
int  lowbit(int n)      { return n & (-n); }
int  popcount(int n)    { return __builtin_popcount(n); }
int  ctz(int n)         { return __builtin_ctz(n); }  // count trailing zeros

// ═══════════════════════════════════════════════
// EXTENDED GCD
// ═══════════════════════════════════════════════
ll extgcd(ll a, ll b, ll &x, ll &y) {
    if (!b) { x = 1; y = 0; return a; }
    ll x1, y1, g = extgcd(b, a%b, x1, y1);
    x = y1; y = x1 - a/b * y1;
    return g;
}

int main() {
    ios_base::sync_with_stdio(false);
    cin.tie(NULL);
    
    precomputeFactorials();
    sieve();
    
    // Test: C(10, 3) = 120
    cout << C(10, 3) << "\n";
    
    // Test: 2^100 mod (10^9+7)
    cout << power(2, 100) << "\n";
    
    // Test: first few primes
    for (int i = 0; i < 10; i++) cout << primes[i] << " ";
    cout << "\n";
    
    return 0;
}

E.8 Number Theory Quick Reference

Divisibility Rules (useful for manual checks)

Integer Square Root

// Safe integer square root (avoids floating point errors)
ll isqrt(ll n) {
    ll x = sqrtl(n);              // floating point approximation
    while (x * x > n) x--;        // correct downward if needed
    while ((x+1) * (x+1) <= n) x++; // correct upward if needed
    return x;
}

Ceiling Division

// Ceiling division: ceil(a/b) for positive integers
ll ceilDiv(ll a, ll b) {
    return (a + b - 1) / b;
    // Or: (a - 1) / b + 1  (same thing for a > 0)
}

❓ FAQ

Q1: When should I use long long?

A: When values might exceed 2 × 10⁹ (roughly the int limit). Typical cases: ① multiplying two large int values (10⁹ × 10⁹ = 10¹⁸); ② summing path weights (N edges, each weight 10⁶, total up to 10¹¹); ③ factorials/combinations (use long long for intermediate calculations even with modular arithmetic). Rule of thumb: use long long whenever there's multiplication in competitive programming code.

Q2: Why use 10⁹ + 7 as the modulus instead of 10⁹?

A: 10⁹ is not prime (= 2⁹ × 5⁹), so Fermat's little theorem can't be used to compute modular inverses. 10⁹ + 7 = 1,000,000,007 is prime, and (10⁹ + 7)² < 2⁶³ (the long long limit), so multiplying two numbers after taking the modulus won't overflow long long.

Q3: How does the bit-manipulation trick in fast exponentiation work?

A: Write the exponent n in binary: n = b_k × 2^k + ... + b_1 × 2 + b_0. Then a^n = a^(b_k × 2^k) × ... × a^(b_1 × 2) × a^b_0. Each loop iteration squares the base (representing a to the power of 2^k), and multiplies into the result when the current bit is 1. This requires only log₂(n) multiplications.

Q4: Why does the Sieve of Eratosthenes start marking from i×i?

A: Multiples 2i, 3i, ..., (i-1)i have already been marked by the smaller primes 2, 3, ..., i-1. For example, 6 = 2×3 was marked by 2; 7×5=35 was marked by 5. Starting from i×i avoids redundant work and optimizes the constant factor.

Q5: Why does n & (n-1) check if n is a power of 2?

A: Powers of 2 have exactly one 1-bit in binary (e.g., 8 = 1000). Subtracting 1 flips the lowest 1-bit to 0 and all lower 0-bits to 1 (e.g., 7 = 0111). So n & (n-1) clears the lowest 1-bit. If n is a power of 2 (only one 1-bit), the result is 0; otherwise it's nonzero.

End of Appendix E — See also: Algorithm Templates | Competitive Programming Tricks