【BZOJ4818】序列计数（SDOI2017）-矩阵优化DP+线性筛质数

测试地址：序列计数
做法：本题需要用到矩阵优化DP+线性筛质数。
我们很快就能想到一个状态定义：令f(i,j)$f(i,j)$为i$i$个数的数列中，数字和对p$p$取模的结果为j$j$的方案数，那么我们很快就能得出状态转移方程。又根据这个状态转移方程，每一个状态都只和上一层的一些状态有关，而且每一层的状态数较少，所以又想到使用矩阵优化，那么就可以做到O(p3logn)$O(p^3\log n)$的时间复杂度了。
然而上面还忽略了一个条件：至少要选择一个质数。对于质数我们显然可以先O(m)$O(m)$线性筛求出，但是我们要不要在状态中多加一维来表示选没选过质数呢？那样比较麻烦，其实我们只要用总方案数减去不选择质数的方案数就行了，这样我们就可以沿用上面的定义做两次DP，然后就解决了这一题。
以下是本人代码：

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long ll;
const ll mod=20170408;
int n,m,p,prime[20000010],primetot=0;
ll ans=0,tot[110];
bool pr[20000010]={0};
struct matrix
{
    ll mat[110][110];
}M[35],S;

void calc_prime()
{
    pr[1]=1;
    for(int i=2;i<=m;i++)
    {
        if (!pr[i]) prime[++primetot]=i;
        for(int j=1;j<=primetot&&i*prime[j]<=m;j++)
        {
            pr[i*prime[j]]=1;
            if (i%prime[j]==0) break;
        }
    }
}

void mult(matrix A,matrix B,matrix &S)
{
    memset(S.mat,0,sizeof(S.mat));
    for(int i=0;i<p;i++)
        for(int j=0;j<p;j++)
            for(int k=0;k<p;k++)
            {
                S.mat[i][j]+=A.mat[i][k]*B.mat[k][j];
                if (S.mat[i][j]>=mod) S.mat[i][j]%=mod;
            }
}

void power(matrix &S,int x)
{
    memset(S.mat,0,sizeof(S.mat));
    for(int i=0;i<p;i++)
        S.mat[i][i]=1;
    int i=0;
    while(x)
    {
        if (x&1) mult(S,M[i],S);
        x>>=1;i++;
    }
}

void build(int type)
{
    memset(tot,0,sizeof(tot));
    for(int i=1;i<=m;i++)
        if (!type||pr[i]) tot[i%p]++;
    for(int i=0;i<p;i++)
        for(int j=0;j<p;j++)
            M[0].mat[i][j]=tot[(i-j+p)%p];
    for(int i=1;i<=32;i++)
        mult(M[i-1],M[i-1],M[i]);
    power(S,n);
    if (!type) ans=(ans+S.mat[0][0])%mod;
    else ans=(ans-S.mat[0][0]+mod)%mod;
}

int main()
{
    scanf("%d%d%d",&n,&m,&p);
    calc_prime();

    build(0);
    build(1);
    printf("%lld",ans);

    return 0;
}