【BZOJ3309】DZY Loves Math（线性筛）

题目：

BZOJ 3309

分析：

首先，经过一番非常套路的莫比乌斯反演（实在懒得写了），我们得到：

[sum_{T=1}^n sum_{d|T}f(d)mu(frac{T}{d})lfloorfrac{n}{T} floorlfloorfrac{m}{T} floor ]

那么，我们现在如果预处理出 (g(n)=sum_{d|n}f(d)mu(frac{n}{d})) 的前缀和，就可以数论分块 (O(sqrt{n})) 处理每次询问了。

(mu) 函数有一个重要的性质：当 (n) 是某个数平方的倍数时（即某个质因子的次数不小于 (2) ）， (mu(n)=0) 。所以，真正有贡献的项只能是 (frac{n}{d}) 的质因子互不相同的项，共有 (2^k) 项，其中 (k) 是 (n) 的质因子种数（相当于每种质因子的次数要么为 (0) ，要么为 (1) ）。

设这 (k) 个质因子中有 (a) 个的次数为 (f(n)) ，则这 (2^k) 项中有 (2^{k-a}) 项的 (f(d)) 为 (f(n)-1) （即这 (a) 个质因子 全部 选入 (frac{n}{d}) 的情况），其余情况为 (f(n)) 。根据 (mu) 的定义，每一项乘上的系数与选了多少个质因子进入 (frac{n}{d}) 有关，奇数为 (-1) ，偶数为 (1) 。而这 (2^{k-a}) 项中（暂定 (k>a) ），选的质因子数为奇数、偶数的项数相等，所以和为 (0) 。同理，剩下 (2^k-2^{k-a}) 项之和也是 (0)

但是，当 (k=a) 时，(f(d)=f(n)-1) 的只有 (2^{k-a}=1) 项，绝对值为 (f(n)-1) 而不是 (0) 。同时，此时 (f(d)=f(n)) 的项也是奇数个，和的绝对值为 (f(n)) 。这两项具体的符号与 (k) 的奇偶性有关（并且一定符号相反）。稍微推一下可得，当 (k) 为奇数， (g(n)=1) ；当 (k) 为偶数， (g(n)=-1) 。

综上可得

[g(n)=egin{cases} 0 (k eq a)\ 1 (k=a且k是奇数)\ -1 (k=a且k是偶数) end{cases}]

代码：

#include <cstdio>
#include <algorithm>
#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 = 1e7 + 10;
	int prime[N], f[N], g[N], num[N], cnt[N], tot[N], pcnt;
	bool mark[N];
	void init()
	{
		f[1] = 0;
		for (int i = 2; i < N; i++)
		{
			if (!mark[i])
				prime[++pcnt] = i, f[i] = num[i] = cnt[i] = tot[i] = 1;
			for (int j = 1; j <= pcnt && (ll)i * prime[j] < N; j++)
			{
				int k = i * prime[j];
				mark[k] = true;
				if (i % prime[j] == 0)
				{
					num[k] = num[i] + 1;
					f[k] = max(f[i], num[k]);
					tot[k] = tot[i];
					if (f[i] == num[k])
						cnt[k] = cnt[i] + 1;
					else if (f[i] > num[k])
						cnt[k] = cnt[i];
					else
						cnt[k] = 1;
					break;
				}
				else
				{
					num[k] = 1;
					f[k] = max(f[i], 1);
					tot[k] = tot[i] + 1;
					if (f[i] == 1)
						cnt[k] = cnt[i] + 1;
					else
						cnt[k] = cnt[i];
				}
			}
		}
		g[1] = 0;
		for (int i = 2; i < N; i++)
			g[i] = (cnt[i] == tot[i] ? ((tot[i] & 1) ? 1 : -1) : 0) + g[i - 1];
	}
	int work()
	{
		int T;
		read(T);
		init();
		while (T--)
		{
			int n, m;
			read(n), read(m);
			if (n > m)
				swap(n, m);
			int pos = 1;
			ll ans = 0;
			while (pos <= n)
			{
				int tmp = min(n / (n / pos), m / (m / pos));
				ans += ll(g[tmp] - g[pos - 1]) * (n / pos) * (m / pos);
				pos = tmp + 1;
			}
			write(ans), putchar('
');
		}
		return 0;
	}
}
int main()
{
#ifdef BlueSpirit
	freopen("3309.in", "r", stdin);
#endif
	return zyt::work();
}