【洛谷5366】[SNOI2017] 遗失的答案(状压DP)

这个我省三年前的省选真题我肝了将近一天 ……

完了我感觉我省选要凉了怎么办在线等急

虽然我还不知道什么时候省选

这题我写了两种做法。

我这么懒的人怎么会写两种做法呢?

因为这题数据范围太!诡!异!了!

下来再说。

题目

洛谷 5366

分析

显然,如果 \(L\) 不是 \(G\) 的倍数、\(X>N\)\(X\) 不是 \(L\) 的因数或 \(X\) 不是 \(G\) 的倍数则无解,否则可以把题转换成求从 \([1,\lfloor\frac{N}{G}\rfloor]\) 中选出若干个数,使得最大公因数是 \(1\) 且最小公倍数是 \(\frac{L}{G}\) 的方案中,有多少种方案选了 \(\frac{X}{G}\) 。以下的 \(L\) 均为原题中的 \(\frac{L}{G}\)\(X\) 均为原题中的 \(\frac{X}{G}\)

\(x=\prod_j p_j^{a_j}\) (称为 \(x\) 的标准分解),其中 \(p_j\) 是质数,那么有:

\[\gcd(\{x_i\})=\prod_j p_j^{\min_i(\{a_{i,j}\})} \]

\[\mathrm{lcm}(\{x_i\})=\prod_j p_j^{\max_i(\{a_{i,j}\})} \]

(这个应该非常显然吧)

因此,在最终选出的数的标准分解中,每个质数的最小指数应该是 \(0\) ,最大指数应该等于 \(L\) 的标准分解中该质数的指数。

不超过 \(10^8\) 的数最多只能有 \(8\) 个不同的质因数,所以考虑设计状态 \(S\) ,每个质数用两位分别表示是否已经选过标准分解中该质数的指数为 \(0\) 和指数与 \(L\) 的指数相等的数。类似地,给每个数确定一个属性表示能满足哪些条件,最终状态就是选出的数的属性或起来。

显然,只有 \(L\) 的不超过 \(2\sqrt{L}\) 个因数是可以选的,因此可以直接暴力枚举 \(L\) 的每个质因子来确定这些因数的属性。(我想了一个多小时怎么给 \(10^8\) 个数确定属性未果后才意识到只有 \(2\sqrt{L}\) 个有用 …… 菜死了。)虽然属性看似有 \(2^{8\times 2}=65536\) 种,但事实上由于一个数不可能同时满足一个质数的两个条件(也就是不存在一个质数的两位都是 \(1\) 的属性),而只有当 \(L\) 的标准分解中这个质数的指数不小于 \(2\) 时才可能存在一个数同时不满足两个条件(即一个质数的两位都是 \(0\) )。因此事实上不同的属性数是很少的。通过如下 Python 代码发现最多可能只有 \(576\) 个。

(为什么说可能呢?因为我真的不会写 Python ,怎么循环都是现查的。我知道下面这个代码很丑很菜看起来像刚学两天编程的人写的,甚至怀疑正确性不太对。)

ans = 0
p = [2, 3, 5, 7, 11, 13, 17, 19]
now = 1
pos = 8
N = 1e8

for i in p :
	now = now * i
i = 0
while i < pos :
	ans = max(ans, 3 ** i * 2 ** (pos - i))
	now = now / p[i]
	now = now * p[i] * p[i]
	while now > N :
		pos = pos - 1
		now = now / p[pos]
	i = i + 1
print(ans)

接下来有两条路可以走。

方法一

这个做法思路清晰解法自然,所以我想到就直接写了。

将各种属性按任意顺序排成一个序列,设 \(\mathrm{pre}_{i,S}\) 表示前 \(i\) 种属性的数组成了集合 \(S\) 的状态,\(\mathrm{suf}_{i,S}\) 表示后 \(i\) 种属性的数组成了集合 \(S\) 的状态。这两个数组可以直接 DP 求出。一种属性贡献的方案数是 \(2^k-1\) ,其中 \(k\) 是这种属性的数的个数,减 \(1\) 是减去这些数都不取的情况。注意特判属性 \(0\) 即使不取也对 \(S\) 没有影响,方案数是 \(2^k\)

考虑预处理钦定第 \(i\) 种属性的某个数必须选的情况的答案。首先把 \(\mathrm{pre}_{i-1}\)\(\mathrm{suf}_{i+1}\) 合并起来。观察发现这一步是一个或卷积,可以用 FWT 优化。然后再考虑第 \(i\) 种属性的没有被钦定的数,一共有 \(2^{k-1}\) 种选法。最后,如果所有没有被钦定的数组成的状态或上第 \(i\) 种属性可成为全 \(1\) 的状态,则可计入答案。询问时求出 \(X\) 的属性然后直接输出该属性被钦定一个数时的答案即可。代码见下文。

这个做法的时间复杂度 \(O(Mp2^{2p})\) ,其中 \(M\) 是属性的种类数,\(p\)\(L\) 的质因数种类数。最坏情况 \(M=576,p=8\) ,算下来大概 \(6\times 10^8\)

虽然网上题解都说常数小可以过,但是我花了一晚上 + 一上午时间使尽浑身解数卡常后洛谷上还!是!过!不!去!

我也不知道某 s 姓选手是怎么带了个火车头就卡进去的。

方法二

这个做法思路更清晰解法更自然,我也不知道我为什么想不到,可能就是太菜了吧。

先介绍两个小技巧

子集反演

\[f_S=\sum_{T\subset S}g_T \]

\[g_S=\sum_{T\subset S}(-1)^{(|S|-|T|)}f(T) \]

等等这不就是容斥原理吗 …… 这都有名字系列.jpg 。

子集求和

\[f_S=\sum_{T\subset S}g_T \]

在已知 \(g\) 的情况下求 \(f\)

当然可以直接暴力枚举子集,时间复杂度是 \(O(3^n)\) 的,其中 \(n\) 是集合大小。但如下代码可以做到 \(O(n2^n)\) :

memcpy(f, g, sizeof(int[1 << n]));
for (int i = 0; i < n; i++)
	for (int j = 0; j < (1 << n); j++)
		if (j & (1 << i))
			f[j] += f[j ^ (1 << i)];

我第一次看到这个技巧的时候震惊了,竟然改变两层循环的顺序就能防止重复计算 ……

简单证明一下。设 \(f_{i,S}\) 表示前 \(i\) 位是 \(S\) 的子集,其他位和 \(S\) 一样的集合的 \(g\) 之和。那么从 \(i-1\) 转移到 \(i\) 时,新增的就是前 \(i\) 位是 \(S\) 的子集,第 \(i\) 位是 \(0\) (即跟 \(S\) 的第 \(i\) 位不一样)且其他位和 \(S\) 一样的集合之和,也就是 \(f_{i,S\setminus\{i\}}\) 。滚动存储一下就是上面的代码。

回到这道题。仍然是预处理钦定第 \(i\) 种属性的某个数必须选的情况的答案。用 \(f_S\) 表示钦定第 \(i\) 种属性的某个数必须选时选出状态 \(S\) 的方案数。由于第 \(i\) 种属性已经钦定选了一个,所以只有当 \(S\) 是该属性的超集时才有意义。

从所有有意义的 \(S\) 中删去属性 \(i\) ,就变成了 \(i\) 的补集的全体子集,因此可以直接套用子集卷积的结论,即(设属性 \(i\)\(P_i\) ):

\[f_{S\cup P_i}=\sum_{T\subset S}(-1)^{(|S|-|T|)}g(T\cup P_i) \]

此时 \(g_S\) 的定义应该是「钦定第 \(i\) 种属性的某个数必须选时选出状态 \(S\) 或其子集的方案数」。那么 \(g_S=2^{k-1}\) ,其中 \(k\) 是属性为 \(S\) 的子集的数的个数(显然钦定的那个数也在其中,所以要减掉),可以用上述子集求和的方法预处理出。

处理询问就和方法一一样了。时间复杂度 \(O(3^{2p})\) ,其中 \(p\)\(L\) 的质因数种类数。

代码

方法一

LOJ 上 1372 ms ,洛谷 TLE 一个点。

#include <algorithm>
#include <cstdio>
#include <cctype>
#include <cstring>
using namespace std;

namespace zyt
{
	template<typename T>
	inline bool read(T &x)
	{
		char c;
		bool f = false;
		x = 0;
		do
			c = getchar();
		while (c != EOF && c != '-' && !isdigit(c));
		if (c == EOF)
			return false;
		if (c == '-')
			f = true, c = getchar();
		do
			x = x * 10 + c - '0', c = getchar();
		while (isdigit(c));
		if (f)
			x = -x;
		return true;
	}
	template<typename T>
	inline void write(T x)
	{
		static char buf[20];
		char *pos = buf;
		if (x < 0)
			putchar('-'), x = -x;
		do
			*pos++ = x % 10 + '0';
		while (x /= 10);
		while (pos > buf)
			putchar(*--pos);
	}
	typedef long long ll;
	const int N = 600, S = 1 << 16, PR = 10, P = 1e9 + 7, INV2 = (P + 1) / 2;
	int n, g, l, q, pow[N], p[PR], pk[PR], pcnt, val[N], id[S], num[N], cnt;
	int pre[N][S], suf[N][S], tot[N][S], ans[N];
	int power(int a, int b)
	{
		int ans = 1;
		while (b)
		{
			if (b & 1)
				ans = (ll)ans * a % P;
			a = (ll)a * a % P;
			b >>= 1;
		}
		return ans;
	}
	void get_prime(int n)
	{
		for (int i = 2; i * i <= n; i++)
			if (n % i == 0)
			{
				p[pcnt] = i, pk[pcnt] = 1;
				while (n % i == 0)
					n /= i, pk[pcnt] *= i;
				++pcnt;
			}
		if (n > 1)
			p[pcnt] = pk[pcnt] = n, pcnt++;
	}
	int get_type(const int n)
	{
		int ans = 0;
		for (int i = 0; i < pcnt; i++)
		{
			if (n % pk[i] == 0)
				ans |= (1 << i);
			else if (n % p[i])
				ans |= (1 << (i + pcnt));
		}
		return ans;
	}
	namespace FWT
	{
		void fwt(int *const a, const int n)
		{
			for (int l = 1; l < n; l <<= 1)
				for (int i = 0; i < n; i += (l << 1))
					for (int k = i; k < i + l; k++)
						a[l + k] = (a[l + k] + a[k]) % P;
		}
		void ifwt(int *const a, const int n)
		{
			for (int l = 1; l < n; l <<= 1)
				for (int i = 0; i < n; i += (l << 1))
					for (int k = i; k < i + l; k++)
						a[l + k] = (a[l + k] - a[k] + P) % P;
		}
	}
	int work()
	{
		using namespace FWT;
		read(n), read(g), read(l), read(q);
		if (l % g)
		{
			while (q--)
				puts("0");
			return 0;
		}
		l /= g, n /= g;
		get_prime(l);
		for (int i = 1; i * i <= l && i <= n; i++)
			if (l % i == 0)
			{
				int t = get_type(i);
				if (!id[t])
					id[t] = cnt, num[cnt] = t, val[cnt] = 1, cnt++;
				val[id[t]] = val[id[t]] * 2 % P;
				if (i * i < l && l / i <= n)
				{
					t = get_type(l / i);
					if (!id[t])
						id[t] = cnt, num[cnt] = t, val[cnt] = 1, cnt++;
					val[id[t]] = val[id[t]] * 2 % P;
				}
			}
		pre[0][0] = 1, pre[0][num[0]] += val[0] - 1;
		fwt(pre[0], 1 << (pcnt * 2));
		for (int i = 1; i < cnt; i++)
		{
			if (!num[i])
				pre[i][0] = val[i];
			else
				pre[i][0] = 1, pre[i][num[i]] = val[i] - 1;
			fwt(pre[i], 1 << (pcnt * 2));
			for (int j = 0; j < (1 << (pcnt * 2)); j++)
				pre[i][j] = (ll)pre[i][j] * pre[i - 1][j] % P;
		}
		suf[cnt - 1][0] = 1, suf[cnt - 1][num[cnt - 1]] += val[cnt - 1] - 1;
		fwt(suf[cnt - 1], 1 << (pcnt * 2));
		for (int i = cnt - 2; i >= 0; i--)
		{
			if (!num[i])
				suf[i][0] = val[i];
			else
				suf[i][0] = 1, suf[i][num[i]] = val[i] - 1;
			fwt(suf[i], 1 << (pcnt * 2));
			for (int j = 0; j < (1 << (pcnt * 2)); j++)
				suf[i][j] = (ll)suf[i][j] * suf[i + 1][j] % P;
		}
		for (int i = 0; i < cnt; i++)
		{
			static int tmp[S];
			if (!i)
				memcpy(tmp, suf[i + 1], sizeof(int[1 << (pcnt * 2)]));
			else if (i == cnt - 1)
				memcpy(tmp, pre[i - 1], sizeof(int[1 << (pcnt * 2)]));
			else
				for (int j = 0; j < (1 << (pcnt * 2)); j++)
					tmp[j] = (ll)pre[i - 1][j] * suf[i + 1][j] % P; 
			ifwt(tmp, 1 << (pcnt * 2));
			for (int j = 0; j < (1 << (pcnt * 2)); j++)
			{
				tot[i][j | num[i]] = (tot[i][j | num[i]] + tmp[j]) % P;
				if ((j | num[i]) == (1 << (pcnt * 2)) - 1)
					ans[i] = (ans[i] + tot[i][j]) % P;
			}
			ans[i] = (ll)ans[i] * val[i] % P * INV2 % P;
		}
		while (q--)
		{
			int x;
			read(x);
			if (x % g || (x / g) > n || l % (x / g))
				write(0);
			else
			{
				int t = get_type(x / g);
				write(ans[id[t]]);
			}
			putchar('\n');
		}
		return 0;
	}
}
int main()
{
	return zyt::work();
}

方法二

LOJ 上 58 ms ,洛谷上 201 ms 。

#include <algorithm>
#include <cstdio>
#include <cctype>
#include <cstring>
using namespace std;

namespace zyt
{
	template<typename T>
	inline bool read(T &x)
	{
		char c;
		bool f = false;
		x = 0;
		do
			c = getchar();
		while (c != EOF && c != '-' && !isdigit(c));
		if (c == EOF)
			return false;
		if (c == '-')
			f = true, c = getchar();
		do
			x = x * 10 + c - '0', c = getchar();
		while (isdigit(c));
		if (f)
			x = -x;
		return true;
	}
	template<typename T>
	inline void write(T x)
	{
		static char buf[20];
		char *pos = buf;
		if (x < 0)
			putchar('-'), x = -x;
		do
			*pos++ = x % 10 + '0';
		while (x /= 10);
		while (pos > buf)
			putchar(*--pos);
	}
	typedef long long ll;
	const int N = 600, S = 1 << 16, PR = 10, P = 1e9 + 7, SQRT= 1e4 + 10;
	int n, g, l, q, p[PR], pk[PR], pcnt, val[S], id[S], num[N], cnt, pow[SQRT], count[S], ans[N];
	void get_prime(int n)
	{
		for (int i = 2; i * i <= n; i++)
			if (n % i == 0)
			{
				p[pcnt] = i, pk[pcnt] = 1;
				while (n % i == 0)
					n /= i, pk[pcnt] *= i;
				++pcnt;
			}
		if (n > 1)
			p[pcnt] = pk[pcnt] = n, pcnt++;
	}
	int get_type(const int n)
	{
		int ans = 0;
		for (int i = 0; i < pcnt; i++)
		{
			if (n % pk[i] == 0)
				ans |= (1 << i);
			else if (n % p[i])
				ans |= (1 << (i + pcnt));
		}
		return ans;
	}
	void init()
	{
		pow[0] = 1;
		for (int i = 1; i < SQRT; i++)
			pow[i] = (ll)pow[i - 1] * 2 % P;
		for (int i = 0; i < S; i++)
			count[i] = count[i >> 1] + (i & 1);
	}
	int work()
	{
		init();
		read(n), read(g), read(l), read(q);
		if (l % g)
		{
			while (q--)
				puts("0");
			return 0;
		}
		l /= g, n /= g;
		get_prime(l);
		for (int i = 1; i * i <= l && i <= n; i++)
			if (l % i == 0)
			{
				int t = get_type(i);
				if (!id[t])
					id[t] = cnt, num[cnt] = t, cnt++;
				++val[t];
				if (i * i < l && l / i <= n)
				{
					t = get_type(l / i);
					if (!id[t])
						id[t] = cnt, num[cnt] = t, cnt++;
					++val[t];
				}
			}
		for (int i = 0; i < pcnt * 2; i++)
			for (int j = 0; j < (1 << (pcnt * 2)); j++)
				if (j & (1 << i))
					val[j] += val[j ^ (1 << i)];
		for (int i = 0; i < cnt; i++)
		{
			int t = ((1 << (pcnt * 2)) - 1) ^ num[i];
			for (int j = t; j; j = ((j - 1) & t))
				ans[i] = ((ll)ans[i] + (count[j | num[i]] & 1 ? -1 : 1) * pow[val[j | num[i]] - 1] + P) % P;
			ans[i] = ((ll)ans[i] + (count[num[i]] & 1 ? -1 : 1) * pow[val[num[i]] - 1] + P) % P;
		}
		while (q--)
		{
			int x;
			read(x);
			if (x % g || (x / g) > n || l % (x / g))
				write(0);
			else
			{
				int t = get_type(x / g);
				write(ans[id[t]]);
			}
			putchar('\n');
		}
		return 0;
	}
}
int main()
{
	return zyt::work();
}
posted @ 2020-04-23 22:12  Inspector_Javert  阅读(140)  评论(0编辑  收藏  举报