组合数的四种求法

前言：acwing 也做了这么久了，是时候开始复习了

0X00 四种方法以及相应的代码

方法 1 迭代法

acwing 求组合数 I

from sys import stdin

def read(n=2):
    if n == 1:
        return int(stdin.readline())
    return map(int, stdin.readline().split())

mod = 10 ** 9 + 7

N = 2010
n = read(1)
comb = [[0] * N for _ in range(N)]

for i in range(N):
    for j in range(i+1):
        if j == 0:
            comb[i][j] = 1
        else:
            comb[i][j] = (comb[i-1][j]+comb[i-1][j-1]) % mod

while n > 0:
    n -= 1
    a, b = read()
    print(comb[a][b])

核心思想是：

方法 2 简单逆元

求组合数 II

from sys import stdin

def read(n=2):
    if n == 1:
        return int(stdin.readline())
    return map(int, stdin.readline().split())

N = 100010
fact = [0] * N
infact = [0] * N

fact[0] = infact[0] = 1

mod = 10 ** 9 + 7

def qmi(a, k, p):
    res = 1
    while k > 0:
        if k & 1 == 1:
            res = res * a % mod
        k >>= 1
        a = a * a % mod
    return res

for i in range(1, N):
    fact[i] = fact[i-1] * i % mod
    infact[i] = infact[i-1] * qmi(i, mod-2, mod) % mod

n = read(1)
while n > 0:
    n -= 1
    a, b = read()
    print(fact[a] * infact[a - b] * infact[b] % mod)
   

逆元的通俗解释是：， 也就是说除以一个数对 p 取余的结果，等于乘以一个数对 p 取余

在模数 p 的情况下, b 存在逆元的充分必要条件是 p 是质数

由于组合数可能很大，所以题目为了不超过 int 往往会让答案对一个质数取余数（比如 1e9 + 7）

我们回到组合数的基本算法：

其中 表示 a 的阶乘， 表示 b 的阶乘的逆元并且

方法3 lucas 算法

求组合数 III

from sys import stdin 

def read(n = 2):
    if n == 1:
        return int(stdin.readline())
    return map(int, stdin.readline().split())

def qmi(a, k, p):
    res = 1
    while k > 0:
        if k & 1:
            res = res * a % p
        k >>= 1
        a = a * a % p
    
    return res
            
    

def C(a, b, p):
    if b > a: return 0
    x, y = 1, 1
    j = 1
    for i in range(a, a-b, -1):
        x = x * i % p
        y = y * j % p
        j += 1
    # print("x, y", x, y)
    return x * qmi(y, p-2, p) % p 
    
def lucas(a, b, p):
    if a < p and b < p:
        return C(a, b, p)
    return C(a % p, b % p, p) * lucas(a // p, b // p, p) % p

n = read(1)
while n > 0:
    n -= 1
    a, b, p = read()
    print(lucas(a, b, p))
    

摘自 oi wiki

除此之外代码中有对 的化简

方法 4 分解质因数

求组合数 IV

N = 5010
cnt = 0
primes = [0] * N
st = [True] * N

def getPrimes(n):
    global cnt
    for i in range(2, n+1):
        if st[i]:
            primes[cnt] = i
            cnt += 1
        j = 0
        while primes[j] <= n // i:
            st[primes[j] * i] = False
            if not (i % primes[j]):
                break
            j += 1

def get(a, p):
    
    res = 0
    while a > 0:
        res += a // p
        a //= p
    
    return res
    
a, b = map(int, input().split())
getPrimes(a)
sums = [0] * cnt
for i in range(cnt):
    p = primes[i]
    sums[i] = get(a, p) - get(a-b, p) - get(b, p)

res = 1
for t, v in zip(sums, primes):
    res *= v ** t
print(res)

再回到这个式子上，我们知道任意一个数字 都可以写成 也就是分解质因数

所以我们只需要求出 的所有质因数以及他出现的次数就可以了

回到一个最简单的例子 假设所有 <= n 的质数存在 ps 中如何计算所有质数出现的次数呢？（每个质数出现的次数至少一次）

0X01 总结

我们来总结一下四种方法：

• 底数比较小用第一种方法
• 有模数底数比较大用第二种
• 有模数，但模数比较小，底数很大用第三种
• 无模数，底数比较大，用第四种