BZOJ4818 [Sdoi2017]序列计数

来源：互联网发布：星际争霸2mac版本编辑：程序博客网时间：2024/06/05 06:27

看数据范围一眼矩乘

f[i][j]表示前i个数，和模p余j

用任意取的方案减去不取质数的方案即得到至少有一个质数的方案

对任取和不取质数分别随便做个转移矩阵快速幂即可

#include<iostream>#include<cstring>#include<ctime>#include<cmath>#include<algorithm>#include<iomanip>#include<cstdlib>#include<cstdio>#include<map>#include<bitset>#include<set>#include<stack>#include<vector>#include<queue>using namespace std;#define MAXN 100#define MAXM 20000010#define ll long long#define eps 1e-8#define MOD 20170408#define INF 1000000000int P;struct mat{ll x[MAXN];mat(){memset(x,0,sizeof(x));}mat operator *(const mat &y){mat z;int i,j;for(i=0;i<P;i++){for(j=0;j<P;j++){z.x[i]+=x[j]*y.x[(i-j+P)%P];if(z.x[i]>1e18){z.x[i]%=MOD;}}z.x[i]%=MOD;}return z;}};int n,m;int p[MAXM],tot;bool np[MAXM];mat Xa,Xn,Ra,Rn;int ca[MAXN],cn[MAXN];int ans;void su(){int i,j;np[1]=1;for(i=2;i<MAXM;i++){if(!np[i]){p[++tot]=i;}for(j=1;(ll)i*p[j]<MAXM;j++){np[i*p[j]]=1;if(!(i%p[j])){break;}}}}int main(){int i,j;su();scanf("%d%d%d",&n,&m,&P);for(i=1;i<=m;i++){Xa.x[i%P]++;if(np[i]){Xn.x[i%P]++;}}Ra.x[0]=Rn.x[0]=1;while(n){if(n&1){Ra=Ra*Xa;Rn=Rn*Xn;}Xa=Xa*Xa;Xn=Xn*Xn;n>>=1;}printf("%d\n",(Ra.x[0]-Rn.x[0]+MOD)%MOD);return 0;}/**/

0 0