bzoj 3451 Tyvj1953 Normal 点分治+FFT

Description

---

对一棵树做随机点分治（分支中心随机的点分治），定义一次分治时间复杂度为树的size，求期望时间复杂度

Solution

---

根据期望的线性性我们知道答案就是每个点被选到概率之和
考虑一个点x成为y的点分树的祖先的概率，就是它在x到y这条链中第一次就被选到的概率，也就是$\frac{1}{dis(x,y)}$
那么答案就是$\sum _{i=1}^n\sum _{j=1}^n\frac{1}{dis(i,j)}$，这里的链是有顺序的
记$w[x]$表示长度为x的边有多少条，那么答案就是${\sum }_{i=1}^n\frac{w[x]}{x}$
求不同边的数量可以点分治然后FFT优化合并的复杂度

注意点分治的不同写法，可以先算当前子树的答案再减去各个子树内部重叠的答案，这样比较好写
当然也可以枚举子树按照深度来慢慢合并，当然时间和代码也得慢慢了。

没有权限号，没有手机，交不了就当精神AC了把

Code

---

#include <stdio.h>
#include <string.h>
#include <algorithm>
#include <math.h>
#define rep(i,st,ed) for (int i=st;i<=ed;++i)

typedef long long LL;
const double pi=acos(-1);
const int INF=0x3f3f3f3f;
const int MOD=998244353;
const int N=2000005;

struct edge {int x,y,next;} e[N*2];
struct com {double x,y;} b[N];

com operator + (com a,com b) {return (com) {a.x+b.x,a.y+b.y};}
com operator - (com a,com b) {return (com) {a.x-b.x,a.y-b.y};}
com operator * (com a,com b) {return (com) {a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};}
com operator / (com a,double b) {return (com) {a.x/b,a.y/b};}

int rv[N],size[N],mx[N];
int ls[N],edCnt,rt,sum;
int v[N];

double w[N];

bool vis[N];

int read() {
	int x=0,v=1; char ch=getchar();
	for (;ch<'0'||ch>'9';v=(ch=='-')?(-1):(v),ch=getchar());
	for (;ch<='9'&&ch>='0';x=x*10+ch-'0',ch=getchar());
	return x*v;
}

void add_edge(int x,int y) {
	e[++edCnt]=(edge) {x,y,ls[x]}; ls[x]=edCnt;
	e[++edCnt]=(edge) {y,x,ls[y]}; ls[y]=edCnt;
}

int ksm(int x,int dep) {
	int res=1;
	for (;dep;dep>>=1) {
		(dep&1)?(res=1LL*res*x%MOD):0;
		x=1LL*x*x%MOD;
	}
	return res;
}

void FFT(com *a,int n,int f) {
	for (int i=0;i<n;++i) if (i<rv[i]) std:: swap(a[i],a[rv[i]]);
	for (int i=1;i<n;i<<=1) {
		com wn=(com) {cos(pi/i),f*sin(pi/i)};
		for (int j=0;j<n;j+=(i<<1)) {
			com w=(com) {1,0};
			for (int k=0;k<i;++k) {
				com u=a[j+k],v=a[j+k+i]*w;
				a[j+k]=u+v,a[j+k+i]=u-v;
				w=w*wn;
			}
		}
	}
	if (f==-1) for (int i=0;i<n;++i) a[i]=a[i]/n;
}

void get_size(int x,int from) {
	size[x]=1; mx[x]=0;
	for (int i=ls[x];i;i=e[i].next) {
		if (e[i].y==from||vis[e[i].y]) continue;
		get_size(e[i].y,x); size[x]+=size[e[i].y];
		mx[x]=std:: max(mx[x],size[e[i].y]);
	}
	mx[x]=std:: max(mx[x],sum-size[x]);
	if (mx[x]<mx[rt]) rt=x;
}

void get_dep(int x,int from,int d) {
	v[++v[0]]=d;
	for (int i=ls[x];i;i=e[i].next) {
		if (e[i].y==from||vis[e[i].y]) continue;
		get_dep(e[i].y,x,d+1);
	}
}

void calc(double xs) {
	int n=0; rep(i,1,v[0]) n=std:: max(n,v[i]);
	int lim,lg; for (lim=1,lg=0;lim<=n*2;) lim<<=1,lg++;
	for (int i=0;i<lim;++i) rv[i]=(rv[i>>1]>>1)|((i&1)<<(lg-1));
	rep(i,1,v[0]) b[v[i]].x++;
	FFT(b,lim,1);
	for (int i=0;i<lim;++i) b[i]=b[i]*b[i];
	FFT(b,lim,-1);
	for (int i=0;i<lim;++i) w[i]=w[i]+b[i].x*xs;
	for (int i=0;i<lim;++i) b[i]=(com) {0,0};
	v[0]=0;
}

void dfs(int x,int from) {
	size[x]=1;
	for (int i=ls[x];i;i=e[i].next) {
		if (e[i].y==from||vis[e[i].y]) continue;
		dfs(e[i].y,x); size[x]+=size[e[i].y];
	}
}

void solve(int x) {
	vis[x]=1;
	get_dep(x,0,0);
	calc(1);
	for (int i=ls[x];i;i=e[i].next) {
		if (vis[e[i].y]) continue;
		// v[v[0]=1]=0;
		get_dep(e[i].y,0,1);
		calc(-1);
		// w[0]++;
	}
	dfs(x,0);
	for (int i=ls[x];i;i=e[i].next) {
		if (vis[e[i].y]) continue;
		rt=0; sum=size[e[i].y];
		get_size(e[i].y,x);
		solve(rt);
	}
}

int main(void) {
	freopen("data.in","r",stdin);
	int n=read();
	rep(i,2,n) add_edge(read()+1,read()+1);
	mx[0]=INF; rt=0; sum=n;
	get_size(1,0);
	vis[rt]=1;
	solve(rt);
	double ans=0;
	for (int i=0;i<=2*n;++i) ans=(ans+w[i]/(i+1));
	printf("%.6lf\n", ans);
	return 0;
}