LOJ6503 「雅礼集训 2018 Day4」Magic

Link

考虑经典二项式反演，设\(f_n\)为恰好\(n\)个魔术对的方案数，\(g_n\)为至少有\(n\)个魔术对的方案数，那么

\[g_i=\sum\limits_{j\ge i}{j\choose i}f_j,f_i=\sum\limits_{j\ge i}(-1)^{j-i}{j\choose i}g_j \]

对于第\(i\)种卡片，假如它在序列中出现了\(j\)个连续段，那么就会有\(a_i-j\)对魔术对，方案数为\(h_{i,j}={a_i-1\choose j-1}\)。
\(h\)合成到\(g\)的过程是一个有标号对象拼接的过程，考虑用EGF表示。

\[\widehat H_i(x)=\sum\limits_{j\ge0}\frac{h_{i,j}x^j}{j!},\widehat G(x)=\prod\limits_{i=1}^m\widehat H_i(x) \]

注意到\(\sum\limits_{i=1}^m(a_i-x_i)=n-\sum\limits_{i=1}^mx_i\)，因此\(g_k=(n-k)![x^k]\widehat G(x)\)。
注意到\(\deg\widehat G=n\)，因此这个连乘积可以分治NTT。

#include<queue>
#include<cctype>
#include<cstdio>
#include<vector>
#include<cstring>
#include<algorithm>
using i64=long long;
using poly=std::vector<i64>;
const int N=1<<18,P=998244353;
int len,a[N],rev[N];i64 w[N],fac[N],ifac[N];
struct cmp{int operator()(const poly&a,const poly&b){return a.size()>b.size();}};
std::priority_queue<poly,std::vector<poly>,cmp>q;
int read(){int x=0,c=getchar();while(isspace(c))c=getchar();while(isdigit(c))(x*=10)+=c&15,c=getchar();return x;}
void inc(i64&a,i64 b){a+=b-P,a+=a>>63&P;}
void dec(i64&a,i64 b){a-=b,a+=a>>63&P;}
i64 C(int n,int m){return fac[n]*ifac[m]%P*ifac[n-m]%P;}
i64 pow(i64 a,int b){i64 r=1;for(;b;b>>=1,a=a*a%P)if(b&1)r=r*a%P;return r;}
int getlen(int n){return 1<<(32-__builtin_clz(n));}
void init(int n)
{
    int lim=1<<(len=32-__builtin_clz(n)),pos=lim/2;i64 g=pow(3,(P-1)/lim);
    for(int i=1;i<lim;++i) rev[i]=(rev[i>>1]>>1)|(i&1? pos:0);
    w[pos]=1;for(int i=pos+1;i<lim;++i) w[i]=g*w[i-1]%P;
    for(int i=pos-1;i;--i) w[i]=w[i<<1];
    fac[0]=1;for(int i=1;i<=n;++i) fac[i]=i*fac[i-1]%P;
    ifac[n]=pow(fac[n],P-2);for(int i=n;i;--i) ifac[i-1]=i*ifac[i]%P;
}
void NTT(i64*a,int lim,int f)
{
    if(!~f) std::reverse(a+1,a+lim);
    for(int i=0,x=len-__builtin_ctz(lim);i<lim;++i) if(i<rev[i]>>x) std::swap(a[i],a[rev[i]>>x]);
    for(int i=1;i<lim;i<<=1) for(int j=0,d=i<<1;j<lim;j+=d) for(int k=0,x;k<i;++k) x=w[i+k]*a[i+j+k]%P,dec(a[i+j+k]=a[j+k],x),inc(a[j+k],x);
    if(!~f) for(int i=0,x=P-(P-1)/lim;i<lim;++i) a[i]=x*a[i]%P;
}
void push(int deg)
{
    poly a(deg+1);
    for(int i=1;i<=deg;++i) a[i]=C(deg-1,i-1)*ifac[i]%P;
    q.push(a);
}
void pop()
{
    static i64 f[N],g[N];
    poly a=q.top();q.pop();poly b=q.top();q.pop();
    int n=a.size(),m=b.size(),lim=getlen(n+m-2);poly c(n+m-1);
    memcpy(f,&a[0],8*n),memcpy(g,&b[0],8*m),memset(f+n,0,8*(lim-n)),memset(g+m,0,8*(lim-m));
    NTT(f,lim,1),NTT(g,lim,1);
    for(int i=0;i<lim;++i) f[i]=f[i]*g[i]%P;
    NTT(f,lim,-1),memcpy(&c[0],f,8*(n+m-1)),q.push(c);
}
int main()
{
    int m=read(),n=read(),k=read();i64 ans=0;init(n);
    for(int i=1;i<=m;++i) a[i]=read(),push(a[i]);
    while(q.size()>1) pop();
    poly f=q.top();q.pop();
    for(int i=k;i<=n-m;++i) ((i^k)&1? dec:inc)(ans,C(i,k)*f[n-i]%P*fac[n-i]%P);
    printf("%lld",ans);
}