题解：P16100 [ICPC 2019 NAIPC] Heaps of Fun

· · 题解

题意简述

树上每个点 u 独立地在区间 [0,b_u] 中均匀随机取一个实数。求所有父节点取值都小于其孩子取值的概率，并在模 10^9+7 意义下输出。

解题思路

先计算满足大小关系的取值区域体积，再除以整个取值区域的体积 \prod_u b_u。因为各点独立且均匀分布，联合概率密度在这个长方体内恒定，这个体积比就是所求概率。相等的边界体积为 0，因此积分时可以使用非严格不等式。

设 F_u(x) 表示在 u 的子树中，所有点满足堆序，并且 u 的取值至少为 x 时的取值区域体积。固定 u 的取值为 t 后，每棵孩子子树的根都需要至少为 t。不同孩子的变量相互独立，对应体积相乘，所以当 0\le x\le b_u 时：

F_u(x)=\int_x^{b_u}\prod_{v\in\operatorname{child}(u)}F_v(t)\mathrm{d}t

其中 \operatorname{child}(u) 表示 u 的孩子集合。叶子的空乘积为 1，因此 F_u(x)=b_u-x。当 x>b_u 时没有合法取值，F_u(x)=0。

这些函数不是全局多项式：例如叶子函数会在 b_u 处从一次函数变为 0。不过，所有变化位置都来自某个点的上界 b_u。在相邻两个不同上界之间，每个孩子函数都是多项式；相乘并积分后，父节点函数仍然是多项式。按子树归纳，F_u 在每个这样的区间内的次数至多是子树大小。

将所有 b_u 排序，并补上左端点 0。按从右往左的顺序处理每个非空区间 [L,R]。开始处理时，已经知道所有 F_u(R)：最右侧边界处全部为 0，后续区间则接收上一区间算出的左端点值。

在当前区间中改用变量 z=R-x，令 f_u(z)=F_u(R-z)。这样，z 从 0 增加到区间长度 h=R-L，积分的常数项恰好是已经求出的右端点值。若 b_u\ge R，则有：

f_u(z)=F_u(R)+\int_0^z\prod_{v\in\operatorname{child}(u)}f_v(w)\mathrm{d}w

这是将 F_u(R-z)-F_u(R) 对应的积分区间从 [R-z,R] 换元得到的，因此不需要重新积分区间右侧已经处理过的部分。

若 b_u<R，由于相邻断点之间没有任何上界，必有 b_u\le L。整个区间内 F_u 都为 0，直接清空该点的多项式和端点值即可。即使某个点在本区间的函数为 0，仍要独立处理其孩子，保证孩子在以后区间所需的右端点值已经算好。

每个区间按树的后序遍历计算。先卷积所有孩子的系数，设乘积为 g(z)=\sum_{i=0}^d g_i z^i。逐项积分后，父节点系数为：

\begin{aligned} f_{u,0} & =F_u(R) \\ f_{u,i+1} & =\frac{g_i}{i+1} \end{aligned}

将得到的多项式代入 z=h，就得到 F_u(L)。处理到左端点 0 后，根节点的值就是满足堆序的总体积。

代码中，f[u][i] 保存当前区间的多项式系数，val[u] 保存上一区间传来的 F_u(R)，并在计算结束后更新为 F_u(L)。ord 保存后序顺序，确保合并时孩子系数已经属于当前区间；g 和 tmp 依次合并孩子多项式。求值采用从高次项到低次项的乘加形式，避免重复计算幂。

所有系数都可以直接在模数下计算。积分产生的除数至多为 n，各个 b_u 也都小于模数，因此它们均有逆元。预处理 1\sim n 的逆元用于积分，最终再乘 \prod_u b_u 的逆元。

每个区间的多项式合并总计 O(n^2)：不同孩子子树之间的次数乘积，可以对应到最低公共祖先处的点对；常数项带来的开销由各点子树大小之和控制，也不超过 O(n^2)。积分和求值同样受子树大小之和限制。区间数至多为 n，因此总时间复杂度为 O(n^3)，空间复杂度为 O(n^2)。

参考代码

#include <bits/stdc++.h>
using namespace std;

using ll=long long;
const int N=305;
const int mod=1000000007;
int cnt;
int a[N],b[N],ord[N],siz[N],f[N][N],g[N],tmp[N],val[N],inv[N];
vector<int> G[N];
ll Pow(ll x,ll y)
{
    x%=mod;
    ll res=1;
    while(y)
    {
        if(y&1)res=res*x%mod;
        x=x*x%mod;
        y>>=1;
    }
    return res;
}
void dfs(int u)
{
    siz[u]=1;
    for(auto v:G[u]){dfs(v);siz[u]+=siz[v];}
    cnt++;
    ord[cnt]=u;
}
int main()
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);
    int n;
    cin>>n;
    int rt=0;
    ll ans=1;
    for(int i=1;i<=n;i++)
    {
        int p;
        cin>>b[i]>>p;
        ans=ans*b[i]%mod;
        a[i]=b[i];
        if(p)G[p].push_back(i);
        else rt=i;
    }
    dfs(rt);
    for(int i=1;i<=n;i++)inv[i]=i==1?1:1ll*(mod-mod/i)*inv[mod%i]%mod;
    sort(a+1,a+n+1);
    for(int i=n;i>=1;i--)
    {
        if(a[i]==a[i-1])continue;
        int h=a[i]-a[i-1];
        for(int j=1;j<=n;j++)
        {
            int u=ord[j];
            fill(f[u],f[u]+siz[u]+1,0);
            if(b[u]<a[i]){val[u]=0;continue;}
            g[0]=1;
            int len=0;
            for(auto v:G[u])
            {
                fill(tmp,tmp+len+siz[v]+1,0);
                for(int j=0;j<=len;j++)
                {
                    for(int k=0;k<=siz[v];k++)tmp[j+k]=(tmp[j+k]+1ll*g[j]*f[v][k])%mod;
                }
                len+=siz[v];
                copy(tmp,tmp+len+1,g);
            }
            f[u][0]=val[u];
            for(int j=0;j<=len;j++)f[u][j+1]=1ll*g[j]*inv[j+1]%mod;
            val[u]=0;
            for(int j=siz[u];j>=0;j--)val[u]=(1ll*val[u]*h+f[u][j])%mod;
        }
    }
    cout<<val[rt]*Pow(ans,mod-2)%mod<<'\n';
    return 0;
}