Lucas' Theorem

Introduction

This article discusses the computation of large binomial coefficients modulo a modulus. The binomial coefficient, also known as the combinatorial number, is defined as:

\[ \binom{n}{k} = \dfrac{n!}{k!(n-k)!}. \]

When the scale is not large, binomial coefficients can be computed using the recurrence formula with a time complexity of \(O(nk)\); or under a large prime modulus \(p>n\), they can be computed in \(O(n)\) time by calculating the factorials of the numerator and denominator. However, when the problem scale is very large (\(n\sim 10^{18}\)), these methods are no longer applicable.

Based on Lucas's theorem and its generalizations, this article discusses a method for computing binomial coefficients when the modulus is not too large (\(m \sim 10^6\)). More precisely, as long as the sum of prime powers in the unique factorization of the modulus \(m=\prod p_i^{e_i}\) (i.e., \(\sum p_i^{e_i}\)) is on the order of \(10^6\), this method can be used, because the preprocessing of the algorithm roughly corresponds to this scale.

Lucas's Theorem

First, consider the case where the modulus is a prime \(p\). In this case, Lucas's theorem applies:

Lucas's Theorem

For a prime \(p\), we have

\[ \binom{n}{k}\equiv \binom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor}\binom{n\bmod p}{k\bmod p}\pmod p. \]

Here, when \(n<k\), the binomial coefficient \(\dbinom{n}{k}\) is defined as \(0\).

Proof using generating functions

Consider the value of \(\displaystyle\binom{p}{n} \bmod p\). Since

\[ \binom{p}{n} = \frac{p!}{n!(p-n)!}, \]

when \(n\neq 0,p\), the denominator contains no factor \(p\), but the numerator contains a factor \(p\), so the fraction must be a multiple of \(p\), and its remainder modulo \(p\) is \(0\); when \(n=0,p\), the fraction is \(1\). Therefore,

\[ \binom{p}{n} \equiv [n=0\lor n=p] \pmod p. \]

Let \(f(x) = ax^n + bx^m\). In general, by the binomial expansion and Fermat's little theorem, we have

\[ \begin{aligned} (f(x))^p &= \left(ax^n + bx^m\right)^p \\ &= \sum_{k=0}^p\binom{p}{k}(ax^n)^k(bx^m)^{p-k}\\ &\equiv a^px^{pn} + b^px^{pm} \\ &\equiv a(x^p)^n+b(x^p)^m\\ &= f(x^p) \pmod p. \end{aligned} \]

The congruence in the third line uses the conclusion mentioned earlier, i.e., only when \(k=0,p\) is the binomial coefficient not a multiple of \(p\).

Using this conclusion, examine the binomial expansion:

\[ \begin{aligned} (1+x)^n &= (1+x)^{p\lfloor n/p\rfloor}(1+x)^{n\bmod p} \\ &\equiv (1+x^p)^{\lfloor n/p\rfloor}(1+x)^{n\bmod p} \pmod p. \end{aligned} \]

On the left-hand side, the coefficient of \(x^k\) is

\[ \binom{n}{k}\bmod p. \]

Now compute the coefficient of \(x^k\) on the right-hand side. The degrees of terms in the first factor must be multiples of \(p\), and the degrees of terms in the second factor must be less than \(p\). The way \(k\) is decomposed into the sum of such two parts is unique, i.e., division with remainder: \(k=p\lfloor k/p\rfloor +(k\bmod p)\). Therefore, the first factor can only contribute its term of degree \(p\lfloor k/p\rfloor\), and the second factor can only contribute its term of degree \(k\bmod p\). Thus, the coefficient of \(x^k\) in the right-hand side is the product of the coefficients of the terms contributed by the two factors:

\[ \binom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor}\binom{n\bmod p}{k\bmod p}\bmod p. \]

Equating the coefficients on both sides yields Lucas's theorem.

Proof using factorial modulo results

Here we provide a proof based on results related to factorial modulo, to facilitate connection with the method in the exLucas section later. The binomial coefficient is known as

\[ \binom{n}{k} = \dfrac{n!}{k!(n-k)!}. \]

Separating the powers of \(p\) from other factors in the factorial \(n!\), we obtain the decomposition:

\[ n! = p^{\nu_p(n!)}(n!)_p. \]

Then the binomial coefficient can be expressed as:

\[ \binom{n}{k} = p^{\nu_p(n!)-\nu_p(k!)-\nu_p((n-k)!)}\dfrac{(n!)_p}{(k!)_p((n-k)!)_p}. \]

The exponent \(\nu_p(n!)\) and the factorial remainder \((n!)_p\bmod p\) have recurrence formulas:

\[ \begin{aligned} \nu_p(n!) &= \lfloor n/p\rfloor+\nu_p( \lfloor n/p\rfloor!),\\ (n!)_p &\equiv (-1)^{\lfloor n/p\rfloor}\cdot (n\bmod p)!\cdot (\lfloor n/p\rfloor!)_p\pmod p. \end{aligned} \]

The former is a corollary of Legendre's formula, and the latter is a corollary of Wilson's theorem.

Substituting the recurrence formulas into the expression for the binomial coefficient and simplifying, we get:

\[ \begin{aligned} \binom{n}{k} &\equiv (-p)^{\lfloor n/p\rfloor-\lfloor k/p\rfloor-\lfloor(n-k)/p\rfloor}\cdot\dfrac{(n\bmod p)!}{(k\bmod p)!((n-k)\bmod p)!} \\ &\quad \cdot p^{\nu_p(\lfloor n/p\rfloor!)-\nu_p(\lfloor k/p\rfloor!)-\nu_p(\lfloor(n-k)/p\rfloor!)}\dfrac{(\lfloor n/p\rfloor!)_p}{(\lfloor k/p\rfloor!)_p(\lfloor(n-k)/p\rfloor!)_p} \pmod p. \end{aligned} \]

Now examine the value of \(\lfloor n/p\rfloor-\lfloor k/p\rfloor-\lfloor(n-k)/p\rfloor\). Since

\[ \begin{aligned} n &= \lfloor n/p\rfloor p + (n\bmod p),\\ k &= \lfloor k/p\rfloor p + (k\bmod p),\\ n-k &= \lfloor (n-k)/p\rfloor p + ((n-k)\bmod p),\\ \end{aligned} \]

subtracting the latter two from the first equation gives

\[ (\lfloor n/p\rfloor-\lfloor k/p\rfloor-\lfloor(n-k)/p\rfloor)p = (k\bmod p)+((n-k)\bmod p)-(n\bmod p). \]

On the right-hand side, the sum of the first two terms is strictly less than \(2p\), and the third term \(n\bmod p\) is exactly the remainder of the sum of the first two terms. Therefore, the right-hand side must be non-negative but less than \(2p\), and it must be a multiple of \(p\), so it can only be \(0\) or \(p\). This implies that \(\lfloor n/p\rfloor-\lfloor k/p\rfloor-\lfloor(n-k)/p\rfloor\) can only be \(0\) or \(1\):

• If it is \(0\), then we also have \((n\bmod p) = (k\bmod p)+((n-k)\bmod p)\). Therefore, the exponent of the first factor in the above expression is \(0\), so that factor equals one; the second factor is \(\dbinom{n\bmod p}{k\bmod p}\); and the third factor, as seen from the earlier expansion, equals \(\dbinom{\lfloor n/p\rfloor}{\lfloor k/p\rfloor}\). In this case, Lucas's formula holds;

• If it is \(1\), then the exponent of the first factor is \(1\), so that factor equals zero, meaning the remainder of the binomial coefficient is zero. At the same time, \(\dbinom{n\bmod p}{k\bmod p}\) on the right-hand side of Lucas's theorem must also be zero, because in this case we necessarily have \((n\bmod p)<(k\bmod p)\); otherwise, we would have

  \[ ((n-k)\bmod p) = p + (n\bmod p) - (k\bmod p) \ge p. \]

  This clearly contradicts the definition of a remainder.

Combining the two cases, we obtain the desired Lucas's theorem. This proof shows that when computing binomial coefficients modulo a prime, the results obtained using Lucas's theorem and the exLucas algorithm are equivalent.

Lucas's theorem states that when the modulus is a prime \(p\), the computation of large binomial coefficients can be reduced to computing smaller binomial coefficients. On the right-hand side, the first binomial coefficient can be recursively computed until \(n,k<p\); the second binomial coefficient can be computed directly or preprocessed in advance. In code form, it looks like:

Illustration

long long Lucas(long long n, long long k, long long p) {
  if (k == 0) return 1;
  return (C(n % p, k % p, p) * Lucas(n / p, k / p, p)) % p;
}

Here, C(n, k, p) is used to compute small-scale binomial coefficients.

Recursion proceeds at most \(O(\log_p n)\) times, so the algorithm's complexity is \(O(f(p)+g(p)\log_p n)\), where \(f(p)\) is the complexity of preprocessing binomial coefficients, and \(g(p)\) is the complexity of computing a single binomial coefficient.

Reference Implementation

The reference implementation provided here preprocesses factorials and their inverses up to \(p\) in \(O(p)\) time, allowing computation of a single binomial coefficient in \(O(1)\) time:

Reference Implementation

#include <iostream>
#include <vector>

class BinomModPrime {
  int p;
  std::vector<int> fa, ifa;

  // Calculate binom(n, k) mod p for n, k < p.
  int calc(int n, int k) {
    if (n < k) return 0;
    long long res = fa[n];
    res = (res * ifa[k]) % p;
    res = (res * ifa[n - k]) % p;
    return res;
  }

 public:
  BinomModPrime(int p) : p(p), fa(p), ifa(p) {
    // Factorials mod p till p.
    fa[0] = 1;
    for (int i = 1; i < p; ++i) {
      fa[i] = (long long)fa[i - 1] * i % p;
    }
    // Inverse of factorials mod p till p.
    ifa[p - 1] = p - 1;  // Wilson's theorem.
    for (int i = p - 1; i; --i) {
      ifa[i - 1] = (long long)ifa[i] * i % p;
    }
  }

  // Calculate binom(n, k) mod p.
  int binomial(long long n, long long k) {
    long long res = 1;
    while (n || k) {
      res = (res * calc(n % p, k % p)) % p;
      n /= p;
      k /= p;
    }
    return res;
  }
};

int main() {
  int t, p;
  std::cin >> t >> p;
  BinomModPrime bm(p);
  for (; t; --t) {
    long long n, k;
    std::cin >> n >> k;
    std::cout << bm.binomial(n, k) << '\n';
  }
  return 0;
}

The time complexity of this implementation is \(O(p+T\log_p n)\), where \(T\) is the number of queries.

exLucas Algorithm

Lucas's theorem requires the modulus \(p\) to be prime. For cases where \(p\) is not prime, the exLucas algorithm is needed. Despite its name, the algorithm does not actually use Lucas's theorem in its operation. Its key step is computing factorials modulo prime powers. The second proof above indicates its connection with Lucas's theorem.

Prime Power Modulus Case

First, consider the case where the modulus is a prime power \(p^\alpha\). Separating the powers of \(p\) from other factors in the factorial \(n!\), we obtain the decomposition:

\[ n! = p^{\nu_p(n!)}(n!)_p. \]

Here, \(\nu_p(n!)\) is the exponent of \(p\) in the prime factorization of \(n!\), and \((n!)_p\) is clearly coprime with \(p\). Therefore, the binomial coefficient can be written as:

\[ \binom{n}{k} = p^{\nu_p(n!)-\nu_p(k!)-\nu_p((n-k)!)}\dfrac{(n!)_p}{(k!)_p((n-k)!)_p}. \]

The terms \(\nu_p(n!)\), etc., can be computed using Legendre's formula, and \((n!)_p\), etc., can be computed using recurrence relations. Since the latter are coprime with \(p^\alpha\), the inverse of the product in the denominator can be computed via the extended Euclidean algorithm. Thus, the problem is solved.

Note that if the exponent \(\nu_p(n!)-\nu_p(k!)-\nu_p((n-k)!)\ge\alpha\), the remainder must be zero, and no further computation is needed.

General Modulus Case

For the case where \(m\) is a general composite number, we first perform prime factorization:

\[ m = p_1^{\alpha_1}p_2^{\alpha_2}\cdots p_s^{\alpha_s}. \]

Then, compute the remainder of the binomial coefficient \(\dbinom{n}{k}\) modulo each \(p_i^{\alpha_i}\), yielding \(s\) congruence equations:

\[ \begin{cases} \dbinom{n}{k} \equiv r_1, &\pmod{p_1^{\alpha_1}}, \\ \dbinom{n}{k} \equiv r_2, &\pmod{p_2^{\alpha_2}}, \\ \quad\quad\cdots\\ \dbinom{n}{k} \equiv r_s, &\pmod{p_s^{\alpha_s}}. \end{cases} \]

Finally, use the Chinese Remainder Theorem to find the remainder modulo \(m\).

Reference Implementation

Finally, we provide a reference implementation for the template problem Binomial Coefficient.

Reference Implementation

#include <iostream>
#include <vector>

// Extended Euclid.
void ex_gcd(int a, int b, int& x, int& y) {
  if (!b) {
    x = 1;
    y = 0;
  } else {
    ex_gcd(b, a % b, y, x);
    y -= a / b * x;
  }
}

// Inverse of a mod m.
int inverse(int a, int m) {
  int x, y;
  ex_gcd(a, m, x, y);
  return (x % m + m) % m;
}

// Coefficient in CRT.
int crt_coeff(int m_i, int m) {
  long long mm = m / m_i;
  mm *= inverse(mm, m_i);
  return mm % m;
}

// Binominal Coefficient Calculator Modulo Prime Power.
class BinomModPrimePower {
  int p, a, pa;
  std::vector<int> f;

  // Obtain multiplicity of p in n!.
  long long nu(long long n) {
    long long count = 0;
    do {
      n /= p;
      count += n;
    } while (n);
    return count;
  }

  // Calculate (n!)_p mod pa.
  long long fact_mod(long long n) {
    bool neg = p != 2 || pa <= 4;
    long long res = 1;
    while (n > 1) {
      if ((n / pa) & neg) res = pa - res;
      res = res * f[n % pa] % pa;
      n /= p;
    }
    return res;
  }

 public:
  BinomModPrimePower(int p, int a, int pa) : p(p), a(a), pa(pa), f(pa) {
    // Pretreatment.
    f[0] = 1;
    for (int i = 1; i < pa; ++i) {
      f[i] = i % p ? (long long)f[i - 1] * i % pa : f[i - 1];
    }
  }

  // Calculate Binom(n, k) mod pa.
  int binomial(long long n, long long k) {
    long long v = nu(n) - nu(n - k) - nu(k);
    if (v >= a) return 0;
    auto res = fact_mod(n - k) * fact_mod(k) % pa;
    res = fact_mod(n) * inverse(res, pa) % pa;
    for (; v; --v) res *= p;
    return res % pa;
  }
};

// Binominal Coefficient Calculator.
class BinomMod {
  int m;
  std::vector<BinomModPrimePower> bp;
  std::vector<long long> crt_m;

 public:
  BinomMod(int n) : m(n) {
    // Factorize.
    for (int p = 2; p * p <= n; ++p) {
      if (n % p == 0) {
        int a = 0, pa = 1;
        for (; n % p == 0; n /= p, ++a, pa *= p);
        bp.emplace_back(p, a, pa);
        crt_m.emplace_back(crt_coeff(pa, m));
      }
    }
    if (n > 1) {
      bp.emplace_back(n, 1, n);
      crt_m.emplace_back(crt_coeff(n, m));
    }
  }

  // Calculate Binom(n, k) mod m.
  int binomial(long long n, long long k) {
    long long res = 0;
    for (size_t i = 0; i != bp.size(); ++i) {
      res = (bp[i].binomial(n, k) * crt_m[i] + res) % m;
    }
    return res;
  }
};

int main() {
  int t, m;
  std::cin >> t >> m;
  BinomMod bm(m);
  for (; t; --t) {
    long long n, k;
    std::cin >> n >> k;
    std::cout << bm.binomial(n, k) << '\n';
  }
  return 0;
}

In preprocessing, this algorithm decomposes the modulus \(m\) into prime powers, then for each \(p^\alpha\), precomputes the product of all natural numbers from \(1\) to \(p^\alpha\) that are not multiples of \(p\), as well as the corresponding coefficients for merging answers via the Chinese Remainder Theorem. The preprocessing time complexity is \(O(\sqrt{m}+\sum_ip_i^{\alpha_i})\). For each query, the complexity is \(O(\log m+\sum_i\log_{p_i}n)\), where the two terms correspond to computing inverses and computing exponents and factorial remainders, respectively.

Exercises

• Luogu3807【Template】Lucas's Theorem
• SDOI2010 Ancient Pig Script & Lucas's Theorem
• Luogu4720【Template】Extended Lucas
• Ceizenpok’s formula

---

本页面最近更新：，更新历史
发现错误？想一起完善？ 在 GitHub 上编辑此页！
本页面贡献者：OI-wiki
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用