bzoj4818 SDOI2017 序列计数

题目链接

思路

先考虑暴力\(dp\)，\(f[i][j]\)表示前\(i\)个数，数字之和模\(p\)余\(j\)的方案数。
我们先不考虑必须有质数这个条件，先统计出全部方案。然后再减去没有质数的方案就行了。
那么就有\(f[i + 1][(j + k) \% p] += f[i][j](1\le k \le m)\)
然后发现这个其实并不需要\(o(m)\)的转移，因为\((j + k) \% p\)对于每个余数来说，最少有\(\lfloor\frac{m}{p}\rfloor\)个。然后有\(m\%p\)个数有\(\lfloor\frac{m}{p}\rfloor + 1\)个.只要对这\(m\%p\)个数单独处理一下就行了。
然后计算没有质数的方案数，全部用合数转移即可。
然后发现这个是可以矩阵优化的。
转移矩阵的第\(i\)行第\(j\)列表示从\(1 + i\)到\(m+i\)中有多少个数字模\(p\)余\(j\)。同样可以计算一下，然后\(o(p^2)\)的列出矩阵。
对于全部用合数转移的情况，转移矩阵可以\(o(p\frac{m}{lnm})\)转移,亲测吸氧可过

代码

/*
* @author: wxyww
* @date:   2019-02-13 16:22:09
* @last modified time: 2019-02-13 20:51:05
*/
#include<cstdio>
#include<iostream>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#include<queue>
#include<vector>
#include<ctime>
using namespace std;
typedef long long ll;
// #define int ll
const int mod = 20170408;
ll read() {
    ll x=0,f=1;char c=getchar();
    while(c<'0'||c>'9') {
        if(c=='-') f=-1;
        c=getchar();
    }
    while(c>='0'&&c<='9') {
        x=x*10+c-'0';
        c=getchar();
    }
    return x*f;
}
int n,m,p;
namespace bf2 {
    const int n = 110;
    struct node {
        int a[n][n],n,m;
        node() {
            memset(a,0,sizeof(a));n = 0,m = 0;
        }
        node(int nn,int mm) {
            memset(a,0,sizeof(a));n = nn;m = mm;
        }
        node(int nn) {
            n = nn;
            memset(a,0,sizeof(a));
            for(int i = 0;i <= nn;++i) a[i][i] = 1;
        }
    };
    node operator * (const node &x,const node &y) {
        int n = x.n,m = y.m;
        node c(n,m);
        int k = x.m;
        for(int i = 0;i <= k;++i) {
            for(int j = 0;j <= n;++j) {
                for(int z = 0;z <= m;++z) {
                    c.a[j][z] += 1ll * x.a[j][i] * y.a[i][z] % mod;
                    c.a[j][z] %= mod;
                }
            }
        }
        return c;
    }
    node operator ^ (node x,int y) {
        node ans(x.n);
        for(;y;y >>= 1,x = x * x)
            if(y & 1) ans = ans * x;
        return ans;
    }
    int dis[2000003],vis[20000003];
    int tot,tmp[2000003];
    int js;
    void eur() {
        vis[1] = 1;
        for(int i = 2;i <= m;++i) {
            if(!vis[i]) dis[++js] = i;
            for(int j = 1;1ll * dis[j] * i <= m;++j) {
                vis[dis[j] * i] = 1;
                if(dis[j] % i == 0) break;
            }
        }
    }
    void main() {
        node c(p - 1,p - 1);
        int k = m % p;
        for(int i = 0;i < p;++i) {
            for(int j = i + 1;j <= i + k;++j) {
                c.a[i][j % p]++;
            }
        }
        for(int i = 0;i < p;++i) for(int j = 0;j < p;++j) c.a[i][j] += m / p;
        node ans(0,p - 1);
        ans.a[0][0] = 1;
        ans = ans * (c ^ n);
        int anss = ans.a[0][0];
        memset(ans.a,0,sizeof(ans.a));
        memset(c.a,0,sizeof(c.a));
        ans.a[0][0] = 1;
        eur();
        for(int i = 0;i < p;++i) {
            for(int j = i + 1;j <= i + k;++j) {
                c.a[i][j % p]++;
            }
        }
        for(int i = 0;i < p;++i) for(int j = 0;j < p;++j) c.a[i][j] += m / p;
        for(int i = 0;i < p;++i)
            for(int j = 1;j <= js;++j)
                c.a[i][(i + dis[j]) % p]--;
        ans = ans * (c ^ n);
        cout<<(anss - ans.a[0][0] + mod) % mod;
    }

}
signed main() {
    n = read(),m = read(),p = read();
    bf2::main();
    return 0;
}