HDU 3944

来源：互联网发布：线切割锥度编程编辑：程序博客网时间：2024/05/01 20:01

原来剑哥讲组合数取模的时候说过Lucas定理，但是不是很明白。总之C(n,k)%p = C(n/p,k/p)*C(n%p,k%p)%p。这里要注意如果k>n, C(n,k)=0。C(n,k)%p可以通过求逆元的方法计算。不过C(n,k)%p = n!/(k!*(n-k)!)%p = n!*((k!*(n-k)!)^(p-2))%p是更好的方法，编码简单一点，时间效率高一点。

还要注意一点，这题比较恶心，每次都初始化阶乘数组肯定会超时。考虑到10000以内只有1229个素数，所以可以先打表。

#include <stdio.h>#include <stdlib.h>#include <math.h>#include <string.h>#include <memory.h>#include <string>#include <vector>#include <map>#include <set>#include <queue>#include <algorithm>#include <iostream>#include <sstream>using namespace std;// Lucas定理求解C(n,k)%p, p是素数// 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)!)^(p-2))%p// pow_mod(a,b,p)求a^b%p的值, 快速幂取模算法// com_mod(n,k,p)求C(n,k)%p的值// lucas(n,k,p)使用lucas定理求C(n,k)%p的值int pow_mod(int a, int b, int p) {    int r = 1, t = a;    while (b) {        if (b & 1) r = r * t % p;        t = t * t % p;        b >>= 1;    }    return r;}int com_mod(int n, int k, int p, int *fac) {    int a, b;    if (k > n) return 0;    a = fac[n];    b = fac[n-k] * fac[k] % p;    return a * pow_mod(b, p-2, p) % p;}int lucas(int n, int k, int p, int *fac) {    int r = 1, a, b;    while (n && k && r) {        a = n % p;        b = k % p;        n /= p;        k /= p;                r = r * com_mod(a, b, p, fac) % p;    }    return r;}int pri[1250], cnt;int fac[1250][10000];bool not_prime[10005];void init() {    int i, j, p;    cnt = 0;    for (i = 2; i < 10000; i++) {        if (not_prime[i]) continue;        pri[cnt++] = i;        for (j = i + i; j < 10000; j += i)            not_prime[j] = 1;    }    for (i = 0; i < cnt; i++) {        fac[i][0] = 1;        p = pri[i];        for (j = 1; j < p; j++)            fac[i][j] = j * fac[i][j-1] % p;    }}int *find(int p) {    int l = 0, r = cnt - 1;    while (l <= r) {        int m = (l + r) / 2;        if (pri[m] == p) return fac[m];        else if (pri[m] < p) l = m + 1;        else r = m - 1;    }    return fac[0];}int main() {    int n, k, p, t = 1;    init();    while (scanf("%d %d %d", &n, &k, &p) != EOF) {        if (k > n - k) k = n - k;        printf("Case #%d: ", t++);        printf("%d\n", (lucas(n+1, k, p, find(p)) + n - k) % p);    }    return 0;}