洛谷 7445 「EZEC-7」线段树

官方题解！

首先令 \(v=V+2\)。

考虑一个结点 \([l,r]\)。
则实际上一次操作对它产生影响的概率即这次操作的区间包含它的概率，即 \(p=\frac{l(n-r+1)}{\binom{n+1}2}\)。
设随机加 \(k\) 次 \([-1,V]\) 中的整数，结果为 \(0\) 的方案数为 \(f_k\)，那么枚举被操作次数，一个结点的贡献即为 \[ \sum\limits_{i=0}^m \binom mi p^i (1-p)^{m-i} \left(1-\frac{f_i}{v^i}\right) \]

观察到这其实是一个关于 \(\frac p{1-p}\) 的多项式乘一个 \((1-p)^m\)，设其等于 \[ (1-p)^m \sum\limits_{i=0}^m c_i \left(\frac p{1-p}\right)^i \]

一个简单的做法是直接多点求值，只要常数不爆炸应该都可以过。
不过存在不多点求值的做法。

将有贡献的结点（非叶结点）任意标号为 \(1,2,\dots,n-1\)，令 \(p_i\) 表示第 \(i\) 个结点的以上提到的 \(p\) 值，令 \(a_i = \frac{p_i}{1-p_i},b_i = (1-p_i)^m\)，则总答案为 \[ \begin{aligned} \sum\limits_{i=1}^{n-1} b_i \sum\limits_{j=0}^m c_j a_i^j &= \sum\limits_{j=0}^m c_j \sum\limits_{i=1}^{n-1} b_i a_i^j \end{aligned} \]

于是问题变为对 \(j=0,1,\dots,m\) 求出 \(\sum b_i a_i^j\)。
发现 \[ \sum\limits_{i=1}^{n-1} b_i a_i^j = [x^j] \sum\limits_{i=1}^{n-1} \frac{b_i}{1-a_ix} \]

可以考虑分治 NTT 实现通分并相加的过程，时间复杂度 \(O(n \log^2 n)\)。

接下来，问题只剩下求出 \(f_k\)。
通过组合数最后可以推得一个 \(O\left(\frac {m^2}V\right)\) 的做法，可以参见出题人的题解，这里略去不谈。
通过大力生成函数处理，显然 \[ f_k = [x^0](x^{-1}+1+x+x^2+\dots+x^{V})^k \]

令 \[ F(x) = x^{-1}+1+x+x^2+\dots+x^{V} = \frac{1-x^v}{x(1-x)} \]

然后问题变成 \[ [x^0] F^k(x) \]

事实上这个问题的形式与一个经典问题极其相似，于是这启发我们考虑拉格朗日反演。
但是由于要提取常数项，所以并不能直接运用经典的拉格朗日反演；并且，\(F\) 显然不存在复合逆。

考虑令 \(G = \frac1F\)，令 \(H\) 为 \(G\) 的复合逆，然后考虑另类拉格朗日反演： \[ [x^0] F^k(x) = [x^0] G^{-k}(x) = [x^k] \frac{xH'(x)}{H(x)} \]

至此，如果能求出 \(H(x)\)，问题即迎刃而解。
注意到 \[ G = \frac{x(1-x)}{1-x^v} \]

代入 \(x = H\) 得 \[ xH^v - H^2 + H - x = 0 \]

使用多项式牛顿迭代求解即可。

时间复杂度 \(O(n \log^2 n)\)（视 \(n,m\) 同阶）。

代码：

#include <cstdio>
#include <vector>
#include <cstring>
#include <utility>
#include <algorithm>
#define add(a,b) (a + b >= mod ? a + b - mod : a + b)
#define dec(a,b) (a < b ? a - b + mod : a - b)
using namespace std;
const int N = 1e5;
const int mod = 998244353;
const int LIM = 128;
int n,m,v,inv,lim;
int a[N + 5],b[N + 5],tot;
int c[N + 5];
int ans;
inline int fpow(int a,int b)
{
    int ret = 1;
    for(;b;b >>= 1)
        (b & 1) && (ret = (long long)ret * a % mod),a = (long long)a * a % mod;
    return ret;
}
namespace Poly
{
    const int LG = 19;
    const int N = 1 << LG + 1;
    const int G = 3;
    int lg2[N + 5];
    int rev[N + 5],fac[N + 5],ifac[N + 5],inv[N + 5];
    int rt[N + 5],irt[N + 5];
    inline void init()
    {
        for(register int i = 2;i <= N;++i)
            lg2[i] = lg2[i >> 1] + 1;
        int w = fpow(G,(mod - 1) / N / 2),wi = fpow(w,mod - 2);
        rt[0] = irt[0] = 1;
        for(register int i = 1;i <= N;++i)
            rt[i] = (long long)rt[i - 1] * w % mod,
            irt[i] = (long long)irt[i - 1] * wi % mod;
        for(register int i = 0;i < N;++i)
            rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << LG),
            i < rev[i] && (swap(rt[i],rt[rev[i]]),swap(irt[i],irt[rev[i]]),1);
        fac[0] = 1;
        for(register int i = 1;i <= N;++i)
            fac[i] = (long long)fac[i - 1] * i % mod;
        ifac[N] = fpow(fac[N],mod - 2);
        for(register int i = N;i;--i)
            ifac[i - 1] = (long long)ifac[i] * i % mod;
        for(register int i = 1;i <= N;++i)
            inv[i] = (long long)ifac[i] * fac[i - 1] % mod;
    }
    struct poly
    {
        vector<int> a;
        inline poly(int x = 0)
        {
            x && (a.push_back(x),1);
        }
        inline poly(const vector<int> &o)
        {
            a = o,shrink();
        }
        inline poly(const poly &o)
        {
            a = o.a,shrink();
        }
        inline void shrink()
        {
            for(;!a.empty() && !a.back();a.pop_back());
        }
        inline int size() const
        {
            return a.size();
        }
        inline void resize(int x)
        {
            a.resize(x);
        }
        inline int operator[](int x) const
        {
            if(x < 0 || x >= size())
                return 0;
            return a[x];
        }
        inline void clear()
        {
            vector<int>().swap(a);
        }
        inline poly rever() const
        {
            return poly(vector<int>(a.rbegin(),a.rend()));
        }
        inline void dit()
        {
            int n = size();
            for(register int i = 1,l = n >> 1;i < n;i <<= 1,l >>= 1)
                for(register int j = 0,w,o = 0;j < i;++j,o += l << 1)
                {
                    w = rt[i | j];
                    for(register int k = o,t;k < o + l;++k)
                        t = (long long)a[k + l] * w % mod,
                        a[k + l] = dec(a[k],t),
                        a[k] = add(a[k],t);
                }
        }
        inline void dif()
        {
            int n = size();
            for(register int i = n >> 1,l = 1;i;i >>= 1,l <<= 1)
                for(register int j = 0,w,o = 0;j < i;++j,o += l << 1)
                {
                    w = irt[i | j];
                    for(int k = o,t;k < o + l;++k)
                        t = a[k + l],
                        a[k + l] = dec(a[k],t),a[k + l] = (long long)a[k + l] * w % mod,
                        a[k] = add(a[k],t);
                }
            for(register int i = 0;i < n;++i)
                a[i] = (long long)a[i] * inv[n] % mod;
        }
        inline void ntt(int type = 1)
        {
            type == 1 ? dit() : dif();
        }
        friend inline poly operator+(const poly &a,const poly &b)
        {
            vector<int> ret(max(a.size(),b.size()));
            for(register int i = 0;i < ret.size();++i)
                ret[i] = add(a[i],b[i]);
            return poly(ret);
        }
        friend inline poly operator-(const poly &a,const poly &b)
        {
            vector<int> ret(max(a.size(),b.size()));
            for(register int i = 0;i < ret.size();++i)
                ret[i] = dec(a[i],b[i]);
            return poly(ret);
        }
        friend inline poly operator*(poly a,poly b)
        {
            if(a.a.empty() || b.a.empty())
                return poly();
            int n = a.size(),m = b.size();
            int lim = 1,tot = n + m - 1;
            for(;lim < tot;lim <<= 1);
            a.resize(lim),b.resize(lim);
            for(register int i = n;i < lim;++i)
                a.a[i] = 0;
            for(register int i = m;i < lim;++i)
                b.a[i] = 0;
            a.ntt(),b.ntt();
            for(register int i = 0;i < lim;++i)
                a.a[i] = (long long)a[i] * b[i] % mod;
            a.ntt(-1),a.shrink();
            return a;
        }
        poly &operator+=(const poly &o)
        {
            resize(max(size(),o.size()));
            for(register int i = 0;i < o.size();++i)
                a[i] = add(a[i],o[i]);
            return *this;
        }
        poly &operator-=(const poly &o)
        {
            resize(max(size(),o.size()));
            for(register int i = 0;i < o.size();++i)
                a[i] = dec(a[i],o[i]);
            return *this;
        }
        poly &operator*=(poly o)
        {
            return (*this) = (*this) * o;
        }
        poly deriv() const
        {
            if(a.empty())
                return poly();
            vector<int> ret(size() - 1);
            for(register int i = 0;i < size() - 1;++i)
                ret[i] = (long long)(i + 1) * a[i + 1] % mod;
            return poly(ret);
        }
        poly integ() const
        {
            if(a.empty())
                return poly();
            vector<int> ret(size() + 1);
            for(register int i = 0;i < size();++i)
                ret[i + 1] = (long long)a[i] * inv[i + 1] % mod;
            return poly(ret);
        }
        inline poly modxn(int n) const
        {
            if(a.empty())
                return poly();
            n = min(n,size());
            return poly(vector<int>(a.begin(),a.begin() + n));
        }
        inline poly inver(int m) const
        {
            poly ret(fpow(a[0],mod - 2)),f,g;
            for(register int k = 1;k < m;)
            {
                k <<= 1,f.resize(k),g.resize(k);
                for(register int i = 0;i < k;++i)
                    f.a[i] = (*this)[i],g.a[i] = ret[i];
                f.ntt(),g.ntt();
                for(register int i = 0;i < k;++i)
                    f.a[i] = (long long)f[i] * g[i] % mod;
                f.ntt(-1);
                for(register int i = 0;i < (k >> 1);++i)
                    f.a[i] = 0;
                f.ntt();
                for(register int i = 0;i < k;++i)
                    f.a[i] = (long long)f[i] * g[i] % mod;
                f.ntt(-1);
                ret.resize(k);
                for(register int i = (k >> 1);i < k;++i)
                    ret.a[i] = dec(0,f[i]);
            }
            return ret.modxn(m);
        }
        inline pair<poly,poly> div(poly o) const
        {
            if(size() < o.size())
                return make_pair(poly(),*this);
            poly f,g;
            f = (rever().modxn(size() - o.size() + 1) * o.rever().inver(size() - o.size() + 1)).modxn(size() - o.size() + 1).rever();
            g = (modxn(o.size() - 1) - o.modxn(o.size() - 1) * f.modxn(o.size() - 1)).modxn(o.size() - 1);
            return make_pair(f,g);
        }
        inline poly log(int m) const
        {
            return (deriv() * inver(m)).integ().modxn(m);
        }
        inline poly exp(int m) const
        {
            poly ret(1);
            for(register int k = 1;k < m;)
                k <<= 1,ret = (ret * (1 - ret.log(k) + modxn(k))).modxn(k);
            return ret.modxn(m);
        }
        inline poly pow(int m,int k1,int k2 = -1) const
        {
            if(a.empty())
                return poly();
            if(k2 == -1)
                k2 = k1;
            int t = 0;
            for(;t < size() && !a[t];++t);
            if((long long)t * k1 >= m)
                return poly();
            poly ret;
            ret.resize(m);
            int u = fpow(a[t],mod - 2),v = fpow(a[t],k2);
            for(register int i = 0;i < m - t * k1;++i)
                ret.a[i] = (long long)operator[](i + t) * u % mod;
            ret = ret.log(m - t * k1);
            for(register int i = 0;i < ret.size();++i)
                ret.a[i] = (long long)ret[i] * k1 % mod;
            ret = ret.exp(m - t * k1),t *= k1,ret.resize(m);
            for(register int i = m - 1;i >= t;--i)
                ret.a[i] = (long long)ret[i - t] * v % mod;
            for(register int i = 0;i < t;++i)
                ret.a[i] = 0;
            return ret;
        }
    };
}
using Poly::init;
using Poly::poly;
inline int C(int n,int m)
{
    return n < m ? 0 : (long long)Poly::fac[n] * Poly::ifac[m] % mod * Poly::ifac[n - m] % mod;
}
poly f,g;
pair<poly,poly> res;
inline poly solve(int m)
{
    poly ret,f2,fv1,fv,t,t1;
    for(register int k = 1;k < m;)
    {
        k <<= 1,f2.resize(k << 1);
        for(register int i = 0;i < (k << 1);++i)
            f2.a[i] = ret[i];
        f2.ntt();
        fv1 = ret.pow(k,v - 1),fv = fv1,fv.resize(k << 1),fv.ntt();
        for(register int i = 0;i < (k << 1);++i)
            fv.a[i] = (long long)fv[i] * f2[i] % mod,
            f2.a[i] = (long long)f2[i] * f2[i] % mod;
        f2.ntt(-1),f2.resize(k),fv.ntt(-1),fv.resize(k);
        t.resize(k),t1.resize(k);
        for(register int i = 0;i < k;++i)
            t.a[i] = add(fv[i - 1],ret[i]),
            t.a[i] = dec(t[i],f2[i]);
        t.a[1] = dec(t[1],1);
        for(register int i = 0;i < k;++i)
            t1.a[i] = ((long long)v * fv1[i - 1] % mod - 2 * ret[i] + 2 * mod) % mod;
        t1.a[0] = add(t1[0],1);
        ret -= (t * t1.inver(k)).modxn(k);
    }
    return ret.modxn(m);
}
inline pair<poly,poly> operator+(pair<poly,poly> a,pair<poly,poly> b)
{
    int alen = a.first.size(),blen = b.first.size();
    int n = a.second.size(),m = b.second.size();
    int lim = 1,tot = n + m - 1;
    for(;lim < tot;lim <<= 1);
    a.first.resize(lim),a.second.resize(lim),b.first.resize(lim),b.second.resize(lim);
    a.first.ntt(),a.second.ntt(),b.first.ntt(),b.second.ntt();
    for(register int i = 0;i < lim;++i)
        a.first.a[i] = ((long long)a.first[i] * b.second[i] + (long long)a.second[i] * b.first[i]) % mod,
        a.second.a[i] = (long long)a.second[i] * b.second[i] % mod;
    a.first.ntt(-1),a.second.ntt(-1),
    a.first.shrink(),a.second.shrink();
    return a;
}
inline pair<poly,poly> solve(int l,int r)
{
    if(l == r)
    {
        poly t;
        t.resize(2),t.a[0] = 1,t.a[1] = dec(0,a[l]);
        return make_pair(poly(b[l]),t);
    }
    int mid = l + r >> 1;
    return solve(l,mid) + solve(mid + 1,r);
}
void build(int l,int r)
{
    if(l == r)
        return ;
    ++tot;
    int p = (long long)l * (n - r + 1) % mod * inv % mod;
    a[tot] = (long long)p * fpow((1 - p + mod) % mod,mod - 2) % mod,
    b[tot] = fpow((1 - p + mod) % mod,m);
    int mid = l + r >> 1;
    build(l,mid),build(mid + 1,r);
}
int main()
{
    init();
    scanf("%d%d%d",&n,&m,&v),v += 2,lim = LIM,inv = fpow(C(n + 1,2),mod - 2);
    if(n == 1)
    {
        puts("0");
        return 0;
    }
    if(v <= lim)
    {
        f = solve(m + 2),f.resize(m + 2),g = f.deriv();
        for(register int i = 0;i <= m;++i)
            f.a[i] = f[i + 1];
        f.resize(m + 1),f = (g * f.inver(m + 1)).modxn(m + 1);
    }
    else
    {
        f.resize(m + 1),f.a[0] = 1;
        for(register int i = 1;i <= m;++i)
            for(register int j = 0;j * v <= i;++j)
                f.a[i] = (f[i] + (long long)C(i,j) * (j & 1 ? mod - 1 : 1) % mod * C(2 * i - j * v - 1,i - 1)) % mod;
    }
    for(register int i = 0;i <= m;++i)
        c[i] = (long long)C(m,i) * (fpow(v,i) - f[i] + mod) % mod * fpow(v,(mod - i - 1) % (mod - 1)) % mod;
    build(1,n),res = solve(1,tot),res.first = (res.first.modxn(m + 1) * res.second.inver(m + 1)).modxn(m + 1);
    for(register int i = 0;i <= m;++i)
        ans = (ans + (long long)c[i] * res.first[i]) % mod;
    printf("%d\n",ans);
}