Tyvj1953:Normal (点分治+FFT)

来源:互联网 发布:广电网络的宽带怎么样 编辑:程序博客网 时间:2024/06/07 07:37

题目传送门:http://tyvj.joyoi.cn/p/1953


题目分析:好神的一道题,%了一波dalao的题解才会做,然后我发现我对期望一无所知QAQ……

这题的难点其实既不在于点分,也不在于FFT,而是对期望的转化。首先由于期望的线性性,我们可以单独计算每一个点x对期望时间的贡献。而点x每被计算一次,就意味着存在一个分治中心y(x可以等于y),使得x到y路径上的其它点都还没有被作为分治中心。由于x到y路径上的每个点,作为这条路径第一个分治中心的概率是相等的,所以如果设dis(x,y)表示x到y路径上的点数,那么有:

ans=x=1ny=1n1dis(x,y)

又由于dis(x,y)在分母处,所以不能直接用裸的点分治求答案。不妨记num[x]为长度是x的路径条数,然后用点分治+FFT算num数组,最后再统计答案。注意点分治要用容斥的写法,即算完当前连通块的答案再递归子连通块减去非法答案。如果直接将所有子树逐个用FFT卷起来,是没有办法保证时间复杂度的。我不会告诉你tututu用先卷深度小的子树,再卷深度大的子树的方法卡过了


CODE:

#include<iostream>#include<string>#include<cstring>#include<cmath>#include<cstdio>#include<cstdlib>#include<stdio.h>#include<algorithm>using namespace std;const int maxn=300100;const double pi=acos(-1.0);struct Complex{    double X,Y;    Complex (double a=0.0,double b=0.0) : X(a),Y(b) {}} ;Complex operator+(Complex a,Complex b){return Complex(a.X+b.X,a.Y+b.Y);}Complex operator-(Complex a,Complex b){return Complex(a.X-b.X,a.Y-b.Y);}Complex operator*(Complex a,Complex b){return Complex(a.X*b.X-a.Y*b.Y,a.X*b.Y+a.Y*b.X);}int Rev[maxn];Complex A[maxn];struct edge{    int obj;    edge *Next;} e[maxn<<1];edge *head[maxn];int cur=-1;int Size[maxn];int Max[maxn];int now[maxn];int tail;int val[maxn];bool vis[maxn];int num[maxn];int n,N,Lg;void Add(int x,int y){    cur++;    e[cur].obj=y;    e[cur].Next=head[x];    head[x]=e+cur;}void Dfs1(int node,int fa){    now[++tail]=node;    Size[node]=1;    Max[node]=0;    for (edge *p=head[node]; p; p=p->Next)    {        int to=p->obj;        if ( !vis[to] && to!=fa )        {            Dfs1(to,node);            Size[node]+=Size[to];            Max[node]=max(Max[node],Size[to]);        }    }}void Dfs2(int node,int fa,int dep){    val[dep]++;    for (edge *p=head[node]; p; p=p->Next)    {        int to=p->obj;        if ( !vis[to] && to!=fa ) Dfs2(to,node,dep+1);    }}void DFT(Complex *a,double f){    for (int i=0; i<N; i++)        if (i<Rev[i]) swap(a[i],a[ Rev[i] ]);    for (int len=2; len<=N; len<<=1)    {        int mid=len>>1;        double ang=2.0*pi/((double)len);        Complex e( cos(ang) , f*sin(ang) );        for (Complex *p=a; p!=a+N; p+=len)        {            Complex wn(1.0,0.0);            for (int i=0; i<mid; i++)            {                Complex temp=wn*p[mid+i];                p[mid+i]=p[i]-temp;                p[i]=p[i]+temp;                wn=wn*e;            }        }    }}void FFT(int m,int f,int x){    N=1,Lg=0;    while (N<2*m) N<<=1,Lg++;    for (int i=0; i<N; i++)    {        Rev[i]=0;        for (int j=0; j<Lg; j++)            if ( i&(1<<j) ) Rev[i]|=( 1<<(Lg-j-1) );    }    for (int i=0; i<m; i++) A[i]=Complex((double)val[i+1],0.0);    for (int i=m; i<N; i++) A[i]=Complex(0.0,0.0);    DFT(A,1.0);    for (int i=0; i<N; i++) A[i]=A[i]*A[i];    DFT(A,-1.0);    for (int i=0; i<N; i++) A[i].X/=((double)N);    for (int i=0; i<N; i++) num[i+1+x]+=( f*(int)floor( A[i].X+0.5 ) );}void Dec(int node){    tail=0;    Dfs1(node,node);    for (int i=1; i<=tail; i++) val[i]=0;    Dfs2(node,node,1);    FFT(tail,-1,2);}void Solve(int node){    tail=0;    Dfs1(node,node);    if (tail==1)    {        num[1]++;        return;    }    int root=node,sum=Size[node];    for (int i=2; i<=tail; i++)    {        int x=now[i];        if ( max(Max[x],sum-Size[x])<max(Max[root],sum-Size[root]) ) root=x;    }    for (int i=1; i<=sum; i++) val[i]=0;    Dfs2(root,root,1);    FFT(sum,1,0);    vis[root]=true;    for (edge *p=head[root]; p; p=p->Next)    {        int to=p->obj;        if (!vis[to]) Dec(to),Solve(to);    }}int main(){    freopen("1953.in","r",stdin);    freopen("1953.out","w",stdout);    scanf("%d",&n);    for (int i=1; i<=n; i++) head[i]=NULL;    for (int i=1; i<n; i++)    {        int x,y;        scanf("%d%d",&x,&y);        x++,y++;        Add(x,y);        Add(y,x);    }    Solve(1);    double ans=0.0;    for (int i=1; i<=n; i++) ans+=( (double)num[i]/(double)i );    printf("%.4lf\n",ans);    return 0;}//FFT重载运算符的时候要仔细检查!!!
原创粉丝点击