Bzoj4818 [Sdoi2017]序列计数
Submit: 504 Solved: 328
Description
Alice想要得到一个长度为n的序列,序列中的数都是不超过m的正整数,而且这n个数的和是p的倍数。Alice还希望
,这n个数中,至少有一个数是质数。Alice想知道,有多少个序列满足她的要求。
Input
一行三个数,n,m,p。
1<=n<=10^9,1<=m<=2×10^7,1<=p<=100
Output
一行一个数,满足Alice的要求的序列数量,答案对20170408取模。
Sample Input
3 5 3
Sample Output
33
HINT
Source
DP 矩阵乘法
n和m都很大,p那么小,想来就是矩阵加速DP了
预处理出从一个值转移到另一个值的方案数,存到矩阵里,然后大力自乘n次。
答案等于全部方案数-不用质数的方案数
1 /*by SilverN*/ 2 #include<iostream> 3 #include<algorithm> 4 #include<cstring> 5 #include<cstdio> 6 #include<cmath> 7 #define LL long long 8 using namespace std; 9 const int mod=20170408; 10 const int mxm=20000005; 11 const int mxn=105; 12 int read(){ 13 int x=0,f=1;char ch=getchar(); 14 while(ch<'0' || ch>'9'){if(ch=='-')f=-1;ch=getchar();} 15 while(ch>='0' && ch<='9'){x=x*10+ch-'0';ch=getchar();} 16 return x*f; 17 } 18 namespace shaker{ 19 int pri[mxm],cnt=0; 20 bool vis[mxm]; 21 void init(int mxm){ 22 vis[1]=1; 23 for(int i=2;i<mxm;i++){ 24 if(!vis[i]) pri[++cnt]=i; 25 for(int j=1;j<=cnt && (LL)pri[j]*i<(LL)mxm;j++){ 26 vis[pri[j]*i]=1; 27 if(i%pri[j]==0)break; 28 } 29 } 30 return; 31 } 32 } 33 // 34 int n,m,p; 35 int w[mxn]; 36 struct Mat{ 37 int x[mxn][mxn]; 38 Mat operator * (const Mat &b)const{ 39 Mat res; 40 for(int i=0;i<p;i++) 41 for(int j=0;j<p;j++){ 42 res.x[i][j]=0; 43 for(int k=0;k<p;k++){ 44 res.x[i][j]=(res.x[i][j]+(LL)x[i][k]*b.x[k][j])%mod; 45 } 46 } 47 return res; 48 } 49 void clear(){memset(x,0,sizeof x);return;} 50 }mp; 51 LL ans=0; 52 void ksm(int k){ 53 Mat res=mp; 54 while(k){ 55 if(k&1)mp=mp*res; 56 res=res*res; 57 k>>=1; 58 } 59 return; 60 } 61 void get1(){ 62 mp.clear(); 63 for(int i=0;i<p;i++){//枚举一个%p==i的数 64 for(int j=0;j<p;j++){//加一个%p==j的数 65 int ed=(i+j)%p; 66 mp.x[i][ed]+=w[j]; 67 if(mp.x[i][ed]>=mod)mp.x[i][ed]-=mod; 68 } 69 } 70 ksm(n-1); 71 return; 72 } 73 int main(){ 74 // freopen("count.in","r",stdin); 75 // freopen("count.out","w",stdout); 76 using namespace shaker; 77 int i,j; 78 n=read();m=read();p=read(); 79 init(m+1); 80 for(register int i=1;i<=m;i++) ++w[i%p]; 81 get1(); 82 ans=mp.x[0][0]; 83 for(register int i=1;i<=cnt && pri[i]<=m;i++) --w[pri[i]%p]; 84 get1(); 85 ans=(ans-mp.x[0][0]+mod)%mod; 86 printf("%lld\n",ans); 87 return 0; 88 }
本文为博主原创文章,转载请注明出处。