BZOJ2142 禮物 【擴展Lucas】
題目
一年一度的聖誕節快要來到了。每年的聖誕節小E都會收到許多禮物，當然他也會送出許多禮物。不同的人物在小E
心目中的重要性不同，在小E心中分量越重的人，收到的禮物會越多。小E從商店中購買了n件禮物，打算送給m個人
，其中送給第i個人禮物數量為wi。請你幫忙計算出送禮物的方案數（兩個方案被認為是不同的，當且僅當存在某
個人在這兩種方案中收到的禮物不同）。由於方案數可能會很大，你只需要輸出模P後的結果。
輸入格式
輸入的第一行包含一個正整數P，表示模；
第二行包含兩個整整數n和m，分別表示小E從商店購買的禮物數和接受禮物的人數；
以下m行每行僅包含一個正整數wi，表示小E要送給第i個人的禮物數量。
輸出格式
若不存在可行方案，則輸出“Impossible”，否則輸出一個整數，表示模P後的方案數。
輸入樣例
100
4 2
1
2
輸出樣例
12
提示
【樣例說明】
下面是對樣例1的說明。
以“/”分割，“/”前後分別表示送給第一個人和第二個人的禮物編號。12種方案詳情如下：
1/23 1/24 1/34
2/13 2/14 2/34
3/12 3/14 3/24
4/12 4/13 4/23
【數據規模和約定】
設P=p1^c1 * p2^c2 * p3^c3 * … *pt ^ ct，pi為質數。
對於100%的數據，1≤n≤109，1≤m≤5，1≤pi^ci≤10^5。
題解
式子很簡單，記\(sum[i]\)為w[i]前綴和：
\[ans = {n \choose sum[m]} \prod\limits {sum[m] - sum[i - 1] \choose w[i]}\]
重點在於計算\(C_{n}^{m} \mod P\)，其中\(n,m \le 10^9\)且\(P = p_1^{k_1}*p_2^{k_2}*p_3^{k_3}.....\)，其中每一個\(p_i^{k_i} \le 10^5\)
擴展Lucas
對於質數\(P \le 10^5\)，我們可以用Lucas定理計算出
\[C_{n}^{m} = \prod\limits_{i = 1} C_{n \mod P^i}^{m \mod P^i}\]
但對於合數\(P\)，Lucas定理就不再適用了
於是我們使用中國剩余定理
\[\left\{
\begin{array}{c}
x\equiv c_1\pmod {m_1}\\
x\equiv c_2\pmod {m_2} \\
x\equiv c_3\pmod {m_3}\...\x\equiv c_k\pmod {m_k}
\end{array}
\right.\]
顯然\(x\)就是答案，用中國剩余定理合並出\(x\)
我們只需要快速計算\(C_{n}^{m} \mod p_i^{k_i}\)
我們只需要快速計算\(n! \mod p_i^{k_i}\)
因為\((a + P) \equiv a \pmod P\)
而
\[n! = \prod\limits_{i = 1}^{n} i\]
所以我們對\(n\)個數按\(P\)進行分組並提取出其中\(p_i\)的倍數，假使有\(t\)個
可以得出
\[n! = p_i^{t} * (\prod\limits_{x = 1}^{p_i^{k_i}} x [x \mod p_i \ne 0])^{\frac{n}{p_i^{k_i}}} * (\lfloor \frac{n}{p_i} \rfloor)!\]
左邊是\(O(p_i^{k_i})\)的，右邊遞歸\(\lfloor \frac{n}{p_i} \rfloor\)
我們先提取出\(n,m,n - m\)的\(p_i\)，使其結果必定與\(p_i\)
其中\(n!\)中\(p\)的個數為\(x=\lfloor{n\over p}\rfloor+\lfloor{n\over p^2}\rfloor+\lfloor{n\over p^3}\rfloor+...\)
最後結合逆元計算出\(\frac{n!}{m!(n-m)!}\)再乘上\(p_i^{\sum t}\)就行了
#include<iostream>
#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#define LL long long int
#define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
#define REP(i,n) for (int i = 1; i <= (n); i++)
#define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<‘ ‘; puts("");
using namespace std;
const int maxn = 100005,maxm = 100005,INF = 1000000000;
inline int read(){
int out = 0,flag = 1; char c = getchar();
while (c < 48 || c > 57){if (c == ‘-‘) flag = -1; c = getchar();}
while (c >= 48 && c <= 57){out = (out << 3) + (out << 1) + c - 48; c = getchar();}
return out * flag;
}
LL sum;
int md,n,m,a[10];
int p[maxn],pk[maxn],pi,Ans[maxn];
void Sp(){
int x = md;
for (int i = 2; i * i <= x; i++){
if (x % i == 0){
p[++pi] = i; pk[pi] = 1;
while (x % i == 0) x /= i,pk[pi] *= i;
}
}
if (x - 1) ++pi,p[pi] = pk[pi] = x;
}
int qpow(int a,int b,int md){
int ans = 1;
for (; b; b >>= 1,a = 1ll * a * a % md)
if (b & 1) ans = 1ll * ans * a % md;
return ans;
}
void exgcd(int a,int b,int& d,int& x,int& y){
if (!b) {d = a; x = 1; y = 0;}
else exgcd(b,a % b,d,y,x),y -= (a / b) * x;
}
int inv(int a,int P){
int d,x,y; exgcd(a,P,d,x,y);
return (x % P + P) % P;
}
int Fac(int n,int P,int Pi){
if (!n) return 1;
int ans = 1;
if (n / P){
for (int i = 2; i < P; i++)
if (i % Pi) ans = 1ll * ans * i % P;
ans = qpow(ans,n / P,P);
}
int E = n % P;
for (int i = 2; i <= E; i++)
if (i % Pi) ans = 1ll * ans * i % P;
return 1ll * ans * Fac(n / Pi,P,Pi) % P;
}
int C(int n,int m,int P,int pi){
if (m > n) return 0;
int a = Fac(n,P,pi),b = Fac(m,P,pi),c = Fac(n - m,P,pi),k = 0,ans;
for (int i = n; i; i /= pi) k += i / pi;
for (int i = m; i; i /= pi) k -= i / pi;
for (int i = n - m; i; i /= pi) k -= i / pi;
ans = 1ll * a * inv(b,P) % P * inv(c,P) % P * qpow(pi,k,P) % P;
return 1ll * ans * (md / P) % md * inv(md / P,P) % md;
}
int exlucas(int n,int m){
int ans = 0;
for (int i = 1; i <= pi; i++){
ans = (ans + C(n,m,pk[i],p[i])) % md;
}
return ans;
}
int main(){
md = read(); n = read(); m = read();
for (int i = 1; i <= m; i++){
sum += (a[i] = read());
if (sum > n) {puts("Impossible"); return 0;}
}
Sp();
int ans = exlucas(n,sum);
for (int i = 1; i <= m; i++)
ans = 1ll * ans * exlucas(sum,a[i]) % md,sum -= a[i];
printf("%d\n",ans);
return 0;
}BZOJ2142 禮物 【擴展Lucas】