HDU 3944 组合数学+数论

/********************************************************************************
    很好的组合数学题，跟数论结合，Lucas定理的应用，关键是这个转化：
	C(n, k) % p = C(n % p, k % p) * C(n / p, k / p) % p
	C(n, k) % p = n! / (k! * (n - k)!) % p = n! * (k! * (n - k)!)^-1 % p
	现在要求(k! * (n - k)!) % p的逆元，因为k!*(n - k)!不可能包含素因子p，所以逆元必然存在。
注意求得的逆元可能小于0！
	然后dfs下就完事，但是因为每次求阶乘显然会超时，所以事先打表处理～(PS:此段出处Racebug)
	切数论题还是要小心。。。整数域跟素数域都很容易坑爹的～
********************************************************************************/
#include <iostream>
#include <algorithm>
#include <cstdlib>
#include <cstring>
#include <utility>
#include <cstdio>
#include <memory>
#include <string>
#include <vector>
#include <cmath>
#include <ctime>
#include <queue>
#include <stack>
#include <map>
#include <set>
using namespace std;

//typedef long long LL;
typedef __int64 LL;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
typedef map<int, int>::iterator MI;
typedef vector<int>::iterator VI;
typedef set<int>::iterator SI;

const int INF_INT = 0x3f3f3f3f;
const double oo = 10e9;
const double eps = 10e-7;
const double PI = acos(-1.0);

const int MAXN = 10004;

int pcnt, prm[MAXN], dp[1500][MAXN];
bool is[MAXN];

inline bool scan_d(int &num)
{
	char in;
	in = getchar();
	if(EOF == in) 
	{
		return false;
	}
	while(in < '0' || in > '9')
	{
		in = getchar();
	}
	num = in - '0';
	while(in = getchar(), in >= '0' && in <= '9')
	{
		num *= 10, num += in - '0';
	}
	return true;
}
void preprocess()
{
	pcnt = 0;
	memset(is, true, sizeof(is));
	is[0] = is[1] = false;
	prm[pcnt++] = 2;
	for (int i = 4; i < MAXN; i += 2)
	{
		is[i] = false;
	}
	int ind;
	for (ind = 3; ind * ind <= MAXN; ind += 2)
	{
		if (is[ind])
		{
			prm[pcnt++] = ind;
			for (int j = ind * ind, k = ind * 2; j < MAXN; j += k)
			{
				is[j] = false;
			}
		}
	}
	while (ind < MAXN)
	{
		if (is[ind])
		{
			prm[pcnt++] = ind;
		}
		ind += 2;
	}
	for (int i = 0; i < pcnt; ++i)
	{
		dp[i][0] = 1;
		for (int j = 1; j < MAXN; ++j)
		{
			dp[i][j] = (dp[i][j - 1] * j) % prm[i];
		}
	}
	return ;
}
int ext_gcd(int a, int b, int &x, int &y)
{
	if (0 == b)
	{
		x = 1;
		y = 0;
		return a;
	}
	int d = ext_gcd(b, a % b, x, y);
	int t = x;
	x = y;
	y = t - a / b * y;
	return d;
}
int bisearch(int p)
{
	int low = 0, high = pcnt - 1, mid = (low + high) >> 1;
	int res = mid;
	while (low <= high)
	{
		mid = (low + high) >> 1;
		if (p == prm[mid])
		{
			res = mid;
			break ;
		}
		if (p < prm[mid])
		{
			high = mid - 1;
		}
		else
		{
			low = mid + 1;
		}
	}
	return res;
}
int Lucas_dfs(int n, int k, int id)
{
	int p = prm[id];
	int _n, _k, buf, tmp, _tmp, y;
	_n = n / p;
	_k = k / p;
	n %= p;
	k %= p;
	if (n < k || _n < _k)
	{
		return 0;
	}
	tmp = (dp[id][k] * dp[id][n - k]) % p;
	ext_gcd(tmp, p, _tmp, y);
	if (_tmp < 0)
	{
		_tmp += p;
	}
	buf = (dp[id][n] * _tmp) % p;
	if (0 == _n && 0 == _k)
	{
		return buf;
	}
	return buf * Lucas_dfs(_n, _k, id);
}
void ace()
{
	int cas = 0;
	int n, k, p;
	int id;
	while (scan_d(n))
	{
		scan_d(k);
		scan_d(p);
		if (k > n - k)
		{
			k = n - k;
		}
		id = bisearch(p);
		printf("Case #%d: %d\n", ++cas, (Lucas_dfs(n + 1, k, id) + n - k) % p);
	}
	return ;
}
int main()
{
	preprocess();
	ace();
	return 0;
}