组合数实现及卢卡斯定理证明

最新推荐文章于 2024-02-20 01:19:51 发布

原创最新推荐文章于 2024-02-20 01:19:51 发布 · 293 阅读

2 ·

CC 4.0 BY-SA版权

文章标签：

#学习 #算法 #c++

本文介绍了三种计算组合数的方法：暴力实现、预处理阶乘及逆元、以及利用卢卡斯定理优化。重点讲解了卢卡斯定理在大数情况下的高效计算，并提供了相应的代码实现和证明过程。

一、通过定义实现

$Cnm=n!(n−m)!C_n^m = \frac{n!}{(n-m)!}$

1、暴力实现
适用于 $n$ ， $m$ 非常小的时候，不然会爆longlong
优势：代码量小
可以加取模

时间复杂度 $O (m)$

#define ll long long
ll C(int a, int b)
{
	if (a < 0 || b < 0 || a < b) return 0;
    ll res = 1;
    for (int i = 1 , j = a; i <= b; i ++ , j--) 
    	res = res * j / i ;
    return res;
}

2、预处理阶乘及逆元实现

$O (n) 预处理， O (1) 查询$

const int mod = 1e9 + 7;
const int N = 1e6 + 7;
ll jc[N],fjc[N];
ll ksm(ll a, ll b)
{
    ll res = 1;
    while (b)
    {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}
ll C(int a, int b)
{
	if (a < 0 || b < 0 || b > a) return 0;
    if (b == 0 || a == b) return 1;
    return jc[a] * fjc[a - b] % mod * fjc[b] % mod;
}
void init()
{
    jc[0] = 1;
    for (int i = 1; i < N; i++) 
    	jc[i] = jc[i - 1] * i % mod;
    fjc[N - 1] = ksm(jc[N - 1], mod - 2);
    for (int i = N - 2; i; i--)  
    	fjc[i] = fjc[i + 1] * (i + 1) % mod;
}

二、通过递推式实现

$C_n^m=C_{n-1}^m+C_{n-1}^{m-1}$
$C_i^0=C_i^i=1$

$O (n m) 预处理， O (1) 查询$

const int N = 3000
int C[N][N];
void init()
{
	C[0][0] = 1;
	for(int i = 1 ; i < N ; i ++)
	{
		C[i][0] = C[0][i] = 1;
		for(int j = 1 ; j < i  ; j ++)
		{
			C[i][j] = (C[i-1][j] + C[i-1][j-1]) % mod;
		}
	}
}

三、卢卡斯定理

定理内容： $p)C_n^m \equiv C_{n\%p}^{m\%p} * C_{n/p}^{m/p}(mod\ p)$

适用范围： $n$ ， $m$ 很大，但模数不大且是素数时

时间复杂度 $O(log_p^m+p)$
注意更换模数之后要重新初始化

const int N = 1e6 + 10;
int mod;
ll jc[N],fjc[N];
ll ksm(ll a, ll b)
{
    ll res = 1;
    while (b)
    {
        if (b & 1) res = res * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return res;
}
ll C(int a, int b)
{
	if (a < 0 || b < 0 || a < b) return 0;
    if (b == 0 || a == b) return 1;
    return jc[a] * fjc[a - b] % mod * fjc[b] % mod;
}
void init()	//使用预处理阶乘方法求组合数，注意上界为mod而不是之前的N
{
    jc[0] = 1;
    for (int i = 1; i < mod; i++)
    	jc[i] = jc[i - 1] * i % mod;
    fjc[mod - 1] = ksm(jc[mod - 1], mod - 2);
    for (int i = mod - 2; i; i--)
    	fjc[i] = fjc[i + 1] * (i + 1) % mod;
}

int LC(ll a , ll b)
{
    if(max(a,b) < mod) return C(a,b);
    return C(a%mod,b%mod) * LC(a/mod , b/mod) % mod;
}

证明：
$p)C_n^m \equiv C_{n\%p}^{m\%p} * C_{n/p}^{m/p}(mod\ p)$

下面过程中将省略 $p)(mod\ p)$ ，注意同余号 $“≡”“\equiv”$ 和等号 $“ = ”$ 的区别。

$1+x)^p$
$\sum _{i=0}^{p}C_p^ix^i$
$1+x^p+\sum _{i=1}^{p-1}C_p^ix^i$

注意 $Cpi=p!i!(p−i)!=(p−1)!∗pi!(p−i)!C_p^i = \frac{p!}{i!(p-i)!} = \frac{(p-1)!*p}{i!(p-i)!}$ ，
由于 $p$ 是素数所以分子中的 $p$ 是不会被约分掉的，
所以 $0 < i < p$ 时，有 $C_p^i\%p = 0$

所以得到 $(1+x)p≡1+xp(1+x)^p \equiv 1+x^p$

将 $n$ ， $m$ 转化为 $p$ 进制数
即设： $n=∑i=0kni∗pin=\sum_{i=0}^kn_i*p^i$ ， $m=∑i=0kmi∗pim=\sum_{i=0}^km_i*p^i$

$(1+x)n=∏i=0k(1+x)ni∗pi=∏i=0k((1+x)pi)ni(1+x)^n=\prod_{i=0}^{k}(1+x)^{n_i*p_i}=\prod_{i=0}^{k}((1+x)^{p_i})^{n_i}$

$(1+x)n≡∏i=0k(1+xpi)ni(1+x)^n\equiv \prod_{i=0}^{k}(1+x^{p_i})^{n_i}$

左式中 $x^m$ 的系数为 $C_n^m$

右式中 $x^{p_im_i}$ 的系数为 $C_{n_i}^{m_i}$
而 $∏i=0kxpimi=xm\prod_{i=0}^{k}x^{p_im_i} = x^m$

所以 $Cnm≡∏i=0kCnimiC_{n}^{m} \equiv \prod_{i=0}^{k}C_{n_i}^{m_i}$ ，这个式子也就是卢卡斯定理的非递归形式。

继续推导：
$∏i=0kCnimi\prod_{i=0}^{k}C_{n_i}^{m_i}$
$=Cn0m0∏i=1kCnimi=C_{n_0}^{m_0}\prod_{i=1}^{k}C_{n_i}^{m_i}$
根据进制转换的性质可以得到
上式 $C_{n\%p}^{m\%p} * C_{n/p}^{m/p}$
证毕。