Lucas 适合大组合数取模。
Lucas定理是用来求 c(n,m) mod p,p为素数的值。
时间复杂度O(logp(n)*p):
表达式C(n,m)%p=C(n/p,m/p)*C(n%p,m%p)
重新证明一遍:
令n=sp+q , m=tp+r .(q ,r ≤p)
原式子=C(sp+q,tp+r)%p
(1+x)^nΞ(1+x)^(sp+q)Ξ(1+x)^(sp)*(1+x)^(q)(mod p)
对于(1+x)^(sp) ,p为素数,二项式展开易得
((1+x)^p)^s = (1^p+x^p)^s = (1+x^p)^s
原式Ξ(1+x^p)^s*(1+x)^q(mod p)
(modp)
由ip+j=tp+r 。解得i = t , j = r 。
故当且仅当i = t , j =r 时有,x^(tp+r)
二项式系数为C(sp+q,tp+r)
即C(s,t)*x^(t*p)*C(q,r)*x^(r) (mod p)
得证
如此在求大组合数C(n,m),而p 又较小<1e5时十分便利。
不断取模取余的过程实际就是将n(10进制)写成p进制的形式。
实现时,涉及到分数模运算。求逆元;
分母^(p-1)Ξ1(mod P),或exgcd 分母即可。
51nod1120
#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long ll;
ll n;
const ll mod = 10007 ;
ll exgcd(ll a , ll b , ll &x , ll &y ){
if(!b) {
x = 1 ;
y = 0 ;
return a ;
}
else {
ll r = exgcd(b,a%b,x,y) ;
ll t = x ;
x = y ;
y = t - a/b*y ;
return r ;
}
}
ll cal_C(ll n , ll m){
ll ans = 1;
for(ll i = n ; i >= n-m+1 ; i--){
ans *= i ;
ans %= mod ;
}
return ans ;
}
ll lucas(ll n , ll m){
ll ans = 1 ;
while(n&&m){
ll s = n%mod , t = m%mod ;
if(s<t) return 0;
ll x ,y;
ll A = cal_C(s,t);
ll B = cal_C(t,t);
// cout << "B is " << B << endl ;
exgcd(B,mod,x,y);
// cout << "A is " << A << endl ;
// cout << "x is " << x << endl ;
if(x<0) x+=mod ;
n/=mod , m/=mod ;
// cout << " Bx + Mod*y = " << B*x+mod*y << endl ;
// cout << "Ax " << A*x%mod << endl ;
ans = ((ans*A%mod)*x)%mod;
}
// cout << "ans is " << ans << endl ;
return ans ;
}
//ll check(ll n , ll m ) {
// ll t = 1;
// for(ll i = n ; i >= n-m+1 ; --i){
// t *= i ;
t%= mod ;
// }
// ll k = 1 ;
// for(ll i = 1 ; i <= m ; ++i) {
// k*= i ;
k%= mod ;
// }
// cout << "t is " << t << endl << "k is " << k << endl ;
// return t/k/(m+1)%mod ;
//}
int main(){ while(~scanf("%lld",&n)){
n-- ;
// cout << "check is " << check(2*n,n) << endl ;
ll x ,y ;
exgcd(n+1,mod,x,y) ;
if(x<0)x+=mod ;
printf("%lld\n",lucas(2*n,n)*x*2%mod) ;
}
return 0 ;
}