【ybt金牌导航5-2-5】【UOJ#33-UR #2】树上GCD

树上GCD

题目链接：ybt金牌导航5-2-5 / UOJ#33-UR #2

题目大意

给你一个以 1 为根的树，问你有对于每个 1~n-1 的 i，有多少个点对，它们分别到它们 LCA 的距离的 GCD 是 i。

思路

(color{white}{吐槽：这是一道神（e）奇（xin）的题目，整个过程操作都非常的神（xuan）奇（xue）。})

首先我们看到这些路径长度之类的，自然想到点分治。

但你会发现统计 GCD 不是很知道要怎么搞，而这里有一个方法，就是容斥。
首先你求出 GCD 是 (x) 的倍数的个数（这个其实好求，你只需要让两个距离都是 (x) 的倍数就行），然后我们进行容斥，倒序枚举数，然后再枚举它的倍数，那后面是都处理好的，就是 GCD 是 (ix) 的个数，那你要求 GCD 是 (x) 的，那就是用 GCD 是 (x) 倍数的减去 GCD 是 (2x,3x,4x,...) 的个数，然后这么 (O(nlogn)) 搞一遍，就好了。

那你考虑怎么搞之前，你会发现一个很特殊的东西，就是 (gcd(0,x)=x)，那这是什么情况呢？就是其中一个是 (LCA)，就是说在树上是一条链。
那你容易想到可以单独统计这种情况，先求出有多少个点到根节点的距离是 (x)，那容易想到以这些点为链的下面，链的上面可以在这个链的任意一个位置，那就会有 (x,x-1,x-2,...,1) 这么多种情况的距离。
那你就倒叙枚举后缀和加一下就行了。

那接着我们就来搞一个点分治的部分。
然后你会发现它是一个有根树，我们不能像在无根树上面那么乱搞。
那我们考虑搞有根树点分治，要怎么搞呢？我们继续分情况讨论。
假设你现在处理一个点，你要枚举它子树嘛，那如果你选的两个点都在原来树的子树，那就是像无根树那样搞。（第一种情况）

至于怎么搞，我们要四个数组 (num_i,sum_i,Tnum_i,Tsum_i)。
(num_i) 是有多少个深度为 (i) 的点。（对于当前枚举的子树）(Tnum_i) 则是对于你之前枚举过的子树。（这样防止算重，(Tsum_i) 也一样）
(sum_i) 是有多少个深度是 (i) 倍数的点。
这两个都可以跑一遍子树得到。
那容易想到，我们就是每次把 (sum_i imes Tsum_i) 统计如 (ans_i) 中，然后把它们加进 (Tnum_i,Tsum_i) 里。

但如果一个不是呢？也就是说，它在原来的树中是这个树的父亲。（第二种情况）
那你就直接暴力跳这些点，跳到你之前的当前子树根节点，找到它子树的点们，然后进行配对。
跟配对就是跟在原来树的所有子树中的点配对。
那前面我们把一开始求出的 (Tnum_i) 保留，然后算出现在你跳到点的子树的 (sum_i)，然后再配对。
为什么不能用之前的 (Tsum_i) 呢？因为你现在跳到的点到根节点是有距离的，那就是说它一开始会间隔一段的距离，那就跟之前间隔 (0) 不同了。那你就要重新求一次。
那容易想到这个复杂度会裂开。

那你考虑记忆化，因为间隔的次数 (x) 顶多是你枚举的求 GCD 为 (x) 倍数的 (x)。
但容易看到 (x) 是 (10^5) 级别，会超空间。
那你考虑每次你枚举 (x) 之后，你后面加时每次加复杂度是 (O(frac{md}{x}))（(md) 是子树的最深深度），那你容易想到分块。
对于前 (sqrt{n}) 我们记忆化，对于 (>sqrt{n}) 我们就直接暴力求，后面这一部分求就是 (O(md log md)) 了，就可以接受了。

关于具体的实现可以看看代码。

代码

#include<queue>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#define INF 0x3f3f3f3f3f3f3f3f
#define ll long long

using namespace std;

struct node {
	int to, nxt;
}e[400001];
int n, fa[200001], le[200001], KK, SUM;
int deg[200001], root, sz[200001], lstdeg;
int f[200001], num[200001], sum[200001];
int Tnum[200001], Tsum[200001], maxdeg;
int rem[1001][1001], UP, tmp, noww, lim;
ll ans1[200001], ans2[200001], tmpp;
bool in[200001];
queue <int> q;

int read() {
	int re = 0;
	char c = getchar();
	while (c < '0' || c > '9') c = getchar();
	while (c >= '0' && c <= '9') {
		re = (re << 3) + (re << 1) + c - '0';
		c = getchar();
	}
	return re;
}

void write(ll x) {
	if (x > 9) write(x / 10);
	putchar(x % 10 + '0');
}

void add(int x, int y) {
	e[++KK] = (node){y, le[x]}; le[x] = KK;
	e[++KK] = (node){x, le[y]}; le[y] = KK;
}

void work_line() {//处理其中一个到 LCA 是 0 距离的情况
	for (int i = 1; i <= n; i++) {
		deg[i] = deg[fa[i]] + 1;
		ans1[deg[i] - 1]++;//因为你根节点深度是 1，那距离是就深度减一
	}
}

void get_root(int now, int father) {//点分治找重心
	sz[now] = 1;
	f[now] = 0;
	for (int i = le[now]; i; i = e[i].nxt)
		if (!in[e[i].to] && e[i].to != father) {
			get_root(e[i].to, now);
			sz[now] += sz[e[i].to];
			f[now] = max(f[now], sz[e[i].to]);
		}
	f[now] = max(f[now], SUM - sz[now]);
	if (f[now] < f[root]) root = now;
} 

void get_road(int now, int father, int deg, int &maxdeg) {//把一个点的子树上的点对应的路径找出来
	maxdeg = max(maxdeg, deg);
	num[deg]++;
	for (int i = le[now]; i; i = e[i].nxt)
		if (!in[e[i].to] && e[i].to != father)
			get_road(e[i].to, now, deg + 1, maxdeg);
}

void get_case1(int now) {//第一种情况
	for (int i = le[now]; i; i = e[i].nxt)
		if (e[i].to != fa[now] && !in[e[i].to]) {
			maxdeg = 0;
			get_road(e[i].to, now, 1, maxdeg);
			lstdeg = max(lstdeg, maxdeg);
			for (int j = 1; j <= maxdeg; j++) {
				for (int k = j; k <= maxdeg; k += j)
					sum[j] += num[k];
				ans2[j] += 1ll * sum[j] * Tsum[j];
				
				Tsum[j] += sum[j]; Tnum[j] += num[j];
				sum[j] = num[j] = 0;
			}
		}
}

void get_case2(int now, int father, int cnt) {//第二种情况
	maxdeg = 0;
	for (int i = le[now]; i; i = e[i].nxt)
		if (e[i].to != father && !in[e[i].to] && e[i].to != fa[now]) {
			get_road(e[i].to, now, 1, maxdeg);
			lstdeg = max(lstdeg, maxdeg);
		}
	
	for (int i = 1; i <= maxdeg; i++)
		for (int j = i; j <= maxdeg; j += i)
			sum[i] += num[j];
	
	lim = min(UP, maxdeg);
	for (int i = 1; i <= lim; i++) {//小于你分块的地方就记忆化
		if (rem[i][cnt % i] == -1) {
			rem[i][cnt % i] = 0;
			q.push(i); q.push(cnt % i);
			for (int j = (i - (cnt % i)) % i; j <= lstdeg; j += i)
				rem[i][cnt % i] += Tnum[j];
		}
		ans2[i] += 1ll * rem[i][cnt % i] * sum[i];
	}
	
	for (int i = lim + 1; i <= maxdeg; i++) {//大于的就直接暴力搞
		tmpp = 0;
		for (int j = (i - (cnt % i)) % i; j <= lstdeg; j += i)
			tmpp += Tnum[j];
		ans2[i] += 1ll * tmpp * sum[i];
	}
	
	for (int i = 1; i <= maxdeg; i++)
		num[i] = sum[i] = 0;
}

void slove(int now) {//点分治的递归
	lstdeg = 0;
	in[now] = 1;
	get_case1(now);
	Tnum[0] = 1;
	tmp = 0;
	noww = now;
	while (fa[noww] && !in[fa[noww]]) {//往上跳看要从那个点的其它子树起跳
		get_case2(fa[noww], noww, ++tmp);
		noww = fa[noww];
	}
	
	for (int i = 0; i <= lstdeg; i++)//初始化
		Tnum[i] = Tsum[i] = 0;
	while (!q.empty()) {
		int a = q.front();
		q.pop();
		int b = q.front();
		q.pop();
		rem[a][b] = -1;
	}
	
	for (int i = le[now]; i; i = e[i].nxt)
		if (!in[e[i].to]) {//点分治继续递归
			root = 0;
			SUM = sz[e[i].to];
			get_root(e[i].to, 0);
			slove(root);
		}
}

void work_dfz() {
	root = 0;
	SUM = n;
	get_root(1, 0);
	slove(root);
}

void rongchi() {//容斥得到真正要输出的答案
	for (int i = n - 1; i >= 1; i--)
		ans1[i] += ans1[i + 1];
	
	for (int i = n - 1; i >= 1; i--)
		for (int j = i + i; j <= n - 1; j += i)
			ans2[i] -= ans2[j];
}

int main() {
//	freopen("read.txt", "r", stdin);
//	freopen("write.txt", "w", stdout);
	
	n = read();
	for (int i = 2; i <= n; i++) {
		fa[i] = read();
		add(fa[i], i);
	}
	
	work_line();
	
	f[0] = INF;
	UP = 1000;//分块（反正不能全部记忆化是因为空间不够，那就能开多大开多大呗）
	memset(rem, -1, sizeof(rem));
	work_dfz();
	
	rongchi();
	
	for (int i = 1; i <= n - 1; i++) {
		write(ans1[i] + ans2[i]);
		putchar('
');
	}
	
	return 0;
}