【洛谷5398】[Ynoi2018]GOSICK（二次离线莫队）

题目：

洛谷 5398

当我刚学莫队的时候，他们告诉我莫队能解决几乎所有区间问题；

现在，当我发现一个区间问题似乎难以用我所了解的莫队解决的时候，他们就把这题的正解叫做 XXX 莫队。——题记

（以上皆为瞎扯，纯属虚构，请勿当真）

分析：

先转化一下题目：如果允许每次询问都暴力把区间扫一遍，那么每扫到一个数 (i) ，就统计已经扫过的部分中有多少个 (j) 满足 (a_j) 是 (a_i) 的因数（即取数对 ((i,j)) ）或倍数 （即取数对 ((j,i)) ）。注意，因为 ((i,j)) 和 ((j,i)) 是不同的，如果 (a_i=a_j) ，那么答案要加 (2) 。

考虑莫队。每次往当前区间加入或删除一个元素的时候，需要统计当前区间有多少个元素是该元素的因数或倍数（为了和题目中的「询问」区分，我们不妨把这个操作简称为一个数对一个区间「统计」），然后给答案加上或减去这个数量。这样单次修改的复杂度是 (O(n)) ，总时间复杂度高达 (O(n^2sqrt{n})) 。 GG

黑科技来了 ：由于每次当前区间两端点的移动是已知的，所以我们把两端点移动时出现的「统计」离线下来处理（是谓「二次离线」）。

如何离线呢？我们发现每次都是一个数 (a) 对一个区间 ([l,r]) 统计，而这个操作又可以拆分为对 ([1,l)) 和 ([1,r]) 两个前缀统计。

然后开始大力分类讨论 ……

第一，当左端点从 (l) 左移到 (l') 时，是每个 (iin[l',l)) 对 ([i,r]) 统计，即对 ([1,i)) 和 ([1,r]) 统计；

第二，当左端点从 (l) 右移到 (l') 时，是每个 (iin[l,l')) 对 ([i,r]) 统计，即对 ([1,i)) 和 ([1,r]) 统计；

第三，当右端点从 (r) 左移到 (r') 时，是每个 (iin(r',r]) 对 ([l,i]) 统计，即对 ([1,l)) 和 ([1,i]) 统计；

第四，当右端点从 (r) 右移到 (r') 时，是每个 (iin(r,r']) 对 ([l,i]) 统计，即对 ([1,l)) 和 ([1,i]) 统计。

可以看出统计分为两大类。第一类是 (i) 对于 ([1,i)) （或 ([1,i]) ，但很明显这两个差距很小）统计，第二类是一个区间中的 (i) 对一个固定的前缀分别统计（如第一种情况中是所有 (iin[l',l)) 对固定的前缀 ([1,r]) 统计）。

第一类可以预处理。第二类可以按前缀从小到大排序，暴力扫当前前缀对应的所有统计（根据莫队的复杂度证明，端点移动距离之和是 (O(nsqrt{n})) ）。

事实上这两类统计可以抽象成同一个问题：维护一个结构，支持「插入一个数」和「给定 (k) ，查询已经插入的数中有多少个是 (k) 的因数或倍数」。插入次数是 (O(n)) ，查询次数是 (O(nsqrt{n})) 。

既然查询次数远比插入次数多，自然考虑维护每个数的答案。插入一个数时，因为可以在 (O(sqrt{n})) 时间内遍历一个数的所有因数，所以直接暴力修改所有因数的答案即可。同理，当 (k) 大于一个阈值，比如 (S) ，倍数也可以暴力修改，单次时间复杂度不超过 (O(frac{n}{S})) ，据说此处选择 (S=32) 。

那么当 (kleq S) 怎么办呢？这里有一个很巧妙的做法。预处理出每个数能否被 ([1,S]) 中的某个数整除（这个结果下称「状态」），并将结果压成一个 (S) 位的二进制数，并开一个大小为 (2^S) 的桶记录「这个状态的数的答案应该统一加上多少」。当插入 (k(kleq S)) 时，暴力将所有 (k) 这一位上是 (1) 的状态的答案加 (1) 。查询时不光查询前一段话中维护的答案，也要加上这一段话中对应状态的答案。由于 (2^{32}) 太大了，因此把因数分成 (8) 个一组，每个数对应 (4) 个状态。这样插入时只需要修改 (2^{8-1}) 个状态，查询时需要查询这个数对应的 (4) 个状态。

还有一个小小的细节。由于默认任意两个相同的数会给答案贡献 (2) （原因在文首说了），如果 (i) 对 ([l,r]) 统计，而 (iin[l,r]) ，那么 ((i,i)) 这个数对就被算了两次。没关系，我们只需要放任它算两次，最后给所有询问的答案减去区间长度即可。也由于这个原因，(i) 对 ([1,i]) 统计的答案应该是 (i) 对 ([1,i)) 统计的答案加 (2) 。

代码：

注意每次算出的是答案相对于上次的变化量，所以最后要按照处理询问的顺序做一遍前缀和。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cctype>
#include <vector>
#include <cmath>
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 = 1e5 + 10, M = 1e5, Q = N, SS = 8, S = 32;
	const bool ADD = true, SUB = false;
	int n, m, block, blnum, belong[N], arr[N], fac[N][4], facnum[4][1 << SS], now[N], f[N];
	ll ans[Q];
	struct _ask
	{
		int l, r, id;
		bool operator < (const _ask &b) const
		{
			return belong[l] == belong[b.l] ? r < b.r : belong[l] < belong[b.l];
		}
	}ask[Q];
	struct node
	{
		int l, r, id;
		bool type;
		node(const int _l, const int _r, const int _id, const bool _t)
			: l(_l), r(_r), id(_id), type(_t) {}
	};
	vector<node> v[N];
	void insert(const int a)
	{
		for (int i = 1; i * i <= a; i++)
			if (a % i == 0)
			{
				++now[i];
				if (i * i < a)
					++now[a / i];
			}
		if (a > S)
		{
			for (int i = a; i <= M; i += a)
				if (i > S)
					++now[i];
		}
		else
		{
			int bl = (a - 1) / SS, pos = (a - 1) % SS;
			for (int i = 0; i < (1 << SS); i++)
				if (i & (1 << pos))
					++facnum[bl][i];
		}
	}
	int query(const int a)
	{
		int ans = now[a];
		for (int i = 0; i < 4; i++)
			ans += facnum[i][fac[a][i]];
		return ans;
	}
	int work()
	{
		read(n), read(m);
		block = sqrt(n), blnum = ceil(double(n) / double(block));
		for (int i = 1; i <= M; i++)
			for (int j = 0; j < 4; j++)
				for (int k = j * SS + 1; k <= (j + 1) * SS; k++)
					if (i % k == 0)
						fac[i][j] |= (1 << (k - j * SS - 1));
		for (int i = 1; i <= n; i++)
			read(arr[i]), belong[i] = (i - 1) / block + 1;
		for (int i = 1; i <= m; i++)
			read(ask[i].l), read(ask[i].r), ask[i].id = i;
		sort(ask + 1, ask + m + 1);
		for (int i = 1; i <= n; i++)
			f[i] = query(arr[i]), insert(arr[i]);
		for (int i = 1, l = 1, r = 1; i <= m; i++)
		{
			if (ask[i].l < l)
			{
				v[r].push_back(node(ask[i].l, l - 1, ask[i].id, ADD));
				while (ask[i].l < l)
					ans[ask[i].id] -= f[--l];
			}
			if (ask[i].l > l)
			{
				v[r].push_back(node(l, ask[i].l - 1, ask[i].id, SUB));
				while (ask[i].l > l)
					ans[ask[i].id] += f[l++];
			}
			if (ask[i].r > r)
			{
				v[l - 1].push_back(node(r + 1, ask[i].r, ask[i].id, SUB));
				while (ask[i].r > r)
					ans[ask[i].id] += f[++r] + 2;
			}
			if (ask[i].r < r)
			{
				v[l - 1].push_back(node(ask[i].r + 1, r, ask[i].id, ADD));
				while (ask[i].r < r)
					ans[ask[i].id] -= f[r--] + 2;
			}
		}
		memset(facnum, 0, sizeof(facnum));
		memset(now, 0, sizeof(now));
		for (int i = 1; i <= n; i++)
		{
			insert(arr[i]);
			for (vector<node>::iterator it = v[i].begin(); it != v[i].end(); it++)
				for (int j = it->l; j <= it->r; j++)
					ans[it->id] += query(arr[j]) * (it->type == ADD ? 1 : -1);
		}
		ans[0] = 2;
		for (int i = 1; i <= m; i++)
			ans[ask[i].id] += ans[ask[i - 1].id];
		for (int i = 1; i <= m; i++)
			ans[ask[i].id] -= ask[i].r - ask[i].l + 1;
		for (int i = 1; i <= m; i++)
			write(ans[i]), putchar('
');
		return 0;
	}
}
int main()
{
#ifdef BlueSpirit
	freopen("5398.in", "r", stdin);
#endif
	return zyt::work();
}