[洛谷 P4389]付公主的背包

danihao123 posted @ 2018年9月28日 21:42 in 题解 with tags 洛谷 NTT FFT 多项式牛顿迭代多项式exp 多项式ln , 1994 阅读

转载请注明出处：http://danihao123.is-programmer.com/

付公主有一个大小为\(10^5\)的背包。然你有\(n\)种类型的物品，其中第\(i\)种体积为\(V_i\)（是正整数），数量有\(10^5\)件。然后给出一正整数\(m\)，对任意\(i=1\ldots m\)求出包里恰好装了\(i\)体积的物品的方案数，答案对\(998244353\)取模。

\(1\leq n\leq 10^5, m\leq 10^5, V_i\leq m\)。

首先，由于每种物品都有\(10^5\)件，我们可以认为每种物品都是无限件的。

然后很显然答案的生成函数就是：\(F(x)=\prod_{i = 1}^n\frac{1}{1-x^{V_i}}\)。

但这玩意暴力NTT算一定没前途，用分治NTT搞复杂度也不对（因为\(V_i\)的和可能非常大）。然后我们对单个物品的生成函数取一下对数，观察到：

\[\ln(\frac{1}{1 - x^c}) = -\ln(1 - x^c) = \sum_{i = 1}^{+\infty}\frac{x^i}{i}\]

最后一步是麦克劳林展开搞出来的……

然后对数可以把加法拆成乘法，然后\(\ln F(x)\)就是每一种物品的对数加起来，把相同体积的和一块然后就是一个经典的调和数复杂度的东西了……

但是我们只知道了\(\ln F(x)\)，我们最后要用的是\(F(x)\)，这样再求个指数函数就行啦~

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <functional>
#include <utility>
typedef long long ll;
const ll ha = 998244353LL;
ll pow_mod(ll a, ll b) {
  ll ans = 1, res = a;
  while(b) {
    if(1LL & b) ans = ans * res % ha;
    res = res * res % ha; b >>= 1;
  }
  return ans;
}
ll inv(ll x) {
  return pow_mod(x, ha - 2LL);
}
 
const int maxn = 400050;
ll invp[maxn];
void process() {
  invp[1] = 1;
  for(int i = 2; i < maxn; i ++) {
    invp[i] = (ha - (ha / (ll)i));
    invp[i] = invp[i] * invp[ha % (ll)i] % ha;
  }
}
 
int flip(int bi, int x) {
  int ans = 0;
  for(int j = 0; j < bi; j ++) {
    if((1 << j) & x) {
      ans |= (1 << (bi - j - 1));
    }
  }
  return ans;
}
void NTT(ll *A, int bi, bool flag = false) {
  int n = 1 << bi;
  for(int i = 0; i < n; i ++) {
    int v = flip(bi, i);
    if(v < i) std::swap(A[v], A[i]);
  }
  for(int L = 1; L < n; L <<= 1) {
    ll xi_n = pow_mod(3LL, (ha - 1LL) / (ll(L << 1)));
    if(flag) xi_n = inv(xi_n);
    for(int i = 0; i < n; i += (L << 1)) {
      ll xi = 1LL;
      for(int j = i; j < i + L; j ++) {
        ll x = A[j], y = A[j + L];
        ll tmp = y * xi % ha;
        A[j] = (x + tmp) % ha;
        A[j + L] = (x - tmp + ha) % ha;
        xi = xi * xi_n % ha;
      }
    }
  }
  if(flag) {
    ll inv_n = inv(n);
    for(int i = 0; i < n; i ++) {
      A[i] = A[i] * inv_n % ha;
    }
  }
}
 
void poly_inv(int mod, ll *B, ll *BB) {
  if(mod == 1) {
    BB[0] = inv(B[0]);
  } else {
    poly_inv((mod + 1) >> 1, B, BB);
    int bi = 0, sz = 1;
    while(sz <= ((mod * 2) + 1)) {
      bi ++; sz <<= 1;
    }
    static ll tmp[maxn];
    std::copy(B, B + mod, tmp);
    std::fill(tmp + mod, tmp + sz, 0LL);
    NTT(tmp, bi); NTT(BB, bi);
    for(int i = 0; i < sz; i ++) {
      tmp[i] = (tmp[i] * BB[i]) % ha;
      tmp[i] = (tmp[i] * (ha - 1LL)) % ha;
      tmp[i] = (tmp[i] + 2LL) % ha;
      tmp[i] = (tmp[i] * BB[i]) % ha;
    }
    NTT(tmp, bi, true);
    std::copy(tmp, tmp + mod, BB);
    std::fill(BB + mod, BB + sz, 0LL);
  }
}
void poly_ln(int mod, ll *B, ll *BB) {
  ll *tmp = (ll*)calloc(maxn, sizeof(ll));
  poly_inv(mod, B, BB);
  for(int i = 0; i < mod - 1; i ++) tmp[i] = B[i + 1] * (ll(i + 1)) % ha;
   
  int bi = 0, sz = 1;
  while(sz <= ((mod * 2) + 1)) {
    bi ++; sz <<= 1;
  }
  NTT(tmp, bi); NTT(BB, bi);
  for(int i = 0; i < sz; i ++) BB[i] = BB[i] * tmp[i] % ha;
  NTT(BB, bi, true);
  for(int i = mod - 1; i > 0; i --) BB[i] = BB[i - 1] * invp[i] % ha;
  BB[0] = 0; std::fill(BB + mod, BB + sz, 0LL);
  free(tmp);
}
void poly_exp(int mod, ll *B, ll *BB) {
  if(mod == 1) {
    BB[0] = 1;
  } else {
    poly_exp((mod + 1) >> 1, B, BB);
    int bi = 0, sz = 1;
    while(sz <= ((mod * 2) + 1)) {
      bi ++; sz <<= 1;
    }
     
    std::fill(BB + mod, BB + sz, 0LL);
    ll *tmp = (ll*)calloc(maxn, sizeof(ll));
    poly_ln(mod, BB, tmp);
    for(int i = 0; i < mod; i ++) {
      tmp[i] = (B[i] - tmp[i] + ha) % ha;
    }
    tmp[0] = (tmp[0] + 1LL) % ha;
    NTT(tmp, bi); NTT(BB, bi);
    for(int i = 0; i < sz; i ++) BB[i] = tmp[i] * BB[i] % ha;
    NTT(BB, bi, true);
     
    std::fill(BB + mod, BB + sz, 0LL);
    free(tmp);
  }
}
 
ll A[maxn], B[maxn];
int cnt[maxn];
int main() {
  process();
  int n, m; scanf("%d%d", &n, &m);
  for(int i = 1; i <= n; i ++) {
    int x; scanf("%d", &x);
    cnt[x] ++;
  }
  int N = 100000;
  for(int i = 1; i <= N; i ++) {
    if(!cnt[i]) continue;
    for(int j = 1; j * i <= N; j ++) {
      ll delta = (ll)cnt[i] * invp[j] % ha;
      A[j * i] = (A[j * i] + delta) % ha;
    }
  }
  poly_exp(N + 1, A, B);
  for(int i = 1; i <= m; i ++) {
    printf("%lld\n", B[i]);
  }
  return 0;
}