[Luogu P3702] [BZOJ 4818] [SDOI2017]序列计数

最新推荐文章于 2025-08-15 12:18:37 发布

LPA20020220

最新推荐文章于 2025-08-15 12:18:37 发布

阅读量221

点赞数

CC 4.0 BY-SA版权

分类专栏：容斥矩阵文章标签：容斥矩阵

本文链接：https://blog.youkuaiyun.com/LPA20020220/article/details/85059226

矩阵同时被 2 个专栏收录

10 篇文章

订阅专栏

容斥

6 篇文章

订阅专栏

本文深入解析洛谷和BZOJ上一道关于寻找满足特定条件的序列数量的问题，通过矩阵快速幂的方法解决序列求和需为某数倍数且包含至少一个质数的挑战。

摘要生成于 C知道，由 DeepSeek-R1 满血版支持，前往体验 >

题目描述

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

Alice还希望，这 $n$ 个数中，至少有一个数是质数。

Alice想知道，有多少个序列满足她的要求。

输入输出格式

输入格式：

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

输出格式：

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

输入输出样例

输入样例#1：

3 5 3

输出样例#1：

说明

对 $20\%$ 的数据， $1\leq n,m\leq100$

对 $50\%$ 的数据， $1\leq m \leq 100$

对 $80\%$ 的数据， $1\leq m\leq 10^6$

对 $100\%$ 的数据， $1\leq n \leq 10^9,1\leq m \leq 2\times 10^7,1\leq p\leq 100$

解题分析

保证有质数这个条件很难弄，我们就用所有的方案减去没有质数的方案就可以了。

设 $d p [i] [j]$ 表示前 $i$ 个数选取总和模 $p$ 余 $j$ 的方案数之和，显然有 $dp[i][j]=\sum_{k=0}^{p-1}dp[i-1][(j-k+p)\ mod\ p]*cnt[k]$ ，其中 $c n t [k]$ 为范围内模 $p$ 余 $k$ 的数的数量。

然后我们发现，似乎这玩意的转移很套路？那我们可以写成矩阵的形式，下图是 $p = 3$ 的一个样例：

$\left[ \begin{matrix} dp[0][0] & dp[0][1] & dp[0][2] \end{matrix} \right] \times \left[ \begin{matrix} cnt[0] & cnt[1] & cnt[2] \\ cnt[2] & cnt[0] & cnt[1] \\ cnt[1] & cnt[2] & cnt[0] \end{matrix} \right] = \left[ \begin{matrix} dp[1][0] & dp[1][1] & dp[1][2] \end{matrix} \right]$
就可以中间矩阵快速幂搞搞，然后就 $A$ 了。

代码如下：

#include <cstdio>
#include <cmath>
#include <cstdlib>
#include <cctype>
#include <algorithm>
#include <cstring>
#define R register
#define IN inline
#define W while
#define gc getchar()
#define ll long long
#define MOD 20170408ll
#define MX 20005000
struct Matrix
{
	ll mat[105][105];
	Matrix () {std::memset(mat, 0, sizeof(mat));}
} unit, st1, st2, tran1, tran2;
int pcnt, n, m, p;
int pri[MX], cnt1[105], cnt2[105];
bool npr[MX];
IN Matrix operator * (const Matrix &x, const Matrix &y)
{
	Matrix ret;
	R int i, j, k;
	for (i = 0; i < p; ++i)
	for (j = 0; j < p; ++j)
	for (k = 0; k < p; ++k)
	ret.mat[i][j] += x.mat[i][k] * y.mat[k][j];
	for (i = 0; i < p; ++i)
	for (j = 0; j < p; ++j)
	ret.mat[i][j] %= MOD;
	return ret;
}
IN Matrix fpow(Matrix base, R int tim)
{
	Matrix ret = unit;
	W (tim)
	{
		if (tim & 1) ret = ret * base;
		base = base * base, tim >>= 1;
	}
	return ret;
}
void pre()
{
	cnt2[1 % p]++;
	R int i, j, tar;
	for (i = 2; i <= m; ++i)
	{
		if (!npr[i]) pri[++pcnt] = i;
		else cnt2[i % p]++;
		for (j = 1; j <= pcnt; ++j)
		{
			tar = i * pri[j];
			if (tar > m) break;
			npr[tar] = true;
			if (!(i % pri[j])) break;
		}
	}
	int base = m / p, sur = m % p;
	for (i = 0; i < p; ++i) cnt1[i] = base, unit.mat[i][i] = 1;
	for (i = 1; i < p; ++i) cnt1[i] += (i <= sur);
	for (j = 0; j < p; ++j)
	{
		for (i = 0; i < p; ++i)
		{
			tran1.mat[i][j] = cnt1[(j + p - i) % p];
			tran2.mat[i][j] = cnt2[(j + p - i) % p];
		}
	}
	for (i = 0; i < p; ++i) st1.mat[0][i] = cnt1[i], st2.mat[0][i] = cnt2[i];
}
int main(void)
{
	scanf("%d%d%d", &n, &m, &p); pre();
	printf("%lld", ((st1 * fpow(tran1, n - 1)).mat[0][0] - (st2 * fpow(tran2, n - 1)).mat[0][0] + MOD) % MOD);
}