「CF1349F2」Slime and Sequences (Hard Version)
「CF1349F2」Slime and Sequences (Hard Version)

2020 年 05 月 14 日

定义一个排列 pp 是好的当且仅当对于每个 k<max{p}k < \max\{p\},存在 1i<jn1 \leq i < j \leq n 使得 ai=k1a_i = k-1aj=ka_j = k

定义 fa(k)f_a(k) 为序列 aa 中数值 kk 的出现次数,假设所有合法序列集合为 SS,对于每个 k[1;n]k \in [1;n],求

(aSfa(k))mod998244353\left( \sum_{a \in S} f_a(k) \right) \bmod 998244353

n105n \leq 10^5

题解

考虑由排列 pp 生成序列,在 pip_ipi+1p_{i+1} 之间填入大于或小于号,则 apia_{p_i} 的值为 p1...ip_{1...i} 中的小于号个数 +1+1,不难验证合法序列集合和排列集合构成双射。

注意到这是个欧拉数的形式,考虑容斥维护:

ansi=m=0nn!m!ji(mj)!(1)ji(ji){mmj}=jin!(1)ji(ji)m=j+1n{mmj}(mj)!m!\newcommand{\strling}[2]{\left\{\begin{matrix}#1\\#2\end{matrix}\right\}} \begin{aligned} ans_i &= \sum_{m=0}^n \frac {n!} {m!} \sum_{j \ge i} (m-j)! (-1)^{j-i}\binom j i \strling m {m-j} \\ &= \sum_{j \ge i} n! (-1)^{j-i} \binom j i \sum_{m=j+1}^n \strling m {m-j} \frac {(m-j)!} {m!} \end{aligned}

第二类斯特林数的生成函数

{nk}=n![zn](ez1)kk!\newcommand{\strling}[2]{\left\{\begin{matrix}#1\\#2\end{matrix}\right\}} \strling n k = n! [z^n] \frac {(e^z-1)^k} {k!}

hh,先带上 m=im=i 的情况计算,最后再令 h0h_011

hi=m=in{mmi}(mi)!m!=m=in[zm](ez1)mi=j=0ni[zj+i](ez1)j=[zi]j=0ni(ez1z)j\newcommand{\strling}[2]{\left\{\begin{matrix}#1\\#2\end{matrix}\right\}} \begin{aligned} h_i &= \sum_{m=i}^n \strling m {m-i} \frac{(m-i)!} {m!} \\ &= \sum_{m=i}^n [z^m] (e^z-1)^{m-i}\\ &= \sum_{j=0}^{n-i} [z^{j+i}] (e^z-1)^j \\ &= [z^i] \sum_{j=0}^{n-i} \left( \frac{e^z-1} {z} \right)^j \\ \end{aligned}

F(z)=(ez1)/zF(z) = (e^z-1)/z,有

hi=[zi]j=0niFj(z)=[zi]Fni+1(z)1F(z)1=[zi](1F(z)1+Fni+1(z)F(z)1)\begin{aligned} h_i &= [z^i] \sum_{j=0}^{n-i} F^j(z) \\ &= [z^i] \frac{F^{n-i+1}(z)-1} {F(z)-1} \\ &= [z^i] \left( \frac{-1}{F(z)-1} + \frac{F^{n-i+1}(z)}{F(z)-1} \right) \\ \end{aligned}

前半部分容易直接多项式求逆处理,现在考虑后半部分

[zi]Fni+1(z)F(z)1=[zn+1](zF(z))ni+1F(z)1[z^i] \frac{F^{n-i+1}(z)}{F(z)-1} = [z^{n+1}] \frac{(z F(z))^{n-i+1}}{F(z)-1}

(开始在这里卡住了,不求甚解的从题解那里拉了个式子,第二天冷静了一下,为了行文连贯,把补充理解附在后面)

ω(z)=zF(z)\omega(z) = zF(z)φ(z)\varphi(z) 满足 ω(z)φ(ω(z))=z\dfrac{\omega(z)}{\varphi(\omega(z))}=z,则 zF(z)φ(ω(z))=z\dfrac{zF(z)}{\varphi(\omega(z))}=z,即 F=φ(ω)F=\varphi(\omega)

拓展拉格朗日反演

f,g,hf,g,hF[[x]]F[[x]] 上的多项式,已知 f(g(x))=g(f(x))=xf(g(x)) = g(f(x)) = x,则

[xn]h(f(x))=1n[xn1]h(x)(xg(x))n[x^n] h(f(x)) = \frac 1 n [x^{n-1}] h'(x) \left( \frac {x} {g(x)} \right)^n

考虑

f(x)=ω(x)g(x)=xφ(x)h(x)=1(1φ(x))(1ux)\begin{aligned} f(x) &= \omega(x) \\ g(x) &= \frac x {\varphi(x)} \\ h(x) &= \frac 1 {(1-\varphi(x)) (1 - ux)} \\ \end{aligned}

[uni+1zn+1]11φ(ω(z))11uw(z)=[uni+1]1n+1[zn]((11φ(z)11uz)φ(z)n+1)\begin{aligned} & [u^{n-i+1} z^{n+1}] \frac 1 {1-\varphi (\omega(z))} \frac 1 {1-uw(z)} \\ =& [u^{n-i+1}] \frac 1 {n+1} [z^n] \left(\left(\frac 1 {1-\varphi(z)} \frac 1 {1-uz}\right)'\cdot \varphi(z)^{n+1}\right) \\ \end{aligned}

1(1x)k=i=0(i+k1k1)xi\displaystyle \frac 1 {(1-x)^k} = \sum_{i=0}^{\infty} \binom {i+k-1} {k-1} x^i

直接考虑系数组合意义就可以证明。

其中

(11φ(z)11uz)=φ(z)(1uz)(1φ(z))2+u(1φ(z))(1uz)2=φ(z)(1φ(z))2i=0uizi+11φ(z)i=0(i+1)ui+1zi\begin{aligned} & \left(\frac 1 {1-\varphi(z)} \frac 1 {1-uz}\right)' \\ =& \frac {\varphi'(z)} {(1-uz) (1-\varphi(z))^2} + \frac {u} {(1-\varphi(z))(1-uz)^2} \\ =& \frac {\varphi'(z)} {(1-\varphi(z))^2} \sum_{i=0}^\infty u^i z^i + \frac {1} {1 - \varphi(z)} \sum_{i=0}^\infty (i+1) u^{i+1} z^i \\ \end{aligned}

代入到原式中有

[zn]φn+1(z)n+1(zni+1φ(z)(1φ(z))2+(ni+1)zni1φ(z))=[zi1]φn+1(z)φ(z)(n+1)(1φ(z))2+[zi]φn+1(z)(ni+1)(n+1)(1φ(z))\begin{aligned} & \frac {[z^n] \varphi^{n+1}(z)} {n+1} \left( \frac {z^{n-i+1} \varphi'(z)} {(1-\varphi(z))^2} + \frac {(n-i+1) z^{n-i}} {1 - \varphi(z)} \right) \\ =& [z^{i-1}] \frac {\varphi^{n+1}(z) \varphi'(z)} {(n+1) (1-\varphi(z))^2} + [z^i] \frac {\varphi^{n+1}(z) (n-i+1)} {(n+1) (1-\varphi(z))} \\ \end{aligned}

唯一的问题就是 φ(x)\varphi(x) 怎么求了,注意到 ω(x)=ex1\omega(x) = e^x-1,构造得 zln(1+z)\dfrac z {\ln (1+z)}

坑点

  1. φ(z)=zln(1+z)\varphi(z) = \dfrac {z} {\ln(1+z)},注意到 ln(1+z)\ln(1+z) 的常数项是 00,不能直接求逆;
  2. 1φ(z)1-\varphi(z) 常数项也是 00,同样不能直接求逆,只能求出 (1φ(z))/z(1-\varphi(z))/z,上面的式子大概变成
    [zi+1]φn+1(z)φ(z)(n+1)((1φ(z))/z)2+[zi+1]φn+1(z)(ni+1)(n+1)((1φ(z)/z)) [z^{i+1}] \frac {\varphi^{n+1}(z) \varphi'(z)} {(n+1) ((1-\varphi(z))/z)^2} + [z^{i+1}] \frac {\varphi^{n+1}(z) (n-i+1)} {(n+1) ((1-\varphi(z)/z))}
  3. 上面两种情况的处理可能带来更多的多项式长度要求。
  4. h0h_011

补充

现在我们需要对 i[0;n]i \in [0;n],求出

[zn+1](zF(z))ni+1F(z)1[z^{n+1}] \frac{(z F(z))^{n-i+1}}{F(z)-1}

考虑用二元生成函数表示

[zn+1uni+1]1F(z)1k=0(zuF(z))k=[zn+1uni+1]1(F(z)1)(1zuF(z))\begin{aligned} & [z^{n+1} u^{n-i+1}] \frac 1 {F(z)-1} \sum_{k=0}^\infty (z u F(z))^{k} \\ =& [z^{n+1} u^{n-i+1}] \frac 1 {(F(z)-1) (1-zuF(z))} \\ \end{aligned}

如果直接对着这个式子做,是没有办法拉格朗日反演的,因为本质上我们有两个关于 zz 的多项式:F(z)F(z)zF(z)zF(z)。我们设法构造一个关于 FF 的多项式使其满足 φ(F(z))=zF(z)\varphi(F(z)) = z F(z)

形式化的,我们想要构造多项式 φ(z)\varphi(z) 使得

φ(F(z))=zF(z)\varphi(F(z)) = z F(z)

同时,由该式我们也可以直接得到 FF 的复合逆形式:

φ(F(z))F(z)=z\frac {\varphi(F(z))} {F(z)} = z

说回到 φ(x)\varphi(x) 的构造,我们只能利用 F(z)F(z) 本身的性质设法构造出 zz:我们有

F(z)=ez1zF(z) = \frac {e^{z} - 1} {z}

相当于有

ln(1+zF(z))=z\ln(1 + zF(z)) = z

不幸的是,这个式子本身还和 zz 有关,但启发我们改变方向,设 ω(z)=zF(z)\omega(z) = z F(z),构造多项式 φ(z)\varphi(z) 使得

{ω(z)=zF(z)φ(ω(z))=F(z)\begin{cases} \omega(z) = z F(z) \\ \varphi(\omega(z)) = F(z) \\ \end{cases}

φ(z)=zln(1+z)\varphi(z) = \dfrac {z} {\ln(1 + z)}ω(z)φ(ω(z))=z\dfrac {\omega(z)} {\varphi(\omega(z))} = z。即把原式转化为了:

1(1uω(z))(ω(z)ln(1+ω(z))1)\frac 1 {(1 - u \omega(z)) \left(\dfrac {\omega(z)} {\ln(1 + \omega(z))}-1\right)}

代码

#include <bits/stdc++.h>
using namespace std;
const int N = 1e5 + 9, mod = 998244353;
struct z {
  int x;
  z(int x = 0) : x(x) {}
  friend inline z operator*(z a, z b) { return (long long)a.x * b.x % mod; }
  friend inline z operator-(z a, z b) { return (a.x -= b.x) < 0 ? a.x + mod : a.x; }
  friend inline z operator+(z a, z b) { return (a.x += b.x) >= mod ? a.x - mod : a.x; }
} h[N], ans[N], inv[N], fac[N], ifac[N];
inline z fpow(z a, int b) {
  z s = 1;
  for (; b; b >>= 1, a = a * a)
    if (b & 1) s = s * a;
  return s;
}
int len = 1, rev[N << 2];
z w[N << 2];
struct vec : vector<z> {
  using vector<z>::vector;
  inline vec divx() {
    vec res = *this;
    return res.erase(res.begin()), res;
  }
  inline vec setl(size_t len) {
    vec res = *this;
    return res.resize(len), res;
  }
  inline vec fun1() {
    vec res(this->begin() + 1, this->end());
    for (int i = 0; i < res.size(); i++) res[i].x = mod - res[i].x;
    return res;
  }
};
int init(int n) {
  int lim = 1, k = 0;
  while (lim < n) lim <<= 1, ++k;
  for (int i = 0; i < lim; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (k - 1));
  for (; len < lim; len <<= 1) {
    z wn = fpow(3, (mod - 1) / (len << 1));
    w[len] = 1;
    for (int i = 1; i < len; i++) w[i + len] = w[i + len - 1] * wn;
  }
  return lim;
}
void dft(vec &a, int lim) {
  a.resize(lim);
  for (int i = 0; i < lim; i++)
    if (i < rev[i]) swap(a[i], a[rev[i]]);
  for (int i = 0; i < lim; i += 2) {
    z x = a[i], y = a[i + 1] * w[1];
    a[i] = x + y, a[i + 1] = x - y;
  }
  for (int len = 2; len < lim; len <<= 1)
    for (int i = 0; i < lim; i += (len << 1))
      for (int j = 0; j < len; j++) {
        z x = a[i + j], y = a[i + j + len] * w[j + len];
        a[i + j] = x + y, a[i + j + len] = x - y;
      }
}
void idft(vec &a, int lim) {
  z inv = fpow(lim, mod - 2);
  dft(a, lim), reverse(&a[1], &a[lim]);
  for (int i = 0; i < lim; i++) a[i] = a[i] * inv;
}
vec operator+(vec a, const vec &b) {
  a.resize(max(a.size(), b.size()));
  for (int i = 0; i < b.size(); i++) a[i] = a[i] + b[i];
  return a;
}
vec operator-(vec a, const vec &b) {
  a.resize(max(a.size(), b.size()));
  for (int i = 0; i < b.size(); i++) a[i] = a[i] - b[i];
  return a;
}
vec operator*(vec a, vec b) {
  if (a.size() < 20 || b.size() < 20 || (uint64_t)a.size() + b.size() < 400) {
    vec c(a.size() + b.size() - 1);
    for (int i = 0; i < a.size(); i++)
      for (int j = 0; j < b.size(); j++) c[i + j] = c[i + j] + a[i] * b[j];
    return c;
  }
  int len = a.size() + b.size() - 1, lim = init(len);
  dft(a, lim), dft(b, lim);
  for (int i = 0; i < lim; i++) a[i] = a[i] * b[i];
  return idft(a, lim), a.resize(len), a;
}
vec inve(const vec &f, int len = -1) {
  if ((len = ~len ? len : f.size()) == 1) return vec{fpow(f[0], mod - 2)};
  vec a = inve(f, (len + 1) >> 1), b(&f[0], &f[len]);
  int lim = init((len << 1) - 1);
  dft(a, lim), dft(b, lim);
  for (int i = 0; i < lim; i++) a[i] = a[i] * (2 - a[i] * b[i]);
  return idft(a, lim), a.resize(len), a;
}
vec inte(vec a) {
  for (int i = a.size() - 1; i; i--) a[i] = a[i - 1] * ::inv[i];
  return *a.begin() = 0, a;
}
vec deri(vec a) {
  for (int i = 0; i < a.size() - 1; i++) a[i] = a[i + 1] * (i + 1);
  return *a.rbegin() = 0, a;
}
vec ln(const vec &f) { return inte((deri(f) * inve(f)).setl(f.size())); }
vec exp(const vec &f, int len = -1) {
  if ((len = ~len ? len : f.size()) == 1) return vec{1};
  vec a = exp(f, (len + 1) >> 1), b = a;
  b.resize(len), b = ln(b);
  for (int i = 0; i < len; i++) b[i] = 0 - b[i];
  b[0] = b[0] + 1;
  for (int i = 0; i < len; i++) b[i] = b[i] + f[i];
  return (a * b).setl(len);
}
vec pow(vec a, int b) {
  a = ln(a);
  for (int i = 0; i < a.size(); i++) a[i] = a[i] * b;
  return exp(a);
}
int main() {
#ifdef memset0
  freopen("1.in", "r", stdin);
#endif
  cin.tie(0)->sync_with_stdio(0);
  fac[0] = ifac[0] = inv[0] = inv[1] = 1;
  for (int i = 2; i < N; i++) inv[i] = (mod - mod / i) * inv[mod % i];
  for (int i = 1; i < N; i++) fac[i] = fac[i - 1] * i, ifac[i] = ifac[i - 1] * inv[i];
  int n, len;
  cin >> n, len = n + 5;
  vec Finv = inve(vec(ifac + 1, ifac + len + 1).fun1());
  for (int i = 0; i <= n; i++) h[i] = Finv[i + 1];
  h[0] = h[0] - 1;
  vec P(len);
  for (int i = 0; i < len; i++) P[i] = (i & 1 ? mod - 1 : 1) * inv[i + 1];
  P = inve(P);
  vec Ppow = pow(P, n + 1);
  vec Pinv = inve(P.fun1());
  vec lpart = (Ppow * Pinv).setl(len);
  vec rpart = ((deri(P) * Pinv).setl(len) * lpart).setl(len);
  z tmp = fpow(n + 1, mod - 2);
  for (int i = 0; i <= n; i++) h[i] = h[i] - tmp * (lpart[i + 1] * (n - i + 1) + rpart[i + 1]);
  vec f(n + 1), g(n + 1);
  for (int i = 0; i <= n; i++) f[i] = ((n - i) & 1 ? mod - 1 : 1) * ifac[n - i];
  for (int i = 0; i <= n; i++) g[i] = fac[i] * h[i];
  f = f * g;
  for (int i = 0; i < n; i++) ans[i] = f[i + n] * fac[n] * ifac[i];
  for (int i = 0; i < n; i++) cout << ans[i].x << " \n"[i + 1 == n];
}

评论

TABLE OF CONTENTS

题解
坑点
补充
代码