【未知来源】Randomized Binary Search Tree

题意

　　求 \(n\) 个点的 Treap 深度为 \(h=0,1,2,\cdots,n\) 的概率。
　　Treap 是一个随机二叉树，每个节点有权值和优先级，权值和优先级都是 \([0,1]\) 中的随机实数。niubi 的是，由于随机的实数精度足够高，你可以近似认为任意两个权值、任意两个优先级相同的概率是 \(0\)。
　　\(n\le 30000\)

题解

　　又是神题，我他吗都做不来
　　官方题解大概是这样，但本蒟蒻完全没看懂，于是只好向 scb 大佬请教了另一种思考方法（得到的 dp 式子一样）。

　　考虑 Treap 的构造方式：先随机选择 \(n\) 个权值，然后从空树开始，每次插入一个节点。插入操作如下：先随机一个优先级 \(p\)，无视优先级，按照二叉搜索树的方式插入这个节点，然后考虑优先级，一直把这个节点往上旋转，直到满足优先级条件。
　　若我们事先确定了每个权值的优先级，那么把权值按优先级从大到小排序，不难发现 Treap 的加点就变成了每次给树加一个叶子。由于 \(n\) 个数的大小关系不变时 Treap 的形态也不变，我们可以把 \(n\) 个权值离散化成 \(1\) 到 \(n\) 这 \(n\) 个整数，问题是完全等价的。

　　权值序列的每个数都不是确定的，而是在 \([0,1]\) 任取一个实数，为什么可以离散化成 \(1\) 到 \(n\) 这 \(n\) 个整数？如何证明每种 \(1\) 到 \(n\) 的排列对应的原权值序列的数量相同？（不然离散化后算的概率不一样啊）
　　遗憾的是，这个需要微积分等高数知识，过程也比较复杂，本蒟蒻不会简单证法。目前把这当成常识记住就好了。

　　现在优先级已经没用了，我们只需要考虑每次给树加一个权值为 \([1,n]\) 内整数的叶子，这棵树要满足二叉搜索树的性质（即任意点的权值小于其左儿子，大于其右儿子）。求每种树高的出现概率。
　　这就跟普通的求方案数类似，设 \(dp(i,j)\) 表示权值为 \(1,2,\cdots j\) 的点构成深度不大于 \(i\) 的树的概率，则我们枚举根的权值 \(k\)，其左子树的权值集合为 \(\{1,2,\cdots,k-1\}\)，概率就是 \(dp(i-1,k-1)\)；其右子树的权值集合为 \(\{k+1,k+2,\cdots,j\}\)，等价于集合 \(\{1,2,\cdots,j-k\}\)，概率是 \(dp(i-1,j-k)\)。把所有 \(k\) 对应的概率求平均值，就得到了 \(dp(i,j)\)。\[dp(i,j)=\frac{1}{j} \sum\limits_{k=1}^j f(i-1,k-1)\times f(i-1,j-k)\]
　　\(O(n^3)\) 转移可得 \(40\) 分。

　　然后发现由于权值随机，而且 Treap 本身就是在随机优先级时树高维持在 \(\log n\) 级别的数据结构，故期望树高为 \(O(\log n)\)，概率都会集中在这附近。实测大概只需要算到 \(dp(50,)\) 即可满足精度要求。\(O(50n^2)\) 可得 \(50\) 分。

　　然后发现上式显然是个卷积形式，把 \(dp(i)\) 看成生成函数，转移就是 \(f(i)\) 自己卷自己，\(\text{FFT}\) 即可。复杂度 \(O(50n\log n)\)，可得 \(100\) 分。

#include<bits/stdc++.h>
#define N 131075
using namespace std;
inline int read(){
    int x=0; bool f=1; char c=getchar();
    for(;!isdigit(c); c=getchar()) if(c=='-') f=0;
    for(; isdigit(c); c=getchar()) x=(x<<3)+(x<<1)+(c^'0');
    if(f) return x;
    return 0-x;
}
int n;
const double PI = acos(-1);
struct cp{
    double r,i;
    cp(){r=i=0;}
    cp(double _r, double _i){r=_r, i=_i;}
    friend cp operator + (cp a, cp b){return cp(a.r+b.r, a.i+b.i);}
    friend cp operator - (cp a, cp b){return cp(a.r-b.r, a.i-b.i);}
    friend cp operator * (cp a, cp b){return cp(a.r*b.r-a.i*b.i, a.r*b.i+a.i*b.r);}
    friend cp operator / (cp a, double b){return cp(a.r/b, a.i/b);}
}dp[N];
struct Poly{
    int n,bit,r[N];
    void init(int x){
        for(n=1,bit=0; n<x; n<<=1,++bit);
        for(int i=1; i<n; ++i) r[i]=(r[i>>1]>>1)|((i&1)<<(bit-1));
    }
    void dft(cp *a, int f){
        for(int i=0; i<n; ++i) if(i<r[i]) swap(a[i],a[r[i]]);
        cp wn,w,x,y;
        for(int i=1; i<n; i<<=1){
            wn=cp(cos(PI/i),sin(f*PI/i));
            for(int j=0; j<n; j+=i<<1){
                w=cp(1,0);
                for(int k=0; k<i; ++k,w=w*wn)
                    x=a[j+k], y=w*a[j+i+k],
                    a[j+k]=x+y, a[j+i+k]=x-y;
            }
        }
        if(f==-1) for(int i=0; i<n; ++i) a[i]=a[i]/n;
    }
}FFT;
int main(){
    n=read();
    FFT.init(n*2+1);
    dp[0]=cp(1,0); double lst=0;
    for(int scx=1; scx<=min(n,50); ++scx){
        FFT.dft(dp,1);
        for(int i=0; i<FFT.n; ++i) dp[i]=dp[i]*dp[i];
        FFT.dft(dp,-1);
        for(int i=n; i<FFT.n; ++i) dp[i]=cp(0,0);
        for(int i=n; i>0; --i) dp[i]=dp[i-1]/i; dp[0]=cp(1,0);
        printf("%.10lf\n",dp[n].r-lst);
        lst=dp[n].r;
    }
    for(int i=min(n,50)+1; i<=n; ++i) printf("%.10lf\n",0);
    return 0;
}

转载于:https://www.cnblogs.com/scx2015noip-as-php/p/11542993.html