組合數取模問題為求$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 }
