题目
题目链接:http://acm.hdu.edu.cn/showproblem.php?pid=5446
题目来源:2015 长春网络赛 模板题 不过吃透了对数论理解很有帮助
简要题意: M=∏ki=1pi 求 CmnmodM
数据范围: 1⩽m⩽n⩽1018;1⩽k⩽10;prime pi⩽105
题解
对于每个 pi 用Lucas去求出 Cmnmodpi 然后用中国剩余定理合并。
题目的代码量较大,由于是很经典的题目,所以最后变成了模板题。
求出 Cmnmodpi 然后去用中国剩余定理合并得到 Ans ,可以保证 Cmn≡AnsmodM 。
而显然结果的解集为 {x∣x=Ans+kM∣k∈Z}
可以保证在 [0,M) 中只有一个数,这个数就一定是结果了。
实现
灵活运用Lucas,组合数取模,求逆,中国剩余定理就能搞出来。
每个质数可以预处理阶乘表来提速,不打表for(1…m)大概是800ms,打表之后的可以到15ms。
代码
#include <iostream>
#include <cstdio>
#include <cmath>
#include <algorithm>
#include <cstring>
#include <stack>
#include <queue>
#include <string>
#include <vector>
#include <set>
#include <map>
#define pb push_back
#define mp make_pair
#define all(x) (x).begin(),(x).end()
#define sz(x) ((int)(x).size())
#define fi first
#define se second
using namespace std;
typedef long long LL;
typedef vector<int> VI;
typedef pair<int,int> PII;
LL powmod(LL a,LL b, LL MOD) {LL res=1;a%=MOD;for(;b;b>>=1){if(b&1)res=res*a%MOD;a=a*a%MOD;}return res;}
// head
LL MOD;
const int N = 1E5+5;
LL comb[15];
LL p[15];
LL M[15];
LL inv[15];
LL fac[N];
inline LL mulmod(LL a, LL b, LL p) {
LL res = 0;
while (b) {
if (b & 1) res = (res+a)%p;
b >>= 1;
a = (a+a)%p;
}
return res;
}
LL getInv(LL a, LL b) {
return powmod(a, b-2, b);
}
LL C(LL n, LL m, LL p) {
if(m > n) return 0;
return fac[n]*getInv(fac[n-m], p)*getInv(fac[m], p)%p;
}
LL Lucas(LL n, LL m, LL p) {
if (m == 0) return 1;
return C(n%p, m%p, p) * Lucas(n/p, m/p, p) % p;
}
LL CRT(LL *mod, LL *inv, LL *M, LL len, LL sum) {
for (int i = 0; i < len; i++) {
M[i] = sum/mod[i];
inv[i] = getInv(M[i], mod[i]);
}
return sum;
}
LL CRTcal(LL *inv, LL *M, LL len, LL *r, LL sum) {
LL ans = 0;
for (int i = 0; i < len; i++) {
LL temp = mulmod(inv[i], M[i], sum);
temp = mulmod(temp, r[i], sum);
ans = (ans + temp)%sum;
}
return ans;
}
int main() {
LL m, n, k, t;
scanf("%I64d", &t);
while (t--) {
scanf("%I64d%I64d%I64d", &n, &m, &k);
fac[0] = MOD = 1;
for (int i = 0; i < k; i++) {
scanf("%I64d", p+i);
for (int j = 1; j < p[i]; j++) {
fac[j] = fac[j-1]*j%p[i];
}
comb[i] = Lucas(n, m, p[i]);
MOD *= p[i];
}
CRT(p, inv, M, k, MOD);
printf("%I64d\n", CRTcal(inv, M, k, comb, MOD));
}
return 0;
}