4338: BJOI2015 糖果
Time Limit: 2 Sec Memory Limit: 256 MBSubmit: 200 Solved: 93[][][]Description
Alice 正在教她的弟弟 Bob 学数学。
每天，Alice 画一个N行M 列的表格，要求 Bob在格子里填数。
Bob已经学会了自然数1到K的写法。因此他在每个格子里填1 ~ K之间的整数。
Alice 告诉 Bob，如果 Bob 填写完表格的 N*M 个数以后，每行的数从第 1 列到第 M
列单调不减，并且任意两行至少有一列的数不同，而且以前 Bob 没有填写过相同的表格，
那么Alice 就给Bob吃一颗糖果。
Bob想知道，如果每天填写一遍表格，最多能吃到多少颗糖果。
答案模P输出。
Input
第一行，四个整数依次是N, M, K, P。
Output
输出一行，一个整数，表示答案模P 后的结果。
Sample Input
【样例输入1】 1 3 3 10 【样例输入2】 2 2 2 10
Sample Output
【样例输出1】 0 【样例输出2】 6
HINT
【样例解释】
样例1。表格只有一行。每个格子可以填写1 ~ 3。有10种填写方法，依次为1 1 1，
1 1 2，1 1 3，1 2 2，1 2 3，1 3 3，2 2 2，2 2 3，2 3 3，3 3 3。
样例2。表格有两行。有6 种填写方法，依次为 1 1/1 2, 1 1/2 2, 1 2/1 1, 1 2/2
2, 2 2/1 1, 2 2/1 2。
【数据规模与约定】
100% 的数据中，1 ≤ N, M ≤ 10^5，1 ≤ P, K ≤ 2*10^9.
没错，这又是扩展卢卡斯+中国剩余定理，是不是已经被我做烂了233333
很容易发现的是一行的方案是可重组合，C(K+M-1，M)；
然后因为有N行，行之间不能相同，所以是排列，P(C(K+M-1，M),N)
求一下这个就行了。。。。
然后会发现的是，P可能是一个大质数或者含有一个大质数，所以我特判了一下，如果P中含有大于M的次数为1的质因子的话，那么就直接用组合数的某一行递推公式来算C;否则就用扩展卢卡斯。
至于算P方法都是一样的2333
还有我为什么成了这个题的rank1了，，，(害怕
#include#define ll long long#define maxn 100005using namespace std;int N,M,K,P,MOD;int d[15],D[15];int mo,phi[15];int ans[15],num;int jc[maxn],inv[maxn];inline int add(int x,int y,const int ha){ x+=y; if(x>=ha) return x-ha; else return x;}inline int ksm(int x,int y,const int ha){ int an=1; for(;y;y>>=1,x=x*(ll)x%ha) if(y&1) an=an*(ll)x%ha; return an;}struct node{ int val,tmp; node operator *(const node &U)const{ return (node){val*(ll)U.val%D[mo],tmp+U.tmp}; } node operator /(const node &U)const{ return (node){val*(ll)ksm(U.val,phi[mo]-1,D[mo])%D[mo],tmp-U.tmp}; }};inline void dvd(){ for(int i=2;i*(ll)i<=P;i++) if(!(P%i)){ d[++num]=i,D[num]=1; while(!(P%i)) P/=i,D[num]*=i; phi[num]=D[num]/d[num]*(d[num]-1); if(P==1) break; } if(P!=1) d[++num]=D[num]=P,phi[num]=P-1;}inline node getjc(int x){ node now=(node){1,0}; if(x>=d[mo]) now=now*getjc(x/d[mo]),now.tmp+=x/d[mo]; if(x>=D[mo]) now=now*(node){ksm(jc[D[mo]-1],x/D[mo],D[mo]),0}; now=now*(node){jc[x%D[mo]],0}; return now;}inline node getC(int x,int y){ return getjc(x)/getjc(y)/getjc(x-y);}inline int getP(int x,int y,const int ha){ int now=x; for(int i=2;i<=y;i++){ now=add(now,ha-1,ha); x=x*(ll)now%ha; } return x;}inline void solve(int x){ mo=x,jc[0]=1; const int ha=D[x]; if(d[x]==D[x]&&d[x]>M){ inv[1]=1; for(int i=2;i<=M;i++) inv[i]=-inv[ha%i]*(ll)(ha/i)%ha+ha; ans[x]=1; int now=K+M-1; for(int i=1;i<=M;i++,now=add(now,ha-1,ha)) ans[x]=ans[x]*(ll)now%ha*(ll)inv[i]%ha; ans[x]=getP(ans[x],N,ha); } else{ for(int i=1;i