bzoj 4818: [Sdoi2017]序列计数

利用补集转化思想
带素数方案数=无限制方案数-只给合数方案数
用cnt[i]表示1-m中%p=i的数各个数
得到朴素动态规划
$dp[i][j]=sum_{k=0}^{p-1} dp[i-1][(j-k)%p]*cnt[k]$
复杂度O(np^2)
对于cnt乘法运算,矩乘优化一下
复杂度O(logn*p^2)

/*
 * dp可行方案数-dp没有素数的方案数
 * */
#include<cstdio>
#include<cstring>
#include<algorithm>

#define mod 20170408
inline int  read() {
    int x=0,f=1;
    char c=getchar(); 
    while(c<'0'||c>'9') { if(c=='-')f=-1;c=getchar();} 
    while(c<='9'&&c>='0') x=x*10+c-'0',c=getchar(); 
    return x*f; 
}
#define mod 20170408 
#define LL long long 
int n,m,p,num=0;
struct Matrix { 
    int f[207]; 
    Matrix() {memset(f,0,sizeof f);} 
}a,b,c,d;
Matrix operator * (Matrix &a ,Matrix &b) { 
    Matrix c,d;  
    for(int i=0;i<p;++i)  
        for(int j=0;j<p;++j) { 
            c.f[i+j]=(c.f[i+j]+1LL*a.f[i]*b.f[j])%mod;  
        }
    for(int i=p-1;i>=0;--i) d.f[i]=c.f[i]+c.f[i+p];
    return d;
}
int prime[20000007];bool vis[20000007];
void get_prime() {
    vis[1]=true;
    for(int i=2;i<=m;++i) {
        if(!vis[i]) prime[++num]=i;
        for(int j=1;j<=num&&i*prime[j]<=m;++j) {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)break;
        }
    }
}
int main() { 
    n=read(),m=read(),p=read();
    get_prime();
    for(int i=1;i<=m;++i) {
        a.f[i%p]++; if(vis[i]) b.f[i%p]++;
    }
    for(c.f[0]=d.f[0]=1;n;n>>=1,a=a*a,b=b*b) 
        if(n&1) c=c*a,d=d*b;
    printf("%d
",(c.f[0]-d.f[0]+mod)%mod);
    return 0;
}