洛谷 P4091 [HEOI2016/TJOI2016]求和

在2016年，佳媛姐姐刚刚学习了第二类斯特林数，非常开心。

现在他想计算这样一个函数的值:

[f(n)=sum_{i=0}^nsum_{j=0}^i S(i,j) imes 2^j imes (j!) ]

S(i, j)表示第二类斯特林数，递推公式为:

(S(i, j) = j imes S(i - 1, j) + S(i - 1, j - 1), 1 le j le i - 1)。

边界条件为：(S(i, i) = 1(0 le i), S(i, 0) = 0(1 le i))

你能帮帮他吗?

不会线性做法，只会无脑化式子/kk

[egin{aligned} f(n)&=sum_{i=0}^nsum_{j=0}^iegin{Bmatrix}i\jend{Bmatrix}2^jj!\ &=sum_{i=0}^nsum_{j=0}^negin{Bmatrix}i\jend{Bmatrix}2^jj!\ &=sum_{j=0}^n2^jj!sum_{j=0}^negin{Bmatrix}i\jend{Bmatrix} end{aligned}]

考虑第二类斯特林数的容斥式(egin{Bmatrix}n\mend{Bmatrix}=frac{1}{m!}sum_{k=0}^m(-1)^kegin{pmatrix}m\kend{pmatrix}(m-k)^n)，则有

[egin{aligned} &=sum_{j=0}^n2^jj!sum_{i=0}^nfrac{1}{j!}sum_{k=0}^i(-1)^kegin{pmatrix}j\kend{pmatrix}(j-k)^i \ &=sum_{j=0}^n2^jsum_{i=0}^nsum_{k=0}^i(-1)^kegin{pmatrix}j\kend{pmatrix}(j-k)^i\ &=sum_{j=0}^n2^jsum_{i=0}^nsum_{k=0}^i(-1)^kfrac{j!}{k!(j-k)!}(j-k)^i\ &=sum_{j=0}^n2^jj!sum_{i=0}^nsum_{k=0}^ifrac{(-1)^k}{k!} imesfrac{(j-k)^i}{(j-k)!}\ &=sum_{j=0}^n2^jj!sum_{k=0}^jfrac{(-1)^k}{k!} imesfrac{sum_{i=0}^n(j-k)^i}{(j-k)!} end{aligned}]

设(A(i)=frac{(-1)^i}{i!},B(i)=frac{sum_{j=0}^ni^j}{i!}=frac{i^{n+1}-1}{i!(i-1)})

那么就有

[=sum_{j=0}^n2^jj!sum_{k=0}^jA(k) imes B(j-k) ]

直接ntt就可以了

Code

#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstring>
const int N = 4e5;
const int p = 998244353;
using namespace std;
int n,a[N + 5],b[N + 5],rev[N + 5],maxn,lg,G[N + 5][2],fac[N + 5],ans;
int mypow(int a,int x){int s = 1;for (;x;x & 1 ? s = 1ll * s * a % p : 0,a = 1ll * a * a % p,x >>= 1);return s;}
void ntt(int *a,int typ)
{
    for (int i = 0;i < maxn;i++)
        if (i < rev[i])
            swap(a[i],a[rev[i]]);
    for (int i = 1;i < maxn;i <<= 1)
        for (int j = 0;j < maxn;j += i << 1)
            for (int k = 0;k < i;k++)
            {
                int x = a[j + k],t = 1ll * G[i + k][typ] * a[j + k + i] % p;
                a[j + k] = (x + t) % p;
                a[j + k + i] = (x - t) % p;
            }
    if (!typ)
    {
        int inv = mypow(maxn,p - 2);
        for (int i = 0;i < maxn;i++)
            a[i] = 1ll * a[i] * inv % p;
    }
}
int main()
{
    scanf("%d",&n);
    maxn = 1;
    while (maxn <= n * 2 + 2)
        maxn <<= 1,lg++;
    for (int i = 0;i < maxn;i++)
        rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << lg - 1);
    for (int i = 1;i < maxn;i <<= 1)
    {
        int g1 = mypow(3,(p - 1) / (i << 1)),g2 = mypow(mypow(3,p - 2),(p - 1) / (i << 1));
        G[i][0] = G[i][1] = 1;
        for (int j = 1;j < i;j++)
            G[i + j][1] = 1ll * G[i + j - 1][1] * g1 % p,G[i + j][0] = 1ll * G[i + j - 1][0] * g2 % p;
    }
    a[0] = b[0] = fac[0] = 1;
    for (int i = 1;i <= n;i++)
        fac[i] = 1ll * fac[i - 1] * i % p;
    for (int i = 1;i <= n;i++)
        a[i] = mypow(-1,i) * mypow(fac[i],p - 2),b[i] = 1ll * (mypow(i,n + 1) - 1) * mypow(1ll * fac[i] * (i - 1) % p,p - 2) % p;
    b[1] = n + 1;
    ntt(a,1);
    ntt(b,1);
    for (int i = 0;i < maxn;i++)
        a[i] = 1ll * a[i] * b[i] % p;
    ntt(a,0);
    for (int i = 0;i <= n;i++)
        ans += 1ll * mypow(2,i) * fac[i] % p * a[i] % p,ans %= p;
    cout<<(ans + p) % p<<endl;
    return 0;
}