題目描述
![](https://img.laitimes.com/img/9ZDMuAjOiMmIsIjOiQnIsICdzFWRoRXdvN1LclHdpZXYyd2LcBzNvwVZ2x2bzNXak9CX90TQNNkRrFlQKBTSvwFbslmZvwFMwQzLcVmepNHdu9mZvwFVywUNMZTY18CX052bm9CX9UEROJza65EMVpXTmZEWjZXUYpVd1kmYr50MZV3YyI2cKJDT29GRjBjUIF2LcRHelR3LcJzLctmch1mclRXY39DN2cDOwgDN3EDOxUDM3EDMy8CX0Vmbu4GZzNmLn9Gbi1yZtl2Lc9CX6MHc0RHaiojIsJye.jpg)
題目解析
呵呵,實際這是NOIP題你信嗎?
如果不考慮包含質數的情況,因為p很小,我們不妨建構p*p的矩陣,然後對于(i,j)填有多少個數x使得(i+x)%p=j,最後矩陣快速幂一下就行了。
如果考慮質數我們便可以求出所有方案數再減去不包含質數的方案數。
代碼
#include<cstdio>
#include<iostream>
#include<cstring>
#include<cmath>
#include<cstdlib>
#include<algorithm>
#include<queue>
#include<stack>
using namespace std;
#define MAXN
#define MAXM 20000000
#define INF 0x3f3f3f3f
typedef long long int LL;
template<class T>
void Read(T &x){
x=;char c=getchar();bool flag=;
while(c<'0'||'9'<c){if(c=='-')flag=;c=getchar();}
while('0'<=c&&c<='9'){x=x*+c-'0';c=getchar();}
if(flag)x=-x;
}
const int MOD = ;
struct Matrix{
int n;
LL data[][];
Matrix(){}
Matrix(int _n){n=_n;memset(data,,sizeof(data));}
Matrix operator + (const Matrix &x)const{
Matrix ret=Matrix(n);
for(int i=;i<n;++i)
for(int j=;j<n;++j)
ret.data[i][j]=data[i][j]+x.data[i][j];
return ret;
}
Matrix operator - (const Matrix &x)const{
Matrix ret=Matrix(n);
for(int i=;i<n;++i)
for(int j=;j<n;++j)
ret.data[i][j]=data[i][j]-x.data[i][j];
return ret;
}
Matrix operator * (const Matrix &x)const{
Matrix ret=Matrix(n);
for(int i=;i<n;++i)
for(int j=;j<n;++j)
for(int k=;k<n;++k)
ret.data[i][j]+=data[i][k]*x.data[k][j];
return ret;
}
Matrix operator % (const int mod)const{
Matrix ret=Matrix(n);
for(int i=;i<n;++i)
for(int j=;j<n;++j)
ret.data[i][j]=data[i][j]%mod;
return ret;
}
};
Matrix quick_power(Matrix a,int p,int mod){
Matrix ret=Matrix(a.n);
for(int i=;i<ret.n;++i)ret.data[i][i]=;
while(p){
if(p&)ret=ret*a%mod;
a=a*a%mod;
p>>=;
}
return ret;
}
bool isprime[MAXM+];
int prime[MAXM+],cnt;
void init(int ed){
memset(isprime,,sizeof(isprime));
cnt=;
isprime[]=;
for(int i=;i<=ed;++i){
if(isprime[i])prime[++cnt]=i;
for(int j=;j<=cnt&&prime[j]*i<=ed;++j){
isprime[prime[j]*i]=;
if(i%prime[j]==)break;
}
}
}
int change[];
int main(){
//freopen("count.in","r",stdin);
//freopen("count.out","w",stdout);
int n,m,p;
Read(n),Read(m),Read(p);
init(m);
Matrix origin=Matrix(p);
memset(change,,sizeof(change));
for(int i=;i<=m;++i)++change[(i%p)];
for(int i=;i<p;++i)
for(int j=;j<p;++j)origin.data[i][j]=change[((j-i)%p+p)%p]%MOD;
Matrix ans1=quick_power(origin,n,MOD);
origin=Matrix(p);
memset(change,,sizeof(change));
for(int i=;i<=m;++i)if(!isprime[i])++change[(i%p)];
for(int i=;i<p;++i)
for(int j=;j<p;++j)origin.data[i][j]=change[((j-i)%p+p)%p]%MOD;
Matrix ans2=quick_power(origin,n,MOD);
printf("%I64d\n",(ans1.data[][]-ans2.data[][]+MOD)%MOD);
}
/*
1 2 1
*/