求随机点分治的期望复杂度,每次对一颗大小为\(n\)的子树需要\(O(n)\)的复杂度。
考虑计算每个点期望下被算的次数,根据期望的线性性,最后将每个点的答案加起来就可以了。
计算点u的计算次数可以考虑v对点u的贡献,即在v作为分治重心的时候u在v所在的子树里面。
不难发现如果v对u产生了贡献,那么从u到v的路径上,v必定是第一个选的,路径外的点怎么选没有影响,于是期望贡献为\(\frac{1}{dis(u,v)+1}\)。
答案即\(\sum_{i=1}^{n}\sum_{j=1}^{n}\frac{1}{dis(i,j)+1}\),又转化成了树上路径问题,考虑点分治,计算出每颗子树内的所有路径长度的出现次数,直接FFT优化即可。
/*=======================================
* Author : ylsoi
* Time : 2019.2.13
* Problem : bzoj3451
* E-mail : ylsoi@foxmail.com
* ====================================*/
#include<bits/stdc++.h>
#define REP(i,a,b) for(int i=a,i##_end_=b;i<=i##_end_;++i)
#define DREP(i,a,b) for(int i=a,i##_end_=b;i>=i##_end_;--i)
#define debug(x) cout<<#x<<"="<<x<<" "
#define fi first
#define se second
#define mk make_pair
#define pb push_back
typedef long long ll;
using namespace std;
void File(){
freopen("bzoj3451.in","r",stdin);
freopen("bzoj3451.out","w",stdout);
}
template<typename T>void read(T &_){
_=0; T f=1; char c=getchar();
for(;!isdigit(c);c=getchar())if(c=='-')f=-1;
for(;isdigit(c);c=getchar())_=(_<<1)+(_<<3)+(c^'0');
_*=f;
}
const int maxn=3e4+10;
const int inf=0x3f3f3f3f;
const double pi=acos(-1);
int n;
int beg[maxn],las[maxn<<1],to[maxn<<1],cnte=1;
double ans;
void add(int u,int v){
las[++cnte]=beg[u],beg[u]=cnte,to[cnte]=v;
las[++cnte]=beg[v],beg[v]=cnte,to[cnte]=u;
}
struct cp{
double x,y;
cp(double xx=0,double yy=0){
x=xx,y=yy;
}
};
cp operator + (cp a,cp b){return cp(a.x+b.x,a.y+b.y);}
cp operator - (cp a,cp b){return cp(a.x-b.x,a.y-b.y);}
cp operator * (cp a,cp b){return cp(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
cp operator / (cp a,double b){return cp(a.x/b,a.y/b);}
int lim,cnt,dn[maxn<<2];
cp g[maxn<<2],ig[maxn<<2];
void fft(cp *A,int ty){
REP(i,0,lim-1)if(i<dn[i])swap(A[i],A[dn[i]]);
for(int len=1;len<lim;len<<=1){
cp w= ty==1 ? g[len<<1] : ig[len<<1];
for(int L=0;L<lim;L+=len<<1){
cp wk=cp(1,0);
REP(i,L,L+len-1){
cp u=A[i],v=A[i+len]*wk;
A[i]=u+v;
A[i+len]=u-v;
wk=wk*w;
}
}
}
if(ty==-1)
REP(i,0,lim-1)A[i]=A[i]/lim;
}
int sz[maxn],tot_sz,Min_sz,rt;
bool vis[maxn];
void findrt(int u,int fh){
int Max_sz=0;
sz[u]=1;
for(int i=beg[u];i;i=las[i]){
int v=to[i];
if(v==fh || vis[v])continue;
findrt(v,u);
sz[u]+=sz[v];
Max_sz=max(Max_sz,sz[v]);
}
Max_sz=max(Max_sz,tot_sz-sz[u]);
if(Max_sz<Min_sz){
Min_sz=Max_sz;
rt=u;
}
}
int dis[maxn],cnt_dis;
void get_dis(int u,int fh,int d){
dis[++cnt_dis]=d;
for(int i=beg[u];i;i=las[i]){
int v=to[i];
if(v==fh || vis[v])continue;
get_dis(v,u,d+1);
}
}
cp a[maxn<<2];
void solve(int u,int s,int ty){
cnt_dis=0;
get_dis(u,0,0);
int max_dis=0;
REP(i,1,cnt_dis)max_dis=max(max_dis,dis[i]);
lim=1,cnt=0;
while(lim<=max_dis*2)lim<<=1,++cnt;
if(!cnt)cnt=1;
REP(i,0,lim-1){
dn[i]=dn[i>>1]>>1|((i&1)<<(cnt-1));
a[i]=cp(0,0);
}
REP(i,1,cnt_dis)a[dis[i]].x+=1;
fft(a,1);
REP(i,0,lim-1)a[i]=a[i]*a[i];
fft(a,-1);
REP(i,0,lim-1)if(i+s>0)
ans=ans+a[i].x/(i+s)*ty;
}
void divide(int u){
vis[u]=1; findrt(u,0);
solve(u,1,1);
for(int i=beg[u];i;i=las[i]){
int v=to[i];
if(vis[v])continue;
solve(v,3,-1);
tot_sz=sz[v],Min_sz=inf;
findrt(v,0);
divide(rt);
}
}
int main(){
File();
read(n);
int u,v;
REP(i,1,n-1){
read(u),read(v);
add(u+1,v+1);
}
lim=1;
while(lim<=n+n)lim<<=1;
g[lim]=cp(cos(pi*2.0/lim),sin(pi*2.0/lim));
ig[lim]=cp(cos(pi*2.0/lim),-sin(pi*2.0/lim));
for(int i=lim>>1;i;i>>=1){
g[i]=g[i<<1]*g[i<<1];
ig[i]=ig[i<<1]*ig[i<<1];
}
tot_sz=n,Min_sz=inf;
findrt(1,0);
divide(rt);
printf("%.4lf\n",round(ans*1e4)/1e4);
return 0;
}
[bzoj3451]Tyvj1953 Normal——点分治+fft
原文:https://www.cnblogs.com/ylsoi/p/10369642.html