组合数取模
组合数取模问题为求$C_{n}^m % p$的值。根据$n$,$m$,$p$取值不同,方法不同。
在此之前我们先看些前置技能:
同余定理:$a≡b(mod\ m)$
性质:
1.传递性:若$a≡b(mod\ m)$,$b≡c(mod\ m)$,则$a≡c(mod\ m)$;
2.同余式相加:若$a≡b(mod\ m)$,$c≡d(mod\ m)$,则$a±c≡b±d(mod\ m)$;
3.同余式相乘:若$a≡b(mod\ m)$,$c≡d(mod\ m)$,则$ac≡bd(mod\ m)$。
逆元(数论倒数):对于正整数$a$和$m$,如果有$ax≡1(mod\ m)$,那么把这个同余方程中$x$的最小正整数解称为$a mod\ m$的逆元。
为什么叫数论倒数呢,因为没有$mod$操作,那么$x$就相当于$a$的倒数,有$mod$操作时效果上和倒数一样。同时我们把$a$的逆元写做:$inv(a)$。
逆元求解一般用扩展欧几里得算法,如果$m$为素数,那么还可以根据费马小定理得到逆元为$a^{m-2} mod\ m$。
扩展欧几里得算法和费马小定理求解逆元具有局限性,两种算法都要求a和m互素。
1.扩展欧几里得:$a*x + b*y = 1$
解$x$就是$a$关于$b$的逆元,解$y$就是$b$关于$a$的逆元
$a*x \% b + b*y \% b = 1 \% b$
$a*x \% b = 1 \% b$
$a*x = 1 (mod\ b)$
2.费马小定理:$a^{p-1}≡1 (mod\ p)$
两边同时除以$a$,$a^{p-2}≡inv(a)(mod\ p)$;即$inv(a)≡a^{p-2}(mod\ p)$
Lucas定理:
$n=n_kp^k+n_{k-1}p^{k-1}+...+n_1p+n_0$
$m=m_kp^k+m_{k-1}p^{k-1}+...+m_1p+m_0$
得到:$C_n^m=\prod\limits_{i=0}^{k}C_{ni}^{mi}(mod\ p)$
具体代码中实现:$Lucas(n,m,p)=C(n\%p,m\%p)*Lucas(n/p,m/p,p)\%p $
经典例题1:FZU 2020
$1 <= m <= n <= 10^9, m <= 10^4, m < p < 10^9, p$是素数,属于$m$比较小,$n$,$p$比较大的情况。
该题cin输入会比较快,怀疑是scanf读入%lld比较耗时,或者是OJ问题。
解决方案1:Lucas定理
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 using namespace std; 5 6 typedef long long LL; 7 8 LL fast_power(LL a,LL b,LL p){ 9 LL ans=1; 10 a%=p; 11 while(b){ 12 if(b&1) ans=(ans*a)%p; 13 a=(a*a)%p; 14 b>>=1; 15 } 16 return ans; 17 } 18 19 LL C(LL a,LL b,LL p){ 20 if(b>a) return 0; 21 if(a==b) return 1; 22 LL ans1=1,ans2=1; 23 for(LL i=1;i<=b;i++){ 24 ans1=ans1*(a-i+1)%p; 25 ans2=ans2*i%p; 26 } 27 return ans1*fast_power(ans2,p-2,p)%p; 28 } 29 30 LL Lucas(LL a,LL b,LL p){ 31 if(b==0) return 1; 32 return C(a%p,b%p,p)*Lucas(a/p,b/p,p)%p; 33 } 34 35 int main(){ 36 int t; 37 cin>>t; 38 while(t--){ 39 LL n,m,p; 40 cin>>n>>m>>p; 41 cout<<Lucas(n,m,p)<<endl; 42 } 43 return 0; 44 }
解决方案2:暴力+逆元
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 using namespace std; 5 6 typedef long long LL; 7 8 LL fast_power(LL x,LL n,LL mod){ 9 LL ans=1; 10 x%=mod; 11 while(n){ 12 if(n&1) ans=(ans*x)%mod; 13 n>>=1; 14 x=(x*x)%mod; 15 } 16 return ans; 17 } 18 19 int main(){ 20 int t; 21 cin>>t; 22 while(t--){ 23 LL n,m,p; 24 cin>>n>>m>>p; 25 LL ans1=n,ans2=1; 26 for(LL i=1;i<m;i++){ 27 ans1=ans1*(n-i)%p; 28 ans2=ans2*i%p; 29 } 30 ans2=ans2*m%p; 31 cout<<ans1*fast_power(ans2,p-2,p)%p<<endl; 32 } 33 return 0; 34 }
经典例题2:Codeforces 101775A
$1 ≤ N ≤ 10^9,3 ≤ K ≤ 10^5,mod=1000000007$ ,和上题一样都属于$m$比较小,$n$,$p$比较大的情况。
解决方案:通过二项式定理的性质,得到所有情况总和为$2^n-1$,题目要求$C_n^k+C_n^{k+1}+...C_n^n$的值, 我们可以通过用总和减去$C_n^1+C_n^2+...+C_n^{k-1}$,得到最终答案。直接暴力+逆元,在for的过程中,我们用总和减去对应的$C$,注意$+mod$再去减,因为取余,值可能为负。
1 #include <bits/stdc++.h> 2 using namespace std; 3 4 typedef long long LL; 5 const LL mod=1000000007; 6 7 LL fast_power(LL x,LL n){ 8 LL ans=1; 9 x%=mod; 10 while(n){ 11 if(n&1) ans=(ans*x)%mod; 12 n>>=1; 13 x=(x*x)%mod; 14 } 15 return ans; 16 } 17 18 int main(){ 19 int t; 20 scanf("%d",&t); 21 for(int Case=1;Case<=t;Case++){ 22 LL n,k; 23 scanf("%lld%lld",&n,&k); 24 LL res=fast_power(2,n)-1,p=n; 25 for(LL i=1;i<k;i++){ 26 res=(res-p+mod)%mod; 27 p=(p*(n-i))%mod; 28 p=(p*fast_power(i+1,mod-2))%mod; 29 } 30 printf("Case #%d: %lld\n",Case,res); 31 } 32 return 0; 33 }
经典例题3:HDU 4349
题意:求$C_n^0,C_n^1,C_n^2...C_n^n$中有多少个奇数。
官方题解:
本题为Lucas定理推导题,我们分析一下 $C(n,m)\%2$,那么由Lucas定理,我们可以写成二进制的形式观察,比如 $n=1001101$,$m$是从000000到1001101的枚举,我们知道在该定理中$C(0,1)=0$,因此如果$n=1001101$的0对应位置的$m$二进制位为1那么$C(n,m) \% 2==0$,因此$m$对应$n$为0的位置只能填0,而1的位置填0,填1都是1(C(1,0)=C(1,1)=1),不影响结果为奇数,并且保证不会出$n$的范围,因此所有的情况即是$n$中1位置对应$m$位置0,1的枚举,那么结果很明显就是:2^(n中1的个数)。
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 using namespace std; 5 6 int main(){ 7 int n; 8 while(scanf("%d",&n)!=EOF){ 9 int cnt=0; 10 while(n){ 11 if(n&1) cnt++; 12 n>>=1; 13 } 14 printf("%d\n",1<<cnt); 15 } 16 return 0; 17 }
经典例题4:HDU 3944
题意:从$(0,0)$到$(n,k)$位置总和最小的走法,并对$p$取模。$0<=k<=n<10^9,p<10^4$,属于$k$和$n$较大,$p$较小的情况。
解决方案:显然要么从0,0开始先往下走,再往右下角走;要么先往右下角走,再往下走。因为第一列和最右边那一斜列都为1,思路肯定是尽量走多的1,具体比较$n$和$k/2$的大小。
方式1:$n-k+C_{n-k}^0+C_{n-k+1}^1+...C_n^k$, 通过加个为0的$C_{n-k}^{-1}$,推出前式为$n-k+C_{n+1}^k$。
方式2:$k+C_k^k+C_{k+1}^k+C_{k+2}^k+C_{k+3}^k+...C_n^k$,通过加个值为0的$C_k^{k+1}$,推出前式为$k+C_{n+1}^{k+1}$。
Lucas定理计算,但是有100000 组数据,所以我们通过$p$较小进行计算,先处理1e4内素数的阶乘和阶乘的逆元,因为用Lucas的时候,最后肯定都会变成素数的情况进行计算,取模了嘛。
1 #include <cstdio> 2 #include <iostream> 3 #include <algorithm> 4 using namespace std; 5 6 int cnt=0; 7 const int N=1e4+10; 8 typedef long long LL; 9 10 LL prime[N]; 11 LL fac[N][N],inv[N][N]; 12 bool is_prime[N+1]; 13 14 LL fast_power(LL a,LL b,LL p){ 15 LL ans=1; 16 a%=p; 17 while(b){ 18 if(b&1) ans=(ans*a)%p; 19 a=(a*a)%p; 20 b>>=1; 21 } 22 return ans; 23 } 24 25 void init(){ 26 for(int i=0;i<N;i++) is_prime[i]=1; 27 is_prime[0]=is_prime[1]=0; 28 for(int i=2;i<N;i++){ 29 if(is_prime[i]){ 30 prime[cnt++]=i; 31 for(int j=2*i;j<N;j+=i) is_prime[j]=0; 32 } 33 } 34 35 for(int i=0;i<cnt;i++){ 36 fac[prime[i]][0]=inv[prime[i]][0]=1; 37 for(int j=1;j<prime[i];j++){ 38 fac[prime[i]][j]=(fac[prime[i]][j-1]*j)%prime[i]; 39 inv[prime[i]][j]=fast_power(fac[prime[i]][j],prime[i]-2,prime[i]); 40 } 41 } 42 } 43 44 LL C(LL a,LL b,LL p){ 45 if(b>a) return 0; 46 if(a==b) return 1; 47 return fac[p][a]*(inv[p][b]*inv[p][a-b]%p)%p; 48 } 49 50 LL Lucas(LL a,LL b,LL p){ 51 if(b==0) return 1; 52 return C(a%p,b%p,p)*Lucas(a/p,b/p,p)%p; 53 } 54 55 int main(){ 56 int Case=1; 57 LL n,k,p; 58 init(); 59 while(scanf("%lld%lld%lld",&n,&k,&p)!=EOF){ 60 if(2*k<=n) k=n-k; 61 printf("Case #%d: %lld\n",Case,(Lucas(n+1,k+1,p)+k)%p); 62 Case++; 63 } 64 return 0; 65 }
经典例题5:ZOJ 4536
题意:从$n$个数$(1-n)$中选出$m$个两两不相邻的数,求有多少种方案。
解决方案:需要选定$m$个数,也就是说剩余$n-m$个数,$n-m+1$个空位,从这些空位中选取$m$个位置,所以最终答案为$C_{n-m+1}^m$。
1 #include <cstdio> 2 #include <cstring> 3 #include <iostream> 4 #include <algorithm> 5 using namespace std; 6 7 typedef long long LL; 8 9 LL fast_power(LL x,LL n,LL mod){ 10 LL ans=1; 11 x%=mod; 12 while(n){ 13 if(n&1) ans=(ans*x)%mod; 14 n>>=1; 15 x=(x*x)%mod; 16 } 17 return ans; 18 } 19 20 LL C(LL a,LL b,LL p){ 21 if(b>a) return 0; 22 if(a==b) return 1; 23 LL ans1=1,ans2=1; 24 for(LL i=1;i<=b;i++){ 25 ans1=ans1*(a-i+1)%p; 26 ans2=ans2*i%p; 27 } 28 return ans1*fast_power(ans2,p-2,p)%p; 29 } 30 31 LL Lucas(LL a,LL b,LL p){ 32 if(b==0) return 1; 33 return C(a%p,b%p,p)*Lucas(a/p,b/p,p)%p; 34 } 35 36 int main(){ 37 LL n,m,p; 38 while(scanf("%lld%lld%lld",&n,&m,&p)!=EOF){ 39 printf("%lld\n",(Lucas(n-m+1,m,p)%p)); 40 } 41 return 0; 42 }