洛谷3702，递推+矩阵乘法优化

传送门
这题一开始想出来了一个二元一次递推式，有点麻烦。
一看题解，发现可以转为求全集和答案的补集，挺妙的。
一开始筛素数的数组只开到106，RE两个点，开到5∗106就AC了，感觉很玄学，毕竟值域m≤2∗107

#include<bits/stdc++.h>using namespace std;typedef long long ll;const ll mo=20170408;struct matrix{    ll a[110][110];    int n;    void clear(){        memset(a,0,sizeof a);    }    matrix operator*(const matrix&x){        matrix ans;        ans.clear();        for(int i=0;i<n;++i)            for(int j=0;j<n;++j)                for(int k=0;k<n;++k)ans.a[i][k]=(ans.a[i][k]+a[i][j]*x.a[j][k])%mo;        ans.n=n;        return ans;    }}a,b,c,d;matrix pow(matrix a,int b){    if(b==1)return a;    matrix u=pow(a,b>>1);    if(b&1)return u*u*a;        else return u*u;}ll n,m,p,x[110],y[110],pp[5000010],xb,i,j,z,ans;bool bb[20000010];int main(){    cin>>n>>m>>p;    y[0]=m/p;    for(i=1;i<=m%p;++i)y[i]=m/p+1;    for(;i<p;++i)y[i]=m/p;    for(i=2;i<=m;++i){        if(!bb[i])pp[++xb]=i,++x[i%p];        for(j=1;j<=xb;++j){            z=pp[j]*i;            if(z>m)break;            bb[z]=1;            if(i%pp[j]==0)break;        }    }    if(n==1)return printf("%d\n",!bb[p] && p<=m),0;    d.n=b.n=p;    for(i=0;i<p;++i)        for(j=0;j<p;++j){            z=(j-i+p)%p;            b.a[i][j]=y[z]-x[z];            d.a[i][j]=y[z];        }    a=pow(b,n-1);    c=pow(d,n-1);    for(j=0;j<p;++j)        ans=(ans+mo+y[j]*c.a[j][0]%mo-(y[j]-x[j])*a.a[j][0]%mo)%mo;    printf("%lld\n",ans);    return 0;}