记录编号 |
337443 |
评测结果 |
AAAAAAAAAA |
题目名称 |
[HZOI 2016]定约servant |
最终得分 |
100 |
用户昵称 |
kito |
是否通过 |
通过 |
代码语言 |
C++ |
运行时间 |
5.587 s |
提交时间 |
2016-11-04 15:27:16 |
内存使用 |
0.31 MiB |
显示代码纯文本
#include<cstdio>
#include<cmath>
#include<vector>
#include<ctime>
using namespace std;
#define fcl fclose(stdin); fclose(stdout); return 0
#define SUBMIT 2333
int m,p,phi,Siz,n;
vector<int> Mu[15],Inv[15];
int prime[15],pp[15],tot=0,mi[15];
int AA[15];
/***********************Matrix****************************/
struct Matrix{
int a[4][4];
void Init(int k){
for(int i=1;i<3;++i){
for(int j=1;j<3;++j){
a[i][j]=0;
}
}
for(int i=1;i<3;++i) a[i][i]=k;
}
}X,Q,F;
Matrix operator * (const Matrix& A,const Matrix& B){
X.Init(0);
for(int i=1;i<3;++i){
for(int j=1;j<3;++j){
for(int k=1;k<3;++k){
X.a[i][j]=(X.a[i][j]+(A.a[i][k]*1ll*B.a[k][j])%phi)%phi;
}
}
}
return X;
}
Matrix QMatrixPow(Matrix a,int b){
Matrix ans;
ans.Init(1);
while(b){
if(b&1){
ans=ans*a;
}
b>>=1;
a=a*a;
}
return ans;
}
void GetMatrix(){
Q.a[1][1]=0; Q.a[1][2]=1;
Q.a[2][1]=1; Q.a[2][2]=1;
Matrix P=QMatrixPow(Q,n-1);
F.a[1][1]=(P.a[1][1]+P.a[2][1])%phi;
F.a[1][2]=(P.a[1][2]+P.a[2][2])%phi;
Siz=(F.a[1][1]*1ll*F.a[1][2])%phi;
}
/***********************Matrix****************************/
int Qpow(int a,int b,int mod){
int c=1;
while(b){
if(b&1){
c=(c*1ll*a)%mod;
}
b>>=1;
a=(a*1ll*a)%mod;
}
return c;
}
void ExGcd(int a,int b,int& x,int& y){
if(b==0){
x=1; y=0;
return;
}
ExGcd(b,a%b,x,y);
int t=x;
x=y; y=t-(a/b)*y;
}
int Gcd(int a,int b){
return b==0?a:Gcd(b,a%b);
}
bool Divprime(){
bool flag=false;
phi=p;
int x=(int)sqrt(double(p)),mod=p;
for(int i=2;i<=x;++i){
if(mod%i==0){
phi/=i;
phi*=(i-1);
prime[++tot]=i;
pp[tot]=1;
while(mod%i==0){
mi[tot]++;
mod/=i;
pp[tot]*=i;
}
if(mi[tot]>1) flag=true;
}
if(mod==1) break;
}
if(mod!=1){
phi/=mod;
phi*=(mod-1);
prime[++tot]=mod;
pp[tot]=mod;
mi[tot]=1;
}
return flag;
}
void Init(){
int x;
for(int j=1;j<=tot;++j){
if(n>p) x=prime[j]-1;
else x=n;
Mu[j].push_back(1); Inv[j].push_back(1);
for(int i=1;i<=x;++i){
Mu[j].push_back((Mu[j][i-1]*1ll*i)%p);
Inv[j].push_back(0);
}
int y;
ExGcd(Mu[j][x],prime[j],Inv[j][x],y);
Inv[j][x]%=p; Inv[j][x]=(Inv[j][x]+p)%p;
for(int i=x-1;i>=0;--i){
Inv[j][i]=(Inv[j][i+1]*1ll*(i+1))%p;
}
}
}
int C(int a,int b,int i){
if(a<b) return 0;
int x=(Mu[i][a]*1ll*Inv[i][b])%prime[i];
x=(x*1ll*Inv[i][a-b])%prime[i];
return x;
}
int Lucas(int a,int b,int i){
if(b==0) return 1;
return (Lucas(a/prime[i],b/prime[i],i)*1ll*C(a%prime[i],b%prime[i],i))%prime[i];
}
void Union(int j,int i){
int a=pp[j],b=pp[i],c=AA[i]-AA[j],x,y;
ExGcd(a,b,x,y);
x=(x*1ll*c)%b;
x=(x+b)%b;
pp[i]*=pp[j];
AA[i]=((x*1ll*pp[j])%pp[i]+AA[j])%pp[i];
}
int CRT(){
for(int i=2;i<=tot;++i){
Union(i-1,i);
}
return AA[tot];
}
void work1(){
GetMatrix();
//处理模数拆分后每个质因子的幂为1的情况
Init();
for(int i=1;i<=tot;++i){
for(int j=1;j<=m+1;++j){
if(Gcd(j,m)==1){
AA[i]=(AA[i]+Lucas(n,j-1,i))%prime[i];
}
}
}
int x=CRT();
printf("%d\n",Qpow(x,Siz,p));
}
int INV(int a,int mod){
int x,y;
ExGcd(a,mod,x,y);
x%=mod; x=(x+mod)%mod;
return x;
}
int size;
vector<int> STD;
int Fact(int N,int i){
if(N==0) return 1;
int x=Qpow(STD[pp[i]],N/pp[i],pp[i]);
x=(x*1ll*STD[N%pp[i]])%pp[i];
size+=N/prime[i];
return (x*1ll*Fact(N/prime[i],i))%pp[i];
}
void work2(){
GetMatrix();
//处理p的素因子指数大于1的情况
int x,y;
for(int i=1;i<=tot;++i){
STD.clear();
STD.push_back(1);
for(int j=1;j<=pp[i];++j){
if(j%prime[i]){
STD.push_back((STD[j-1]*1ll*j)%pp[i]);
}
else STD.push_back(STD[j-1]);
}
for(int j=1;j<=m+1;++j){
if(Gcd(j,m)==1){
x=1; y=0;
size=0; x=(x*1ll*Fact(n,i))%pp[i]; y+=size;
size=0; x=(x*1ll*INV(Fact(j-1,i),pp[i]))%pp[i]; y-=size;
size=0; x=(x*1ll*INV(Fact(n-j+1,i),pp[i]))%pp[i]; y-=size;
x=(x*1ll*Qpow(prime[i],y,pp[i]))%pp[i];
AA[i]=(AA[i]+x)%pp[i];
}
}
}
x=CRT();
printf("%d",Qpow(x,Siz,p));
}
int main(){
#ifdef SUBMIT
freopen("servant.in","r",stdin);
freopen("servant.out","w",stdout);
#endif
scanf("%d%d%d",&n,&m,&p);
if(!Divprime()) work1();
else work2();
#ifndef SUBMIT
getchar(); getchar();
#endif
fcl;
}