組合數(shù)的四種求法

前言：acwing 也做了這么久了，是時(shí)候開始復(fù)習(xí)了

0X00 四種方法以及相應(yīng)的代碼

方法 1 迭代法

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])

核心思想是： $C_i^{j} = C_{i-1}^{j-1} + C_{i-1}^{j}$

方法 2 簡單逆元

求組合數(shù) 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)

逆元的通俗解釋是： $a * b = a / b\ mod\ p$ ，也就是說除以一個(gè)數(shù)對 p 取余的結(jié)果，等于乘以一個(gè)數(shù)對 p 取余

在模數(shù) p 的情況下, b 存在逆元的充分必要條件是 p 是質(zhì)數(shù)

由于組合數(shù)可能很大，所以題目為了不超過 int 往往會(huì)讓答案對一個(gè)質(zhì)數(shù)取余數(shù)（比如 1e9 + 7）

我們回到組合數(shù)的基本算法：

$C_{a}^衬衬 = \frac{a!}{b!(a-b)!} = fact[a] * infact[a-b] * infact[b]$

其中 $fact[a]$ 表示 a 的階乘， $infact[b]$ 表示 b 的階乘的逆元并且 $fact[0] = infact[0] = 1$

方法3 lucas 算法

求組合數(shù) 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

除此之外代碼中有對 $C_{a}^$ 的化簡

$C_{a}^铐殃 = \frac{a!}{(a-b)!b!} = \frac{a(a-1)...(a-b+1)}{b!}$

方法 4 分解質(zhì)因數(shù)

求組合數(shù) 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)

$C_{a}^ = \frac{a!}{(a-b)!b!}$

再回到這個(gè)式子上跨新，我們知道任意一個(gè)數(shù)字 $n$ 都可以寫成 $n = p_1^{\alpha_1}p_2^{\alpha_2}...p_n^{\alpha_n}$ 也就是分解質(zhì)因數(shù)

所以我們只需要求出 $C_a^富腊$ 的所有質(zhì)因數(shù)以及他出現(xiàn)的次數(shù)就可以了

回到一個(gè)最簡單的例子 $n!$ 假設(shè)所有 <= n 的質(zhì)數(shù)存在 ps 中如何計(jì)算所有質(zhì)數(shù)出現(xiàn)的次數(shù)呢？（每個(gè)質(zhì)數(shù)出現(xiàn)的次數(shù)至少一次）

$\alpha = \lfloor \frac{n}{p} \rfloor + \lfloor \frac{n}{p^2} \rfloor...\lfloor \frac{n}{p^t} \rfloor$

0X01 總結(jié)

我們來總結(jié)一下四種方法：

底數(shù)比較小用第一種方法
有模數(shù)底數(shù)比較大用第二種
有模數(shù)域帐，但模數(shù)比較小赘被，底數(shù)很大用第三種
無模數(shù)是整，底數(shù)比較大，用第四種

?著作權(quán)歸作者所有,轉(zhuǎn)載或內(nèi)容合作請聯(lián)系作者