[CTS2019]珍珠——生成函数+二项式反演

题面

　　LOJ3120

解析

　　考虑如何判断一个序列是否满足条件，设$c_i$表示数$i$出现的次数，则：$$sum_{i=1}^Dleftlfloor frac{c_i}{2} ight floorgeqslant m\ sum_{i=1}^Dfrac{c_i-(c_i\%2)}{2}geqslant m\n-sum_{i=1}^D(c_i\%2)geqslant 2m\ sum_{i=1}^D(c_i\%2)leqslant n-2m$$

　　$n-2m<0$时，输出$0$；$n-2mgeqslant D$时，输出$D^n$

　　设$g_i$表示恰有$i$个数出现次数为奇数，$f_i$表示至少有$i$个数出现次数为奇数，且$f_i=sum_{j=i}inom{j}{i}g_j$

　　如果我们求出$f$那么就可以二项式反演求出$g$，因此现在需要想办法求出$f$

　　考虑生成函数，令$A(x)=sum_{i=0}^{infty}[i\%2==1]x^i=frac{e^x-e^{-x}}{2}$，$B(x)=sum_{i=0}^{infty}x^i=e^x$

　　钦定有$i$个数出现奇数次，剩下的数出现任意次数，则有：$$egin{align*}f_i&=inom{D}{i}n![x^n]A^iB^{D-i}\&=inom{D}{i}n^ie^{(D-i)x}\&=frac{D!}{i!(D-i)!}n!*frac{1}{2^i}[x^n]sum_{j=0}^iinom{i}{j}(-1)^{i-j}e^{jx}e^{-(i-j)x}e^{(D-i)x}\&=frac{D!n!}{i!(D-i)!2^i}*[x^n]sum_{j=0}^iinom{i}{j}(-1)^{i-j}e^{(D-2i+2j)x}\&=frac{D!n!}{i!(D-i)!2^i}sum_{j=0}^ifrac{i!}{j!(i-j)!}*(-1)^{i-j}*frac{(D-2i+2j)^n}{n!}\&=frac{D!}{(D-i)!2^i}sum_{j=0}^ifrac{(-1)^{i-j}(D-2i+2j)^n}{(i-j)!}*frac{1}{j!}end{align*}$$

　　卷积即可求出$f$

　　然后二项式反演：$$egin{align*}g_i&=sum_{j=i}^D(-1)^{j-i}inom{j}{i}f_j\&=frac{1}{i!}sum_{j=i}^Dfrac{(-1)^{j-i}}{(j-i)!}*j!*f_jend{align*}$$

　　最后的答案就是：$$Ans=sum_{i=0}^{n-2m}g_i$$

　　$O(D log D)$

　代码：

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
typedef long long ll;
const int maxn = 200005, mod = 998244353, g = 3;

int add(int x, int y)
{
    return x + y < mod? x + y: x + y - mod;
}

int rdc(int x, int y)
{
    return x - y < 0? x - y + mod: x - y;
}

ll qpow(ll x, int y)
{
    ll ret = 1;
    while(y)
    {
        if(y&1)
            ret = ret * x % mod;
        x = x * x % mod;
        y >>= 1;
    }
    return ret;
}

int D, n, m;
int lim, bit, rev[maxn<<1];
ll fac[maxn], fnv[maxn];
ll ginv, f[maxn<<1], A[maxn<<1], B[maxn<<1];

void init()
{
    ginv = qpow(g, mod - 2);
    fac[0] = 1;
    for(int i = 1; i <= D; ++i)
        fac[i] = fac[i-1] * i % mod;
    fnv[D] = qpow(fac[D], mod - 2);
    for(int i = D - 1; i >= 0; --i)
        fnv[i] = fnv[i+1] * (i + 1) % mod;
}

void NTT_init(int x)
{
    lim = 1;
    bit = 0;
    while(lim <= x)
    {
        lim <<= 1;
        ++ bit;
    }
    for(int i = 1; i < lim; ++i)
        rev[i] = (rev[i>>1] >> 1) | ((i & 1) << (bit - 1));
}

void NTT(ll *x, int y)
{
    for(int i = 1; i < lim; ++i)
        if(i < rev[i])
            swap(x[i], x[rev[i]]);
    ll wn, w, u, v;
    for(int i = 1; i < lim; i <<= 1)
    {
        wn = qpow((y == 1)? g: ginv, (mod - 1) / (i << 1));
        for(int j = 0; j < lim; j += (i << 1))
        {
            w = 1;
            for(int k = 0; k < i; ++k)
            {
                u = x[j+k];
                v = x[j+k+i] * w % mod;
                x[j+k] = add(u, v);
                x[j+k+i] = rdc(u, v);
                w = w * wn % mod;
            }
        }
    }
    if(y == -1)
    {
        ll linv = qpow(lim, mod - 2);
        for(int i = 0; i < lim; ++i)
            x[i] = x[i] * linv % mod;
    }
}

int main()
{
    scanf("%d%d%d", &D, &n, &m);
    if(n - 2 * m < 0)
    {
        printf("0");
        return 0;
    }
    if(n - 2 * m >= D)
    {
        printf("%lld", qpow(D, n));
        return 0;
    }
    init();
    for(int i = 0; i <= D; ++i)
    {
        A[i] = qpow(rdc(D, 2 * i), n) * fnv[i] % mod;
        A[i] = ((i & 1)? rdc(0, A[i]): A[i]);
        B[i] = fnv[i];
    }
    NTT_init(D << 1);
    NTT(A, 1);
    NTT(B, 1);
    for(int i = 0; i < lim; ++i)
        A[i] = A[i] * B[i] % mod;
    NTT(A, -1);
    for(int i = 0; i <= D; ++i)
        A[i] = ((A[i] * fac[D] % mod) * qpow(qpow(2, i), mod - 2) % mod) * fnv[D-i] % mod;
    for(int i = D + 1; i < lim; ++i)
        A[i] = 0;
    for(int i = 0; i <= D; ++i)
    {
        A[i] = A[i] * fac[i] % mod;
        f[D-i] = ((i & 1)? rdc(0, fnv[i]): fnv[i]);
    }
    NTT(A, 1);
    NTT(f, 1);
    for(int i = 0; i < lim; ++i)
        f[i] = f[i] * A[i] % mod;
    NTT(f, -1);
    int ans = 0;
    for(int i = 0; i <= n - 2 * m; ++i)
        ans = add(ans, f[D+i] * fnv[i] % mod);
    printf("%d", ans);
    return 0;
}