【BZOJ4818】[Sdoi2017]序列计数

Description

Alice想要得到一个长度为n的序列，序列中的数都是不超过m的正整数，而且这n个数的和是p的倍数。Alice还希望

，这n个数中，至少有一个数是质数。Alice想知道，有多少个序列满足她的要求。

Input

一行三个数，n,m,p。

1<=n<=10^9,1<=m<=2×10^7,1<=p<=100

Output

一行一个数，满足Alice的要求的序列数量，答案对20170408取模。

Sample Input

3 5 3

Sample Output

题解：至少包含1个质数的数量=总数-不包含质数的数量（这种补集法不是一次两次见到了吧？）

于是我们考虑用DP求解，先快筛出1..m内的质数，1..m内除以P模为j的数的个数，1..m内除以P模为j的合数的个数

然后设f[i][j]表示i个数，总和除以P模j的方案数，g[i][j]表示i个合数，总和除以P模j的方案数，容易得出

f[i+1][(j+k)%P]+=f[i][j]+1..m内除以P模为j的数的个数
g[i+1][(j+k)%P]+=g[i][j]+1..m内除以P模为j的合数的个数

发现时间复杂度O(np)，用矩乘快速幂优化一下就好啦

#include <cstdio>#include <cstring>#include <iostream>#define mod 20170408using namespace std;typedef long long ll;int np[20000010],cnt[110],sum[110],pri[10000010];int n,m,p,tot;typedef struct matrix{	ll v[110][110];}M;M x,ans,emp;ll ans1;M mmul(M a,M b){	M c=emp;	int i,j,k;	for(i=0;i<p;i++)		for(j=0;j<p;j++)			for(k=0;k<p;k++)				c.v[i][j]=(c.v[i][j]+a.v[i][k]*b.v[k][j])%mod;	return c;}void pm(int y){	while(y)	{		if(y&1)	ans=mmul(ans,x);		x=mmul(x,x),y>>=1;	}}int main(){	scanf("%d%d%d",&n,&m,&p);	int i,j;	np[1]=cnt[1]=sum[1]=1;	for(i=2;i<=m;i++)	{		sum[i%p]=(sum[i%p]+1)%mod;		if(!np[i])	pri[++tot]=i;		else	cnt[i%p]=(cnt[i%p]+1)%mod;		for(j=1;j<=tot&&i*pri[j]<=m;j++)		{			np[i*pri[j]]=1;			if(i%pri[j]==0)	break;		}	}	for(i=0;i<p;i++)		for(j=0;j<p;j++)			x.v[i][(i+j)%p]=(x.v[i][(i+j)%p]+sum[j])%mod;	ans.v[0][0]=1;	pm(n);	ans1=ans.v[0][0];	memset(ans.v,0,sizeof(ans.v)),memset(x.v,0,sizeof(x.v));	ans.v[0][0]=1;	for(i=0;i<p;i++)		for(j=0;j<p;j++)			x.v[i][(i+j)%p]=(x.v[i][(i+j)%p]+cnt[j])%mod;	pm(n);	printf("%lld",(ans1-ans.v[0][0]+mod)%mod);	return 0;}

【BZOJ4818】[Sdoi2017]序列计数 DP+矩阵乘法

声明：以上内容来自用户投稿及互联网公开渠道收集整理发布，本网站不拥有所有权，未作人工编辑处理，也不承担相关法律责任，若内容有误或涉及侵权可进行投诉：投诉/举报工作人员会在5个工作日内联系你，一经查实，本站将立刻删除涉嫌侵权内容。

联系
我们

首页 > 代码库 > 【BZOJ4818】[Sdoi2017]序列计数 DP+矩阵乘法