求C((m, n) % p) (p为素数)
表达式:C(n,m)%p=C(n/p,m/p)*C(n%p,m%p)%p
当p<1e6时
#include<algorithm>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long ll;
const int N = 100000;
ll f[N];
void init(int p) {
f[0] = 1;
for (int i=1; i<=p; ++i) f[i] = f[i-1] * i % p;
}
ll pow_mod(ll a, ll x, int p) {
ll ret = 1;
while (x) {
if (x & 1) ret = ret * a % p;
a = a * a % p;
x >>= 1;
}
return ret;
}
ll Lucas(ll n, ll k, int p) { //C (n, k) % p
ll ret = 1;
while (n && k) {
ll nn = n % p, kk = k % p;
if (nn < kk) return 0; //inv (f[kk]) = f[kk] ^ (p - 2) % p
ret = ret * f[nn] * pow_mod (f[kk] * f[nn-kk] % p, p - 2, p) % p;
n /= p, k /= p;
}
return ret;
}
int main(void) {
int p = 10007, n , m;
while(~scanf("%d %d", &n, &m)) {
init (p);
printf ("%lld\n", Lucas (n, m, p));
}
return 0;
}
扩展:当p较大(1e12左右)或者p不是素数时
#include<algorithm>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long ll;
ll POW(ll a, ll b, ll mod)
{
ll ans=1;
while(b)
{
if(b & 1) ans = ans * a % mod;
a = a *a % mod;
b >>= 1;
}
return ans;
}
ll POW(ll a, ll b)
{
ll ans = 1;
while(b)
{
if(b & 1) ans = ans * a;
a = a * a;
b >>= 1;
}
return ans;
}
ll exGcd(ll a, ll b, ll &x, ll &y)
{
ll t, d;
if(!b)
{
x=1;
y=0;
return a;
}
d = exGcd(b, a % b, x, y);
t = x;
x = y;
y = t - a / b * y;
return d;
}
bool modular(ll *a, ll *m, ll k)
{
ll d, t, c, x, y, i;
for(i = 2; i <= k; i++)
{
d = exGcd(m[1], m[i], x, y);
c = a[i]-a[1];
if(c % d) return false;
t = m[i] / d;
x = (c / d * x % t + t) % t;
a [1] =m[1] * x + a[1];
m [1] =m[1] * m[i] / d;
}
return true;
}
ll reverse(ll a, ll b)
{
ll x, y;
exGcd(a, b, x, y);
return (x % b + b) % b;
}
ll C(ll n, ll m, ll mod)
{
if(m > n) return 0;
ll ans = 1, i, a, b;
for(i = 1; i <= m; i++)
{
a = (n + 1 - i) % mod;
b = reverse(i % mod, mod);
ans = ans * a % mod * b % mod;
}
return ans;
}
ll C1(ll n, ll m, ll mod)
{
if(m == 0) return 1;
return C(n % mod, m % mod, mod) * C1(n / mod, m / mod, mod) % mod;
}
ll cal(ll n, ll p, ll t)
{
if(!n) return 1;
ll x = POW(p, t), i, y = n / x, temp = 1;
for(i = 1; i <= x; i++) if(i % p) temp = temp * i%x;
ll ans=POW(temp,y,x);
for(i=y*x+1;i<=n;i++) if(i%p) ans=ans*i%x;
return ans*cal(n/p,p,t)%x;
}
ll C2(ll n, ll m, ll p, ll t)
{
ll x=POW(p, t);
ll a, b, c, ap = 0, bp = 0, cp = 0, temp;
for(temp = n; temp ;temp /= p) ap += temp / p;
for(temp = m; temp; temp /= p) bp += temp / p;
for(temp = n - m; temp; temp /= p) cp += temp / p;
ap = ap - bp - cp;
ll ans = POW(p, ap, x);
a = cal(n, p, t);
b = cal(m, p, t);
c = cal(n - m, p, t);
ans = ans * a % x * reverse(b, x) % x * reverse(c, x) % x;
return ans;
}
//计算C(n,m)%mod
ll Lucas(ll n,ll m,ll mod)
{
ll i, t, cnt = 0;
ll A[205], M[205];
for(i = 2; i * i <= mod; i++)
if(mod % i == 0)
{
t = 0;
while(mod % i == 0)
{
t++;
mod /= i;
}
M[++cnt] = POW(i,t);
if(t==1) A[cnt] = C1(n, m, i);
else A[cnt] = C2(n, m, i, t);
}
if(mod > 1)
{
M[++cnt] = mod;
A[cnt] = C1(n, m, mod);
}
modular(A, M, cnt);
return A[1];
}
int main(void) {
ll mod = 10000006, n , m;
while(~scanf("%lld %lld", &n, &m)) {
printf ("%lld\n", Lucas (n, m, mod));
}
return 0;
}