各種逆元求法 組合數取模 comb (組合數 Lucas)
組合數取模(comb)
【問題描述】
計算C(m,n)mod 9901的值
【輸入格式】
從檔案comb.in中輸入資料。
輸入的第一行包含兩個整數,m和n
【輸出格式】
輸出到檔案comb.out中。
輸出一行,一個整數
【樣例輸入】
2 1
【樣例輸出】
2
【資料規模與約定】
對於 20%的資料,n<=m<=20
對於 40%的資料,n<=m<=2000
對於 100%的資料,n<=m<=20000
題解:組合數取模(comb)
對於 20%資料:m<=20,直接計算最後取模。
對於 40%資料:m<=2000,預處理 1~2000 的逆元,暴力計算。
對於 100%資料:m<=20000,使用 lucas 定理把縮小至 9901 以內,再用上一個方
法計算。
求組合數直接算階乘肯定是不行的,會爆。(不能邊算邊mod,做除法是會出錯)遞推太慢了,捨棄。那麼怎麼辦呢?只好又回去考慮,既然一定要取模,那麼自然就想到一些在取模意義上相等的演算法,要解決除法問題就剩下逆元(在模意義下將除法轉為乘法)了。
逆元求法較多,下面提供幾種:
(取模意義下除法 ->( a 除以 b )mod 一個數)
1.擴充套件歐幾里得
inline long long extend_gcd(long long a,long long b,long long &x,long long &y){
if(a == 0 && b == 0 )
return -1;
if(b == 0){
x = 1; y = 0;
return a;
}
long long d = extend_gcd(b, a % b, y, x);
y -= a / b * x;
return d;
}
inline long long mod_reverse(long long a,long long n){
long long d = extend_gcd(a,n,x,y);
if(d == 1)
return ( x % n + n ) % n;
else
return -1;
}
int main(){
cin >> a >> b;
long long nn = mod_reverse(b, mod);//逆元
cout << a * nn % mod<< endl;
return 0;
}
2.費馬小定理(快速冪)
long long power_mod(long long z,long long y){
long long ans = 1; z %= MOD;
while( y ){
if(y & 1) ans = (ans * z) % MOD;
z = (z * z) % MOD;
y >>= 1;
}
return ans;
}
int main(){
cin >> t;
while( t-- ){
cin >> a >> b;
cout << (a * power_mod(b, MOD - 2)) % MOD <<endl;
}
return 0;
}
3.遞推求階乘逆元
a = mod / x
b = mod % x
(a*x + b) % mod = 0
a*x % mod = -b % mod
-a % mod = b * inv(x) % mod
(mod - a) % mod = b * inv(x) % mod
(mod - a) % mod * inv(b) = inv(x)
long long inva[1000005];//1--1000005的逆元
long long MOD = 1e9 + 7, ans;
int main(){
ans = 1;
inva[1] = 1;
for(int i=2; i<=1000000; i++){
inva[i] = inva[ MOD % i ] * ( -MOD / i );
inva[i] = ( inva[i] % MOD + MOD ) % MOD;
ans *= inva[i];
ans %= MOD;
}
cout << ans;
return 0;
}
這樣一來正常的資料就可以安安穩穩地AC了。
不過這道題,還要在優化一下–>Lucas
顯然有了這個公式,我們就可以把組合數限定在MOD以內。
#include <cstdio>
#include <iostream>
#include <algorithm>
#include <cmath>
#define LL long long
using namespace std;
const int mod = 9901;
LL mub[10010];
LL x, y ;
void init(){//預處理階乘
mub[0] = 1;//注意0! = 1
for(int i=1; i<=10000; i++){
mub[i] = mub[i-1] * i % mod;
}
}
LL work( LL m , LL n ) {//遞迴方法,要爆棧已捨棄
if(m == n) return 1;
if(n == 1) return m;
else return (work( m-1, n-1 ) + work( m-1, n )) % mod;
}
inline LL extend_gcd(LL a, LL b, LL &x, LL &y){//擴充套件歐幾里得求逆元x
if(a == 0 && b == 0)
return -1;
if(b == 0){
x = 1; y = 0;
return a;
}
LL d = extend_gcd(b, a % b, y, x);
y -= a / b * x;
return d;
}
inline LL mod_reverse(LL a, LL n){//規範逆元X
LL d = extend_gcd(a, n, x, y);
if(d == 1)
return ( x % n + n ) % n;
else
return -1;
}
LL solve(LL a, LL b){
if(a > b) return 0;//
LL nn = mod_reverse((mub[a] * mub[(b + mod - a) % mod]) % mod, mod);
return mub[b] * nn % mod;
}
void to_solve(LL a, LL b){//Lucas降資料
if(b < mod){
solve(a, b);
return;
}
cout << solve(a/mod, b/mod) * solve(a%mod, b%mod) << endl;
}
int main(){
freopen("comb.in","r",stdin);
freopen("comb.out","w",stdout);
LL a, b;
cin >> b >> a;
init();
if(b < mod){
printf("%lld", solve(a, b));
return 0;
}
else{
to_solve(a, b);
return 0;
}
}