UOJ33【UR #2】树上GCD【长链剖分，阈值法】

给定 \(n\) 个点的以 \(1\) 为根的树，设 \(d(u,v)\) 是 \(u,v\) 间的距离，\(\text{LCA}(u,v)\) 是 \(u,v\) 的最近公共祖先。

对于两个节点 \(u,v\pod{u<v}\)，设 \(a=\text{LCA}(u,v)\)，定义 \(f(u,v)=\gcd(d(u,a),d(v,a))\)。

\(\forall i\in[1,n-1]\) 求有多少个无序对 \((u,v)\) 满足 \(f(u,v)=i\)。

\(n\le 2\cdot 10^5\)。

---

根据套路显然要容斥，转为 \(\forall i\in[1,n-1]\) 计算 \(i\mid d(u,a),d(v,a)\) 的方案数，然后后缀差分即可。

然而 \(u,v\) 是祖先-后代关系时比较麻烦，所以单独考虑这种情况：即为深度数组的桶的后缀和。

虽然这是有根树，但仍然可以考虑点分治。

设 \(rt\) 是当前分治中心，只考虑 \(u\) 到 \(v\) 的路径经过 \(rt\) 的无序对 \((u,v)\)。

• 若 \(u,v\) 都在 \(rt\) 的子树内（以 \(1\) 为根），那么直接对这棵子树 dfs 出深度数组的桶，然后做狄利克雷后缀和。
• 若 \(u\) 在 \(rt\) 的子树外，\(v\) 在 \(rt\) 的子树内，对 \(u\) 的祖先的每个子树分别 dfs 出深度数组的桶，然后查询 \(rt\) 子树的桶，是查询下标\(\bmod x=y\) 之和的形式，这个东西可以根号分治（阈值法），也可以记忆化（本质相同）

分治的单次计算的时间复杂度是 \(O(n\sqrt n)\)，所以总复杂度是 \(O(n\sqrt n)\)。

也有更简单的方法，考虑处理有根树上距离的问题，就想到长剖，合并的时候计算贡献即可。

然而计算贡献时如果直接暴力就跟普通暴力一样了，所以仍然是根号分治：计算复杂度是 \(O(l/i)\)，其中 \(l\) 是长链长度，所以待定阈值 \(B\)，当 \(i\le B\) 时用 dp 方法每次 \(O(n)\) 计算，当 \(i>B\) 时用上述方法，总复杂度是 \(O(nB+n^2\log n/B)\)，取 \(B=O(\sqrt{n\log n})\) 就得到复杂度是 \(O(n\sqrt{n\log n})\)，但因为常数小所以跑起来时间跟点分治方法差不多。

对于 dp 部分，有一种避免容斥的方法：设 \(f_u\) 表示 \(u\) 的子树中到 \(u\) 的距离为 \(i\) 的倍数的点个数，\(g_u\) 表示 \(u\) 的子树中到 \(u\) 的距离为 \(i\) 的倍数 \(-1\) 的点个数，则有 \(f_u=\sum_{v\in\text{child}(u)}g_v+1\)，\(g_i=\sum_{c\in\text{child}^{(i-1)}(u)}f_v\)（\(\text{child}^{(i)}(u)\) 表示 \(u\) 的 \(i\) 级儿子），在 \(g\rightarrow f\) 的时候计算贡献即可。

对于长剖部分，把数组 reverse 之后用 vector 维护是比较简洁的实现方法。

#include<bits/stdc++.h>
#define PB emplace_back
using namespace std;
typedef long long LL;
const int N = 200003, B = 100;
template<typename T>
void rd(T &x){
	int ch = getchar(); x = 0;
	for(;ch < '0' || ch > '9';ch = getchar());
	for(;ch >= '0' && ch <= '9';ch = getchar()) x = x * 10 + ch - '0';
}
int n, fa[N], dep[N], anc[N], f[N], g[N];
LL ans[N], ans2[N];
vector<int> v[N];
int main(){
	rd(n);
	for(int i = 2;i <= n;++ i){
		rd(fa[i]);
		++ans2[dep[i] = dep[fa[i]] + 1];
	}
	for(int i = n-1;i;-- i) ans2[i] += ans2[i+1];
	for(int i = 1;i <= n;++ i) anc[i] = i;
	for(int d = 1;d <= B;++ d){
		memset(f, 0, n+1<<2);
		memset(g, 0, n+1<<2);
		for(int i = n;i > 1;-- i){
			g[anc[i]] += ++f[i]; anc[i] = fa[anc[i]];
			ans[d] += (LL)f[fa[i]] * g[i];
			f[fa[i]] += g[i];
		}
	}
	for(int i = n;i;-- i){
		v[i].PB(1);
		if(v[i].size() > v[fa[i]].size())
			v[i].swap(v[fa[i]]);
		int l1 = v[i].size(), l2 = v[fa[i]].size();
		for(int d = B+1;d <= l1;++ d){
			int r1 = 0, r2 = 0;
			for(int j = d;j <= l1;j += d) r1 += v[i][l1-j];
			for(int j = d;j <= l2;j += d) r2 += v[fa[i]][l2-j];
			ans[d] += (LL)r1 * r2;
		}
		for(int j = 1;j <= l1;++ j)
			v[fa[i]][l2-j] += v[i][l1-j];
		v[i].resize(0);
	}
	for(int i = n>>1;i;-- i)
		for(int j = i<<1;j <= n;j += i)
			ans[i] -= ans[j];
	for(int i = 1;i < n;++ i) printf("%lld\n", ans[i]+ans2[i]);
}