BZOJ 4555: [Tjoi2016&Heoi2016]求和 (NTT + 第二类斯特林数)

题意

给你一个数 n 求这样一个函数的值 :

f(n)=i=0nj=0i{ij}2jj!

(1n100000)

题解

这个题直接划式子 然后 NTT 就行了qwq

需要知道一个容斥求斯特林数的东西

{nm}=1m!k=0m(1)k(mk)n

这个用组合意义去理解 我们考虑把 n 个物品放进 m 个盒子里中的方案数就是斯特林数(盒子不区分)

然后我们考虑枚举至少有 k 个空的盒子的方案数 那么 n 就可以随便放入剩下的 (mk) 个盒子中去

这个式子我们划一下 就可以得到一个用来 NTT 的式子...

拆一下组合数....

{nm}=1m!k=0m(1)km!k!(mk)!(mk)n

然后再简单整理一下qwq

{nm}=k=0m[(1)kk!][(mk)n(mk)!]

然后这个就是 NTT 的式子了

对于这道题我们也可以这样做qwq

f(n)=i=0nj=0i{ij}×2j×(j!)

如果 j>i{ij} 是为 0 的 (没有方案数) 那么就有

=i=0nj=0n{ij}×2j×(j!)

把之前的那个套进来 卷积形式 我们可以将 jkk 互换

=i=0nj=0n2j×(j!)k=0j[(1)jk(jk)!][kik!]

不难发现只有一个地方与 i 有关 那么我们再放进去

=j=0n2j×(j!)k=0j[(1)jk(jk)!][i=0nkik!]

然后那个是等比数列 我们用等比数列求和公式 就可以直接处理出来了

=j=0n2j×(j!)k=0j[(1)jk(jk)!][kn+11(k1)k!]

右边卷积就可以求出对于每个 j 的取值咯qwq

注意程序中卷积之前 要特判右边等比数列次数 =0,1 的答案 一个是 1 另一个是 n+1

代码

#include <bits/stdc++.h> #define For(i, l, r) for(register int i = (l), i##end = (int)(r); i <= i##end; ++i) #define Fordown(i, r, l) for(register int i = (r), i##end = (int)(l); i >= i##end; --i) #define Set(a, v) memset(a, v, sizeof(a)) using namespace std; inline bool chkmin(int &a, int b) {return b < a ? a = b, 1 : 0;} inline bool chkmax(int &a, int b) {return b > a ? a = b, 1 : 0;} inline int read() { int x = 0, fh = 1; char ch = getchar(); for (; !isdigit(ch); ch = getchar()) if (ch == '-') fh = -1; for (; isdigit(ch); ch = getchar()) x = (x * 10) + (ch ^ 48); return x * fh; } void File() { #ifdef zjp_shadow freopen ("4555.in", "r", stdin); freopen ("4555.out", "w", stdout); #endif } typedef long long ll; const ll Mod = 998244353; ll fpm(ll x, int power) { ll res = 1; for (; power; power >>= 1, (x *= x) %= Mod) if (power & 1) (res *= x) %= Mod; return res; } const int N = 1 << 19; struct Number_Theoretic_Transform { ll pow3[N], invpow3[N], a[N], b[N]; int n, m, rev[N]; void Init(int n1, int n2, ll A[], ll B[]) { For (i, 0, n1) a[i] = A[i]; For (i, 0, n2) b[i] = B[i]; m = n1 + n2; } void NTT(ll P[], int opt) { For (i, 0, n - 1) if (i < rev[i]) swap(P[i], P[rev[i]]); for (int i = 2, p; i <= n; i <<= 1) { p = (i >> 1); ll Wi = (opt == 1) ? pow3[i] : invpow3[i]; for (int j = 0; j < n; j += i) { ll x = 1; for (int k = 0; k < p; ++ k, (x *= Wi) %= Mod) { ll u = P[j + k], v = x * P[j + k + p] % Mod; P[j + k] = (u + v) % Mod; P[j + k + p] = (u - v + Mod) % Mod; } } } if (opt == -1) { ll invn = fpm(n, Mod - 2); For (i, 0, n - 1) (P[i] *= invn) %= Mod; } } void Mult() { int cnt = 0; for (n = 1; n <= m; n <<= 1) ++ cnt; for (int i = 1; i <= n; i <<= 1) pow3[i] = fpm(3, (Mod - 1) / i), invpow3[i] = fpm(pow3[i], Mod - 2); For (i, 1, n) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (cnt - 1)); NTT(a, 1); NTT(b, 1); For (i, 0, n - 1) (a[i] *= b[i]) %= Mod; NTT(a, -1); } } T; int n; ll a[N], b[N], fac[N], ifac[N]; void Init(int maxn) { fac[0] = ifac[0] = 1; For (i, 1, maxn) fac[i] = fac[i - 1] * i % Mod; ifac[maxn] = fpm(fac[maxn], Mod - 2); Fordown (i, maxn - 1, 1) ifac[i] = ifac[i + 1] * (i + 1) % Mod; } ll ans = 0; int main () { File(); n = read(); Init(n); For (i, 0, n) { a[i] = (Mod + ((i & 1) ? -1 : 1) * ifac[i]) % Mod; if (i > 1) b[i] = (fpm(i, n + 1) - 1) * ifac[i] % Mod * fpm(i - 1, Mod - 2) % Mod; else if (i == 1) b[i] = n + 1; else if (i == 0) b[i] = 1; // printf ("a[%d] = %lld; b[%d] = %lld;\n", i, a[i], i, b[i]); } T.Init(n, n, a, b); T.Mult(); For (i, 0, n) (ans += T.a[i] * fpm(2, i) % Mod * fac[i]) %= Mod; printf ("%lld\n", ans); return 0; }

__EOF__

本文作者zjp_shadow
本文链接https://www.cnblogs.com/zjp-shadow/p/8763277.html
关于博主:评论和私信会在第一时间回复。或者直接私信我。
版权声明:本博客所有文章除特别声明外,均采用 BY-NC-SA 许可协议。转载请注明出处!
声援博主:如果您觉得文章对您有帮助,可以点击文章右下角推荐一下。您的鼓励是博主的最大动力!
posted @   zjp_shadow  阅读(374)  评论(0编辑  收藏  举报
编辑推荐:
· 如何编写易于单元测试的代码
· 10年+ .NET Coder 心语,封装的思维:从隐藏、稳定开始理解其本质意义
· .NET Core 中如何实现缓存的预热?
· 从 HTTP 原因短语缺失研究 HTTP/2 和 HTTP/3 的设计差异
· AI与.NET技术实操系列:向量存储与相似性搜索在 .NET 中的实现
阅读排行:
· 地球OL攻略 —— 某应届生求职总结
· 周边上新:园子的第一款马克杯温暖上架
· Open-Sora 2.0 重磅开源!
· 提示词工程——AI应用必不可少的技术
· .NET周刊【3月第1期 2025-03-02】
点击右上角即可分享
微信分享提示