bzoj4818: [Sdoi2017]序列计数（矩阵快速幂优化dp）

传送门
题意简述：问有多少长度为n的序列，序列中的数都是不超过m的正整数，而且这n个数的和是p的倍数，且其中至少有一个数是质数，答案对$20170408$取模$(n\le 1e9,m\le 2e7,p\le 100)$。

---

思路：
首先因为只需要是$p$的倍数，因此可以看成全局和对$p$取模为$0$方案数。
设状态$f_{0/1,i,j}$表示不限制选出的数/选出的数不能是质数，允许选$i$个数，全局和对$p$取模等于$j$方案数，答案就是$f_{0,n,0}-f_{1,n,0}$。
用矩阵快速幂来优化这个$dp$即可。
代码：

#include<bits/stdc++.h>
#define ri register int
using namespace std;
typedef long long ll;
const int mod=20170408,M=2e7+5,N=105;
inline int add(const int&a,const int&b){return a+b>=mod?a+b-mod:a+b;}
inline int dec(const int&a,const int&b){return a>=b?a-b:a-b+mod;}
inline int mul(const int&a,const int&b){return (ll)a*b%mod;}
int n,m,p,s1[105],s2[N],pri[M],tot=0;
bool vis[M];
inline void init(){
	vis[1]=1;
	for(ri i=2;i<=m;++i){
		if(!vis[i])pri[++tot]=i;
		for(ri j=1;j<=tot&&i*pri[j]<=m;++j){
			vis[i*pri[j]]=1;
			if(i%pri[j]==0)break;
		}
	}
}
struct Mat{
	int a[N][N];
	Mat(int v=0){memset(a,0,sizeof(a));for(ri i=0;i<p;++i)a[i][i]=v;}
	friend inline Mat operator*(const Mat&a,const Mat&b){
		Mat ret;
		for(ri i=0;i<p;++i)for(ri k=0;k<p;++k)if(a.a[i][k])
			for(ri j=0;j<p;++j)ret.a[i][j]=add(ret.a[i][j],mul(a.a[i][k],b.a[k][j]));
		return ret;
	}
	friend inline Mat operator^(Mat&a,int p){
		Mat ret(1);
		for(;p;p>>=1,a=a*a)if(p&1)ret=ret*a;
		return ret;
	}
};
int main(){
	cin>>n>>m>>p;
	init();
	Mat f0,f1;
	for(ri i=1;i<=m;++i){
		++s1[i%p];
		if(vis[i])++s2[i%p];
	}
	f0.a[0][0]=f1.a[0][0]=1;
	for(ri i=0;i<p;++i)for(ri j=0;j<p;++j)f0.a[i][(i+j)%p]=s1[j]%mod,f1.a[i][(i+j)%p]=s2[j]%mod;
	f0=f0^n,f1=f1^n;
	cout<<dec(f0.a[0][0],f1.a[0][0]);
	return 0;
}