【Luogu P4389】付公主的背包

Problem

Description

给出 \(n\) 种物品，每种物品体积为 \(v_i\) ，件数可视作无限。给定 \(m\) ，求出用这些商品装满容积分别为 \(1\sim m\) 的这 \(m\) 个背包的方案数。

Input Format

第一行两个正整数 \(n,m\) 。

接下来一行 \(n\) 个数，每个数 \(v_i\) 表示第 \(i\) 种物品的体积。

Output Format

\(m\) 行 \(m\) 个整数，第 \(i\) 行表示对于体积为 \(i\) 的背包有多少种方案，答案对 \(998244353\) 取模。

Range

\(n, m\le 10^5, v_i\le m\)

Algorithm

多项式，生成函数

Mentality

对于一种体积为 \(v_i\) 的物品，它的生成函数为：

\[f_i(x)=\sum_{k = 1} x^{v_ik} \]

可以得到封闭表达式：\(f_i(x)=\frac{1}{1-x^{v_i}}\) 。

显然有答案的生成函数 \(Ans=\prod_{i=1}^n f_i(x)\) 。

当然我们并不能把所有的函数乘起来，转而考虑一下利用 \(ln\) 转化成加法。

\[ln'(f_i(x)) = ln'(\frac{1}{1-x^{v_i}})\\ = (1-x^{v_i})*f_i'(x)\\ = (1-x^{v_i})* \sum v_ik *x^{v_ik-1}\\ = \sum v_i * x^{v_ik-1}\\ \]

这个时候我们再积分回去就得到了：

\[ln(f_i(x))=\sum \frac{1}{k} * x^{v_ik} \]

由于物品体积总数就那么一点，我们直接把相同体积的物品一起加到 \(Ans\) 的对应位置去就好了。

不难看出这部分复杂度为 \(nln(n)\) 。

然后再把 \(Ans\) 给 \(exp\) 回去就好了。

Code

#include <cmath>
#include <cstdio>
#include <iostream>
using namespace std;
#define LL long long
#define inline __inline__ __attribute__((always_inline))
inline LL read() {
  LL x = 0, w = 1;
  char ch = getchar();
  while (!isdigit(ch)) {
    if (ch == '-') w = -1;
    ch = getchar();
  }
  while (isdigit(ch)) {
    x = (x << 3) + (x << 1) + ch - '0';
    ch = getchar();
  }
  return x * w;
}

const int Max_n = 6e5 + 5, mod = 998244353, G = 3, Len = 1 << 19;

int ksm(int a, int b = mod - 2) {
  int res = 1;
  for (; b; b >>= 1, a = 1ll * a * a % mod)
    if (b & 1) res = 1ll * res * a % mod;
  return res;
}

namespace Poly {
int bit, len, rev[Max_n];
int gp[Max_n], invn[Max_n];
struct poly {
  int f[Max_n];
  inline int& operator[](int x) {
    return f[x];
  }
  void Init() {
    for (int i = 0; i < Max_n; i++) f[i] = 0;
  }
  void dft(int t) {
    for (int i = 0; i < len; i++)
      if (rev[i] > i) swap(f[i], f[rev[i]]);
    for (int l = 1; l < len; l <<= 1) {
      int Wn = Len / (l << 1);
      for (int i = 0; i < len; i += l << 1) {
        for (int k = i; k < i + l; k++) {
          int x = f[k], p = (k - i) * Wn;
          int y = 1ll * f[k + l] * (t == -1 ? gp[Len - p] : gp[p]) % mod;
          f[k] = (x + y) % mod, f[k + l] = (x - y + mod) % mod;
        }
      }
    }
    if (t == -1)
      for (int i = 0, Inv = ksm(len); i < len; i++) 
        f[i] = 1ll * f[i] * Inv % mod;
  }
};
void init(int n) {
  len = 1 << (bit = log2(n) + 1);
  for (int i = 0; i < len; i++)
    rev[i] = rev[i >> 1] >> 1 | ((i & 1) << bit - 1);
}
void inv(int n, poly &f, poly &g) {
  static poly F;
  F.Init(), g.Init();
  g[0] = ksm(f[0]);
  for (int deg = 2; deg < (n << 1); deg <<= 1) {
    init(deg * 3);
    for (int i = 0; i < deg; i++) F[i] = f[i];
    for (int i = deg; i < len; i++) F[i] = 0;
    g.dft(1), F.dft(1);
    for (int i = 0; i < len; i++)
      g[i] = 1ll * g[i] * (2 - 1ll * g[i] * F[i] % mod + mod) % mod;
    g.dft(-1);
    for (int i = deg; i < len; i++) g[i] = 0;
  }
}
void ln(int n, poly &f, poly &g) {
  static poly df, Inv;
  df.Init();
  for (int i = 0; i < n - 1; i++) df[i] = 1ll * f[i + 1] * (i + 1) % mod;
  inv(n, f, Inv), g.Init(), init(n * 2);
  df.dft(1), Inv.dft(1);
  for (int i = 0; i < len; i++) g[i] = 1ll * df[i] * Inv[i] % mod;
  g.dft(-1);
  for (int i = n - 1; i; i--) g[i] = 1ll * g[i - 1] * invn[i] % mod;
  g[0] = 0;
}
void exp(int n, poly &f, poly &g) {
  static poly Ln, F;
  g.Init(), F.Init();
  g[0] = 1;
  for (int deg = 2; deg < (n << 1); deg <<= 1) {
    ln(deg, g, Ln), init(deg);
    for (int i = 0; i < deg; i++) F[i] = f[i];
    for (int i = deg; i < len; i++) F[i] = 0;
    g.dft(1), F.dft(1), Ln.dft(1);
    for (int i = 0; i < len; i++) 
      g[i] = 1ll * g[i] * (1 - Ln[i] + F[i] + mod) % mod;
    g.dft(-1);
    for (int i = deg; i < len; i++) g[i] = 0;
  }
}
}
using namespace Poly;

int n, m, exi[Max_n];
poly a, ans;

namespace Input {
void main() {
  n = read(), m = read();
  for (int i = 0; i < n; i++) exi[read()]++;
}
}  // namespace Input

namespace Init {
void main() {
  for (int i = 1; i <= n; i++)
    if (exi[i])
      for (int j = 1; j * i <= m; j++)
        (a[j * i] += 1ll * exi[i] * ksm(j) % mod) %= mod;
  invn[0] = invn[1] = gp[0] = 1;
  int g = ksm(G, (mod - 1) / Len);
  for (int i = 1; i <= Len; i++) gp[i] = 1ll * gp[i - 1] * g % mod;
  for (int i = 2; i < Len; i++) invn[i] = 1ll * (mod - mod / i) * invn[mod % i] % mod;
}
}  // namespace Init

namespace Solve {
void main() {
  exp(m + 1, a, ans);
  for (int i = 1; i <= m; i++) printf("%d\n", ans[i]);
}
}  // namespace Solve

int main() {
#ifndef ONLINE_JUDGE
  freopen("4389.in", "r", stdin);
  freopen("4389.out", "w", stdout);
#endif
  Input::main();
  Init::main();
  Solve::main();
}