题解：[SDOI2017] 序列计数

题目传送门

题意分析

首先，“至少一个数是质数”就很不好判断，考虑容斥掉。答案即 不考虑是否为质数的答案 减去 没有质数的答案。

观察到 \(1\leq p\leq100\)，考虑复杂度与 \(p\) 正相关的算法。

设 \(f_{i,j}\) 表示第 \(1\sim i\) 个数，和模 \(p\) 为 \(j\) 时的不考虑是否为质数答案。没有质数的答案为 \(g_{i,j}\)。答案即：

\[(f_{n,0}-g_{n,0})\bmod20170408 \]

设 \(\textit{cntF}_{i}\) 表示 \(1\sim m\) 中模 \(p\) 余 \(i\) 的数的个数，\(\textit{cntG}_{i}\) 表示 \(1\sim m\) 中模 \(p\) 余 \(i\) 的合数的个数。

则有：

\[f_{i,j}=\sum_{k=0}^{p-1}f_{i-1,k}\cdot \textit{cntF}_{j-k}\\ g_{i,j}=\sum_{k=0}^{p-1}g_{i-1,k}\cdot \textit{cntG}_{j-k}\\ \]

即将两个序列对位相加，满足和模 \(p\) 余 \(0\) 即可。

但是这样复杂度会炸掉（\(\mathcal O\left(np^2\right)\)）。容易发现，序列里的数不定，因此一半的答案和后一半的答案一样，因此可以分治。

即枚举前一半的序列和后一半的序列，只要相加的和模 \(p\) 余 \(0\) 即可：

\[f_{i,j}=\sum_{k=0}^{p-1}f_{\frac i2,k}f_{\frac i2,(j-k)\bmod p}\\ g_{i,j}=\sum_{k=0}^{p-1}g_{\frac i2,k}g_{\frac i2,(j-k)\bmod p}\\ \]

注意 \(i\) 为奇数时需要特判。

---

同时，滚动数组可以将空间复杂度优化至线性。

AC 代码

//#include<bits/stdc++.h>
#include<algorithm> 
#include<iostream>
#include<cstring>
#include<iomanip>
#include<cstdio>
#include<string>
#include<vector>
#include<cmath>
#include<ctime>
#include<deque>
#include<queue>
#include<stack>
#include<list>
using namespace std;
constexpr const int N=1e9,M=2e7,P=20170408;
int vis[M+1],prime[M+1],primeSize;
void pre(){
	for(int i=2;i<=M;i++){
		if(!vis[i]){
			vis[i]=i;
			prime[++primeSize]=i;
		}
		for(int j=1;j<=primeSize&&1ll*i*prime[j]<=M;j++){
			vis[i*prime[j]]=prime[j];
			if(i%prime[j]==0){
				break;
			}
		}
	}
}
int p,f[2][M+1],g[2][M+1],cntF[M+1],cntG[M+1];
void solve(int x,bool mode){
	if(x==1){
		for(int i=0;i<p;i++){
			f[mode][i]=cntF[i];
			g[mode][i]=cntG[i];
		}
		return;
	}
	solve(x>>1,!mode);
	for(int i=0;i<p;i++){
		f[mode][i]=g[mode][i]=0;
	} 
	for(int i=0;i<p;i++){
		for(int j=0;j<p;j++){
			int pl=i-j;
			if(i-j<0){
				pl+=p;
			}
			f[mode][i]=(f[mode][i]+1ll*f[!mode][j]*f[!mode][pl]%P)%P;
			g[mode][i]=(g[mode][i]+1ll*g[!mode][j]*g[!mode][pl]%P)%P;
		}
	}
	if(x&1){
		for(int i=0;i<p;i++){
			f[!mode][i]=f[mode][i];
			g[!mode][i]=g[mode][i];
		}
		for(int i=0;i<p;i++){
			f[mode][i]=g[mode][i]=0;
		}
		for(int i=0;i<p;i++){
			for(int j=0;j<p;j++){
				int pl=i-j;
				if(pl<0){
					pl+=p;
				}
				f[mode][i]=(f[mode][i]+1ll*f[!mode][j]*cntF[pl]%P)%P;
				g[mode][i]=(g[mode][i]+1ll*g[!mode][j]*cntG[pl]%P)%P;
			}
		}
	}
}
int query(int n,int m,int p){
	for(int i=1;i<=m;i++){
		cntF[i%p]++;
		cntG[i%p]+=(vis[i]!=i);
	}
	solve(n,0);
	int ans=(f[0][0]-g[0][0])%P;
	if(ans<0){
		ans+=P;
	}
	return ans;
}
int main(){
	/*freopen("test.in","r",stdin);
	freopen("test.out","w",stdout);*/
	
	ios::sync_with_stdio(false);
	cin.tie(0);cout.tie(0);
	
	pre();
	int n,m;
	cin>>n>>m>>p;
	cout<<query(n,m,p)<<'\n';
	
	cout.flush(); 
	
	/*fclose(stdin);
	fclose(stdout);*/
	return 0;
}