题解 【P7445 「EZEC-7」线段树】

· · 题解

出题人题解。

为了方便把 V2

一个节点 [l,r]\mathrm{Pushdown} 当且仅当所有覆盖它的修改加起来不等于 0

一次操作它被覆盖的期望为 p=\dfrac{l\times(n-r+1)}{\frac{n(n+1)}{2}}

设覆盖 k 次,加起来为 0 的方案数为 f(k)

f(k)=[x^0](x^{-1}+1+x+x^2+\cdots+x^{V-2})^k

一个节点对于答案的贡献就是

S=\sum_{i=0}^{m}p^i(1-p)^{m-i}\binom{m}{i}\dfrac{V^i-f(i)}{V^i}

\dfrac{p}{1-p} 看做参数,可以得到一个多项式。

发现总共只有 n-1 个点(线段树非叶节点个数),多点求值就好了,非常简便。

当然也可以不用这个科技,验题人给了个常数较小的分治FFT做法,稍微推推式子就好了。

现在唯一剩下的问题是计算 k\in [0,m]f(k)

方法一

f(k)=[x^k](1+x+x^2+x^3+\cdots+x^{V-1})^k =[x^k](\frac{1-x^{V}}{1-x})^k (1-x^V)^k=\sum_{i=0}\binom{k}{i}x^{iV}(-1)^i (1-x)^{-k}=\sum_{i=0}\binom{k+i-1}{i}x^{i} f(k)=\sum_{i=0}^{\frac{k}{V}}\binom{k}{i}(-1)^i\binom{2k-Vi-1}{k-1}

复杂度 O(\dfrac{m^2}{V})

方法二

F(x)=1+x+x^2+\cdots+x^{V-1} G(x)=\dfrac{x}{F(x)}

我们要求的是 [x^0]\left(\dfrac{F(x)}{x}\right)^k,即 [x^0]G^{-k}(x)

可以发现 G(x) 没有常数项且一次项系数不为 0 。令 G(x) 的复合逆为 H(x)

根据EI鸽鸽在员交课件里的《另类拉格朗日反演》一节里的一个式子(感谢 EI 的指导。证明可以看验题人题解,大概是普通拉反证明第一步两边 k 次幂):

G,H 互为复合逆,那么

[x^n]G^k=[x^{-k-1}]H'H^{-n-1}

带进去就得到

[x^0]G^{-k}(x)=[x^k]\dfrac{H'(x)}{\frac{H(x)}{x}}

只要算出 H(x) 就好了。

注意到

G(x)=\dfrac{x}{F(x)}=\dfrac{x-x^2}{1-x^V}

相当于要求一个 H(x) 满足

\dfrac{H(x)-H^2(x)}{1-H^{V}(x)}\equiv x

化简一下:

H(x)-H^2(x)-x+xH^V(x)\equiv 0

牛顿迭代即可,复杂度 O(n\log n)n,m 视为同阶)。

综上,使用方法二,复杂度为 O(n\log^2 n)

思维难度应该不是很大,有点板子,放 E 应该差不多吧。

upd:上面那句话是赛前写的,不想删掉/wq。但是原本确实是为了放一个简单一点的 E 才选了这题的,没想到场上最高 30。。。

所有点的时限至少为验题人代码 3 倍,std 2 倍多,并且我的板子不算很快,还有极大的优化空间,我用我写过最慢的板子都能过,所以应该不会有人被卡常了。

#include<bits/stdc++.h>
using namespace std;
#define fi first
#define se second
#define mkp(x,y) make_pair(x,y)
#define pb(x) push_back(x)
#define sz(v) (int)v.size()
typedef long long LL;
typedef double db;
template<class T>bool ckmax(T&x,T y){return x<y?x=y,1:0;}
template<class T>bool ckmin(T&x,T y){return x>y?x=y,1:0;}
#define rep(i,x,y) for(int i=x,i##end=y;i<=i##end;++i)
#define per(i,x,y) for(int i=x,i##end=y;i>=i##end;--i)
inline int read(){
    int x=0,f=1;char ch=getchar();
    while(!isdigit(ch)){if(ch=='-')f=0;ch=getchar();}
    while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    return f?x:-x;
}
const int N=200005;
const int M=N<<2;
int n,m,V,f[M],a[N<<1],tot,IV,c[N],b[N<<1],ans,P[N];
#define mod 998244353
namespace math{
int inv[N],fac[N],ifc[N];
inline int qpow(int n,int k){int res=1;for(;k;k>>=1,n=1ll*n*n%mod)if(k&1)res=1ll*n*res%mod;return res;}
inline void fmod(int&x){x-=mod,x+=x>>31&mod;}
inline int comb(int n,int m){return n<m?0:1ll*fac[n]*ifc[m]%mod*ifc[n-m]%mod;}
void initmath(const int&n=N-1){
    inv[1]=1;for(int i=2;i<=n;++i)inv[i]=1ll*inv[mod%i]*(mod-mod/i)%mod;
    fac[0]=1;for(int i=1;i<=n;++i)fac[i]=1ll*i*fac[i-1]%mod;
    ifc[n]=qpow(fac[n],mod-2);for(int i=n-1;i>=0;--i)ifc[i]=1ll*(i+1)*ifc[i+1]%mod;
}
}
using namespace math;

namespace poly{

int rev[M],lg,lim,rt[M];
void initpoly(const int&n){
    for(lg=0,lim=1;lim<=n;++lg,lim<<=1);
    for(int i=0;i<lim;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1));
    for(int i=1;i<lim;i<<=1){
        const int w=qpow(3,(mod-1)/(i<<1));
        rt[i]=1;
        for(int j=1;j<i;++j)rt[i+j]=1ll*rt[i+j-1]*w%mod;
    }
}
void NTT(int*a,int op){
    if(!op)reverse(a+1,a+lim);
    for(int i=0;i<lim;++i)
        if(i>rev[i])swap(a[i],a[rev[i]]);
    for(int i=1;i<lim;i<<=1){
        for(int j=0;j<lim;j+=i<<1){
            for(int k=0;k<i;++k){
                const int X=a[j+k],Y=1ll*rt[i+k]*a[i+j+k]%mod;
                fmod(a[j+k]=X+Y),fmod(a[i+j+k]=X-Y+mod);
            }
        }
    }
    if(op)return;const int ilim=qpow(lim,mod-2);
    for(int i=0;i<lim;++i)a[i]=1ll*a[i]*ilim%mod;
}
#define clr(a,n) memset(a,0,sizeof(int)*(n))
#define cpy(a,b,n) memcpy(a,b,sizeof(int)*(n))
void poly_mul(int*f,int*g,int*ans,int n,int m){
    static int A[M],B[M];initpoly(n+m);
    cpy(A,f,n),clr(A+n,lim-n),NTT(A,1);
    cpy(B,g,m),clr(B+m,lim-m),NTT(B,1);
    for(int i=0;i<lim;++i)A[i]=1ll*A[i]*B[i]%mod;
    NTT(A,0),cpy(ans,A,n+m-1);
}
void dao(int*g,int*f,int n){
    for(int i=0;i<n-1;++i)g[i]=1ll*(i+1)*f[i+1]%mod;g[n-1]=0;
}
void poly_inv(int*g,int*f,int n){
    static int A[M],B[M];
    if(n==1)return g[0]=qpow(f[0],mod-2),void();
    int m=(n+1)>>1;
    poly_inv(g,f,m),initpoly(n<<1);
    cpy(A,f,n),clr(A+n,lim-n),cpy(B,g,m),clr(B+m,lim-m);
    NTT(A,1),NTT(B,1);
    for(int i=0;i<lim;++i)A[i]=1ll*B[i]*(2-1ll*A[i]*B[i]%mod+mod)%mod;
    NTT(A,0),cpy(g,A,n);
}

#define lc (p<<1)
#define rc (p<<1|1)
int pool[N<<6],*mem=pool,*ev[M];
void eva_mul(int*f,int*g,int*ans,int n,int m){
    static int A[M],B[M];initpoly(n);
    cpy(A,f,n),clr(A+n,lim-n),NTT(A,1);
    cpy(B,g,m),clr(B+m,lim-m),NTT(B,1);
    for(int i=0;i<lim;++i)A[i]=1ll*A[i]*B[i]%mod;
    NTT(A,0),cpy(ans,A+m-1,n-m+1);
}
void eva_solve1(int l,int r,int p,int*a){
    ev[p]=mem,mem+=r-l+1;
    if(r-l==1)return ev[p][0]=mod-a[l],ev[p][1]=1,void();
    int mid=(l+r)>>1;
    eva_solve1(l,mid,lc,a);
    eva_solve1(mid,r,rc,a);
    poly_mul(ev[lc],ev[rc],ev[p],mid-l+1,r-mid+1);
}
void eva_solve2(int l,int r,int p,int*h,int*f){
    if(r-l==1)return h[l]=f[0],void();
    int mid=(l+r)>>1,*al,*ar;
    al=mem,mem+=mid-l,ar=mem,mem+=r-mid;
    eva_mul(f,ev[rc],al,r-l,r-mid+1);
    eva_mul(f,ev[lc],ar,r-l,mid-l+1);
    eva_solve2(l,mid,lc,h,al);
    eva_solve2(mid,r,rc,h,ar);
}
void poly_eva(int*g,int*f,int*a,int n,int m){
    static int A[M];
    n=max(n,m);
    eva_solve1(0,n,1,a);
    reverse(ev[1],ev[1]+n+1);
    poly_inv(A,ev[1],n),reverse(A,A+n);
    poly_mul(A,f,A,n,n),cpy(A,A+n,n);
    eva_solve2(0,n,1,g,A);
    for(int i=0;i<m;++i)fmod(g[i]=1ll*g[i]*a[i]%mod+f[0]);
    for(int i=m;i<n;++i)g[i]=0;
}
void poly_sqa(int*g,int*f,int n){
    initpoly(n<<1);
    NTT(f,1);
    for(int i=0;i<lim;++i)g[i]=1ll*f[i]*f[i]%mod;
    NTT(g,0);
}
void poly_qpow(int*g,int *f,int k,int n){
    static int A[M];
    cpy(A,f,n),g[0]=1;
    for(;k;k>>=1,poly_sqa(A,A,n),clr(A+n,n))
        if(k&1)poly_mul(g,A,g,n,n),clr(g+n,n);
}
void newton(int*g,int n){
    static int A[M],B[M],C[M];
    if(n==1)return g[0]=0,void();
    newton(g,(n+1)>>1);
    clr(A,n),A[0]=1;
    for(int i=0;i<n;++i)A[i]=(A[i]-2ll*g[i]%mod+mod)%mod;
    clr(B,n),poly_qpow(B,g,V-1,n);
    for(int i=1;i<n;++i)A[i]=(A[i]+1ll*V*B[i-1]%mod)%mod;
    clr(C,n),poly_inv(C,A,n);

    poly_mul(B,g,B,n,n),clr(B+n,n);
    for(int i=0;i<n;++i)A[i]=g[i];
    poly_sqa(A,A,n);
    for(int i=0;i<n;++i)fmod(A[i]=g[i]-A[i]+mod);
    fmod(A[1]+=mod-1);
    for(int i=1;i<n;++i)fmod(A[i]+=B[i-1]);
    poly_mul(A,C,A,n,n);

    for(int i=0;i<n;++i)fmod(g[i]+=mod-A[i]);
}
void calcf(int*f,int n){
    static int A[M],B[M];
    if(V>=127){
        f[0]=1;
        for(int k=1;k<=m;++k){
            for(int i=0;i<=k/V;++i){
                int tmp=1ll*comb(k,i)*comb(2*k-V*i-1,k-i*V)%mod;
                fmod(f[k]+=i&1?mod-tmp:tmp);
            }
        }
        return;
    }
    newton(A,n+1),dao(f,A,n+1);
    for(int i=0;i<n;++i) A[i]=A[i+1];
    clr(B,n),poly_inv(B,A,n);
    poly_mul(f,B,f,n,n),clr(f+n,n);
}

}

void solve(int l,int r){
    if(l==r)return;
    int mid=(l+r)>>1;
    int p=1ll*l*(n-r+1)%mod*IV%mod,q=mod+1-p;
    P[tot]=p,a[tot++]=1ll*p*qpow(q,mod-2)%mod;
    solve(l,mid),solve(mid+1,r);
}

signed main(){
    math::initmath();
    n=read(),m=read(),V=read()+2;
    poly::calcf(f,m+1);
    IV=qpow(1ll*n*(n+1)/2%mod,mod-2);
    solve(1,n);
    for(int i=0,j=1;i<=m;++i,j=1ll*j*V%mod)
        c[i]=1ll*comb(m,i)*(j+mod-f[i])%mod*qpow(j,mod-2)%mod;
    poly::poly_eva(b,c,a,m+1,tot);
    for(int i=0;i<tot;++i)fmod(ans+=1ll*b[i]*qpow(mod+1-P[i],m)%mod);
    cout<<ans<<'\n';
}