4818: [Sdoi2017]序列计数

考虑至少含一个质数的方案等于所有方案 $-$不含质数方案
然后考虑如何求，看到$1e9$很容易想到矩阵快速幂，
然后转移也就出来了，对于不含质数的方案跑一边线性筛即可
c++代码如下：

#include<bits/stdc++.h>
#define rep(i,x,y) for(register int i = x ; i <= y ; ++ i)
#define repd(i,x,y) for(register int i = x ; i >= y ; -- i)
using namespace std;
typedef long long ll;
template<typename T>inline void read(T&x)
{
    char c;int sign = 1;x = 0;
    do { c = getchar(); if(c == '-') sign = -1; }while(!isdigit(c));
    do { x = x * 10 + c - '0'; c = getchar(); }while(isdigit(c));
    x *= sign;
}

const int mod = 20170408;
int n,m,p;
struct Matrix
{
    int a[101][101];

    inline void clear() { memset(a,0,sizeof a); }
    Matrix operator * (Matrix x)
    {
        Matrix c;c.clear();
        rep(i,0,p-1) rep(k,0,p-1) rep(j,0,p-1)
            c.a[i][j] = (c.a[i][j] + 1ll*a[i][k]*x.a[k][j] )% mod;
        return c;
    }
    Matrix operator *= (Matrix x)
    {   
        *this = *this * x;
    }   
}a;

inline void quick_pow(Matrix&a,int y)
{
    Matrix ans;
    ans.clear();
    rep(i,0,p-1) ans.a[i][i] = 1;
    while(y)
    {
        if(y&1) ans *= a;
        a *= a;
        y >>= 1;
    }
    a = ans;
}

bool vis[20100000];
int cnt[101],tot,prime[2000000];
int main()
{
    read(n); read(m); read(p);
    rep(i,1,m) ++cnt[i%p];
    rep(i,0,p-1) rep(j,0,p-1)
        a.a[i][j] = cnt[(i+p-j)%p]; 

    quick_pow(a,n);
    int ans1 = a.a[0][0];   
    rep(i,0,p) cnt[i] = 0;
    rep(i,2,m)
    {
        if(!vis[i]) prime[++tot] = i;
        else ++cnt[i%p];
        for(register int j = 1; prime[j] * i <= m;++j)
        {
            vis[prime[j]*i] = 1;
            if(i%prime[j] == 0) break;
        }
    }

    a.clear();++cnt[1%p];
    rep(i,0,p-1) rep(j,0,p-1)
        a.a[i][j] = cnt[(i+p-j)%p]; 

    quick_pow(a,n);

    int ans2 = a.a[0][0];
    cout << (ans1 - ans2 + mod) %mod << endl; 

    return 0;
}