数论 · exLucas 定理

前言

终于看懂了！！！连夜丢掉卡特兰来做了两道 exLucas 的题目。

定理

作为一个没什么用的铺垫：数论 · Lucas 定理

求解 $C_n^m\%p$（不保证 $p$ 为质数）。

证明

1

先把 $p$ 拆分一下：$p=p_1^{a_1}\ast p_2^{a_2}\ast \cdots p_k^{a_k}$。

然后逐一求解 $C_n^m\%p_i^{a_i}$，

最后用 CRT 求解 $\{\begin{array}{l}x\equiv b_1(\mathrm{m}\mathrm{o}\mathrm{d}{p}_1^{a_1})\\ x\equiv b_2(\mathrm{m}\mathrm{o}\mathrm{d}{p}_2^{a_2})\\ \cdots \\ x\equiv b_k(\mathrm{m}\mathrm{o}\mathrm{d}{p}_k^{a_k})\end{array}$ 的最小正整数解 $x$ 即可。

2

假设求解 $C_n^m\%p_i^k$，设 $p_k=p_i^k$。

有 $C_n^m=\frac{n!}{m!(n-m)!}$，

所以就变成求解 $\frac{n!\%p_k}{m!\%p_k(n-m)!\%p_k}$。

但是因为分母部分需要逆元，有可能 $m!\%p_k$ 或 $(n-m)!\%p_k$ 与 $p_k$ 不互质，所以我们要先除去他们中所有的 $p_i$，然后取模完之后再乘回去即可。

假设 $x$ 为 $n!$ 中包含的 $p_i$ 的个数，$y$、$z$ 同理。

则变为求解 $(\frac{\frac{n!}{p_i^x}\%p_k}{\frac{m!}{p_i^y}\%p_k\frac{(n-m)!}{p_i^z}\%p_k}\ast p_i^{x-y-z})\%p_k$

3

考虑求解 $\frac{n!}{p_i^x}\%p_k$。

先举个例子：

设 $n=19,p_i=3,k=2$。

先把 $n!$ 中 $p_i$ 的倍数先提取出来再合并一下：

$n=(1\ast 2\ast 4\ast 5\ast 7\ast 8\ast 10\ast 11\ast 13\ast 14\ast 16\ast 17\ast 19)\ast 3^6\ast (1\ast 2\ast 3\ast 4\ast 5\ast 6)$

按照 2 的后半部分所述，$3^6$ 可以先忽略掉。

而柿子的后半部分 $(1\ast 2\ast 3\ast 4\ast 5\ast 6)$，也就是 $\lfloor \frac{n}{p_i}\rfloor !$ 我们可以用递归再次求解。

对于柿子的前半部分，会发现它对 $p_i^k$ 取余会有一个周期：

$(1\ast 2\ast 4\ast 5\ast 7\ast 8)\equiv (10\ast 11\ast 13\ast 14\ast 16\ast 17)(\mathrm{m}\mathrm{o}\mathrm{d}{p}_i^k)$

所以只用求出第一个周期：$res={\prod }_{i=1}^{p_k}i(p_i\nmid i)$，然后再用快速幂求出 $re{s}^{\lfloor \frac{n}{p_k}\rfloor }$ 即可。

然后就会发先剩下个 $19$，这些数直接暴力乘进 $res$（进行完快速幂之后）即可。

4

关于 2 中 $x,y,z$ 的计算。

小奥中学过试除法（忘了叫啥），可以求解类似“100 的阶乘中有多少个因子 5”的问题。

这里直接除然后求出 $x-y-z$ 即可。

代码如下：

int al = n;
while(al)
	k += al / pi, al /= pi;
al = m;
while(al)
	k -= al / pi, al /= pi;
al = n - m;
while(al)
	k -= al / pi, al /= pi;

最后用 扩展欧几里得 求出 $\frac{m!}{p_i^y}\%p_k\frac{(n-m)!}{p_i^z}\%p_k$ 的逆元即可。

5

关于使用 CRT 合并。

每次求完 $C_n^m\%p_i^k=b_i$，就有了 $x\equiv b_i(\mathrm{m}\mathrm{o}\mathrm{d}{p}_i^{a_i})$。

先给代码：ans += bi * inv(p / pi, pi) % p * (p / pi) % p

bi * inv(p / pi, pi) % p 是为了让它在取模 $p_i$ 的情况下余数为 $b_i$；

* (p / pi) % p 是为了让它在取模其他质数的情况下余数为 0。

详见 CRT。

代码

#include<algorithm>
#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;

#define int long long
#define rep(i, a, b) for(int i = a; i <= b; ++i)

int n, m, p;

inline int power(int x, int t, int mod)
{
	int res = 1;
	while(t >= 1)
	{
		if(t & 1)
			res = res * x % mod;
		t /= 2, x = x * x % mod;
	}
	return res;
}

inline int fac(int nw, int pi, int pk)
{
	if(!nw) return 1;
	int res = 1;
	rep(i, 2, pk)
		if(i % pi) res = res * i % pk;
	res = power(res, nw / pk, pk);
	rep(i, 2, nw % pk)
		if(i % pi) res = res * i % pk;
	return res * fac(nw / pi, pi, pk) % pk;
}

int x, y;

inline void exGcd(int a, int b)
{
	if(!b)
	{
		x = 1, y = 0;
		return;
	}
	exGcd(b, a % b);
	int tmp = x;
	x = y, y = tmp - a / b * y;
	return;
}

inline int inv(int nw, int mod)
{
	exGcd(nw, mod);
	return (x + mod) > mod ? x : x + mod;
}

inline int C(int pi, int pk)
{
	int fn = fac(n, pi, pk), fm = fac(m, pi, pk), fnm = fac(n - m, pi, pk);
	int k = 0;
	int al = n;
	while(al)
		k += al / pi, al /= pi;
	al = m;
	while(al)
		k -= al / pi, al /= pi;
	al = n - m;
	while(al)
		k -= al / pi, al /= pi;
	return fn * inv(fm, pk) % pk * inv(fnm, pk) % pk * power(pi, k, pk) % pk;
}

inline int Crt(int bi, int mi)
{
	return bi * inv(p / mi, mi) % p * (p / mi) % p;
}

inline int exLucas()
{
	int pk, res = 0;
	int lim = ceil(sqrt(p)), tmp = p;
	rep(i, 2, lim)
	{
		if(tmp % i) continue;
		pk = 1;
		while(!(tmp % i)) 
			tmp /= i, pk *= i;
		res = (res + Crt(C(i, pk), pk)) % p;
	} 
	if(tmp > 1)
		res = (res + Crt(C(tmp, tmp), tmp)) % p;
	return res;
}

signed main()
{
	scanf("%lld %lld %lld", &n, &m, &p);
	printf("%lld\n", exLucas());
	return 0;
}

---

——$End$——