CVBB ACM Team

到此差别页面的链接

--- 2020-2021:teams:legal_string:jxm2001:数论_5 [2020/10/16 10:04]
jxm2001
+++ 2020-2021:teams:legal_string:jxm2001:数论_5 [2020/10/16 15:01] (当前版本)
jxm2001
@@ 行 5: / 行 5: @@
 ==== 算法简介 ====
-$O(p+\log_pn)$ 计算 ${n\choose m}\bmod p$ 的算法。
+$O(p+\log_pn)$ 计算 ${n\choose m}\bmod p$ 的算法，$p$ 为素数。
 ==== 算法实现 ====
+$${p\choose i}=\frac pi{p-1\choose i-1}$$
+于是对 $1\le i\le p-1$，有 ${p\choose i}\equiv 0\pmod p$。
+$$(1+x)^p\equiv {p\choose 0}+{p\choose 1}x+\cdots +{p\choose p}x^p\equiv 1+x^p\pmod p$$
+对 ${n\choose m}$，不妨设 $n=k_1p+r_1,m=k_2p+r_2$，于是有
+$$(1+x)^n\equiv (1+x)^{k_1p+r_1}\equiv ((1+x)^p)^{k_1}(1+x)^{r_1}\equiv (1+x^p)^{k_1}(1+x)^{r_1}\pmod p$$
+对比左右两边 $x_m$ 的系数，有
+$${n\choose m}\equiv {k_1\choose k_2}{r_1\choose r_2}\pmod p$$
+于是递归处理，有
+$$\text{Lucas}(n,m)\equiv \text{Lucas}(\lfloor \frac np\rfloor,\lfloor \frac mp\rfloor){n\bmod p\choose m\bmod p}$$
+==== 代码模板 ====
+[[https://www.luogu.com.cn/problem/P3807|洛谷p3807]]
+<hidden 查看代码>
+<code cpp>
+const int MAXN=1e5+5;
+int frac[MAXN],invfrac[MAXN];
+int quick_pow(int a,int b,int p){
+	int ans=1;
+	while(b){
+		if(b&1)ans=1LL*ans*a%p;
+		a=1LL*a*a%p;
+		b>>=1;
+	}
+	return ans;
+}
+int Lucas(int n,int m,int p){
+	if(n<m)return 0;
+	else if(n<p)return 1LL*frac[n]*invfrac[n-m]%p*invfrac[m]%p;
+	else return 1LL*Lucas(n%p,m%p,p)*Lucas(n/p,m/p,p)%p;
+}
+int main()
+{
+	int T=read_int();
+	while(T--){
+		int n=read_int(),m=read_int(),p=read_int();
+		frac[0]=1;
+		_for(i,1,p)frac[i]=1LL*frac[i-1]*i%p;
+		invfrac[p-1]=quick_pow(frac[p-1],p-2,p);
+		for(int i=p-1;i;i--)
+		invfrac[i-1]=1LL*invfrac[i]*i%p;
+		enter(Lucas(n+m,n,p));
+	}
+	return 0;
+}
+</code>
+</hidden>
+===== 拓展卢卡斯 =====
+==== 算法简介 ====
+$O\left(\sum p_i^k+\sum \log_2 n\log_{p_i^k}n\right)$ 计算 ${n\choose m}\bmod p$ 的算法，其中 $p=\prod p_i^k$。
+==== 算法实现 ====
+考虑计算 ${n\choose m}\equiv x_i\pmod {p_i^k}$，然后中国剩余定理合并答案。
+$$
+{n\choose m}\equiv \frac {n!}{m!(n-m)!}\pmod {p^k}
+$$
+分子不一定与模数互质，不能直接使用逆元计算。考虑提取出 $n!$ 所有的 $p$ 因子，记为 $g(n)$，同时记 $f(n)=\frac {n!}{p^{g(n)}}$。代入上式，有
+$$
+\cfrac{\frac {n!}{p^{g(n)}}}{\frac {m!}{p^{g(m)}}\frac {(n-m)!}{p^{g(n-m)}}}p^{g(n)-g(m)-g(n-m)}\equiv \cfrac{f(n)}{f(m)f(n-m)}p^{g(n)-g(m)-g(n-m)}\pmod {p^k}
+$$
+于是有 $(f(n),p^k)=1$，可以参与逆元计算。设 $n=ap^k+b$。
+$$
+n!=\prod_{i=1,p\mid i}^ni \prod_{i=1,p\not\mid i}^ni=(\lfloor\frac np\rfloor)!p^{\lfloor\frac np\rfloor}(\prod_{i=1,p\not\mid i}^{p^k}i)^a\prod_{i=1,p\not\mid i}^b i
+$$
+于是有
+$$
+f(n)=f(\lfloor\frac np\rfloor)(\prod_{i=1,p\not\mid i}^{p^k}i)^a\prod_{i=1,p\not\mid i}^b i
+$$
+$$
+g(n)=g(\lfloor\frac np\rfloor)+\lfloor\frac np\rfloor
+$$
+经过 $O\left(\sum p_i^k\right)$ 预处理后即可 $O\left(\sum \log_2 n\log_{p^k}n\right)$ 递归计算答案。
+==== 代码模板 ====
+[[https://www.luogu.com.cn/problem/P4720|洛谷p4720]]
+<hidden 查看代码>
+<code cpp>
+int quick_pow(int a,int b,int p){
+	int ans=1;
+	while(b){
+		if(b&1)ans=1LL*ans*a%p;
+		a=1LL*a*a%p;
+		b>>=1;
+	}
+	return ans;
+}
+void exgcd(LL a,LL b,LL &tx,LL &ty){
+	if(b==0){
+		tx=1,ty=0;
+		return;
+	}
+	exgcd(b,a%b,ty,tx);
+	ty-=a/b*tx;
+}
+int inv(int x,int mod){
+	LL t1,t2;
+	exgcd(x,mod,t1,t2);
+	return t1%mod;
+}
+namespace Lucas{
+	const int MAXM=30;
+	vector<int> frac[MAXM];
+	int p_cnt,prime[MAXM],y[MAXM],pod[MAXM];
+	void Init(int p,int mod){
+		frac[p_cnt].resize(mod+1);
+		frac[p_cnt][0]=1;
+		_rep(i,1,mod){
+			if(i%p)frac[p_cnt][i]=1LL*frac[p_cnt][i-1]*i%mod;
+			else frac[p_cnt][i]=frac[p_cnt][i-1];
+		}
+		prime[p_cnt]=p;
+		pod[p_cnt]=mod/p*(p-1);
+		y[p_cnt++]=mod;
+	}
+	void init(int mod){
+		p_cnt=0;
+		int t=mod;
+		for(int p=2;p*p<=t;p++){
+			if(t%p==0){
+				int mod2=1;
+				while(t%p==0)t/=p,mod2*=p;
+				Init(p,mod2);
+			}
+		}
+		if(t!=1)Init(t,t);
+	}
+	int f(LL n,int i){
+		int ans=1;
+		while(n){
+			ans=1LL*ans*quick_pow(frac[i][y[i]],(n/y[i])%pod[i],y[i])%y[i]*frac[i][n%y[i]]%y[i];
+			n/=prime[i];
+		}
+		return ans;
+	}
+	int g(LL n,int i){
+		int ans=0;
+		while(n)ans+=(n/=prime[i]);
+		return ans;
+	}
+	int cal(LL n,LL m,int i){
+		int t=1LL*f(n,i)*inv(1LL*f(m,i)*f(n-m,i)%y[i],y[i])%y[i];
+		t=1LL*t*quick_pow(prime[i],g(n,i)-g(m,i)-g(n-m,i),y[i])%y[i];
+		return t;
+	}
+	int C(LL n,LL m){
+		int mod=1,ans=0;
+		_for(i,0,p_cnt)mod*=y[i];
+		_for(i,0,p_cnt){
+			int x=cal(n,m,i);
+			ans=(ans+1LL*x*(mod/y[i])%mod*inv(mod/y[i],y[i]))%mod;
+		}
+		return (ans+mod)%mod;
+	}
+}
+int main()
+{
+	LL n=read_LL(),m=read_LL();
+	int mod=read_int();
+	Lucas::init(mod);
+	enter(Lucas::C(n,m));
+	return 0;
+}
+</code>
+</hidden>

CVBB ACM Team

用户工具

站点工具

页面工具