题目链接:Luogu 4389
付公主有一个可爱的背包,这个背包最多可以装 $10 ^ 5$ 大小的东西。付公主有 $n$ 种商品,她要准备出摊了。每种商品体积为 $V_i$,都有 $10 ^ 5$ 件。
给定 $m$,对于整数 $s \in [1, m]$,请你回答用这些商品恰好装 $s$ 体积的方案数。
数据范围:$1 \le n \le 10 ^ 5$,$1 \le V_i \le m \le 10 ^ 5$。
Solution
显然这就是一个完全背包问题(物品数量极多以至于退化为完全背包)。我们对每种物品建立生成函数,体积为 $v$ 的物品的生成函数为:
$$ F_v(x) = \sum_{i = 0} ^ {\infty} [i \bmod v = 0] x ^ i = \frac{1}{1 - x ^ v} $$
如果我们直接把 $n$ 个生成函数乘起来,复杂度 $\mathcal O(nm \log m)$ 是无法接受的。
考虑求出 $\ln F_v(x)$ 后相加,最后再 $\exp$ 回去。
$$ \ln F_v(x) = G(x) $$
求导得:
$$ \begin{align*} G'(x) = \frac{F'_v(x)}{F_v(x)} & = (1 - x ^ v) \sum_{i = 1} ^ {\infty} iv \cdot x ^ {iv - 1} \\ & =\sum_{i = 1} ^ {\infty} iv \cdot x ^ {iv - 1} - \sum_{i = 1} ^ {\infty} iv \cdot x ^ {(i + 1) v - 1} \\ & = \sum_{i = 1} ^ {\infty} iv \cdot x ^ {iv - 1} - \sum_{i = 1} ^ {\infty} (i - 1)v \cdot x ^ {iv - 1} \\ & = \sum_{i = 1} ^ {\infty} v \cdot x ^ {iv - 1} \end{align*} $$
积分得到:
$$ \begin{align*} G(x) & = \int \left(\sum_{i = 1} ^ {\infty} v \cdot x ^ {iv - 1}\right) \mathrm{d} x \\ & = \sum_{i = 1} ^ {\infty} \frac{1}{i} x ^ {iv} \end{align*} $$
我们可以 $\mathcal O(m \log m)$ 直接求出 $G(x)$,最后再对 $G$ 求 $\exp$ 即可!
时间复杂度:$\mathcal O(m \log m)$。
Code
#include <cstdio>
#include <cmath>
#include <cassert>
#include <algorithm>
#include <vector>
#include <unordered_map>
typedef std::vector<int> Vec;
const int MOD = 998244353, G = 3;
void add(int &x, int y) {
(x += y) >= MOD && (x -= MOD);
}
void sub(int &x, int y) {
(x -= y) < 0 && (x += MOD);
}
int add(int x) {
return x >= MOD ? x - MOD : x;
}
int sub(int x) {
return x < 0 ? x + MOD : x;
}
void mod(int &x) {
x >= MOD && (x -= MOD), x < 0 && (x += MOD);
}
int mul(int x, int y) {
return 1LL * x * y % MOD;
}
int pow(int x, int k) {
int ans = 1;
for (; k; k >>= 1, x = 1LL * x * x % MOD) {
if (k & 1) ans = 1LL * ans * x % MOD;
}
return ans;
}
int inv(int x) {
return pow(x, MOD - 2);
}
int BSGS(int a, int b, int p) {
std::unordered_map<int, int> mp;
int m = ceil(sqrt(p));
for (int i = 0; i <= m; b = 1LL * b * a % p, i++) mp[b] = i;
a = pow(a, m);
for (int i = 0, j = 1; i < m; j = 1LL * j * a % p, i++) {
if (mp.count(j) && i * m >= mp[j]) {
return i * m - mp[j];
}
}
return -1;
}
int degree(int a, int k, int p) {
int x = BSGS(G, a, p);
assert(x >= 0 && x % k == 0);
int r = pow(G, x / k);
return std::min(r, p - r);
}
namespace FFT {
int extend(int x);
void NTT(Vec &A, bool opt);
void DFT(Vec &A);
void IDFT(Vec &A);
int extend(int x) {
int n = 1;
for (; n < x; n <<= 1);
return n;
}
void NTT(Vec &A, bool opt) {
int n = A.size(), k = 0;
for (; (1 << k) < n; k++);
Vec rev(n);
for (int i = 0; i < n; i++) {
rev[i] = rev[i >> 1] >> 1 | (i & 1) << (k - 1);
if (i < rev[i]) std::swap(A[i], A[rev[i]]);
}
for (int l = 2; l <= n; l <<= 1) {
int m = l >> 1, w = pow(G, (MOD - 1) / l);
if (opt) w = inv(w);
for (int j = 0; j < n; j += l) {
int wk = 1;
for (int i = 0; i < m; i++, wk = 1LL * wk * w % MOD) {
int p = A[i + j], q = 1LL * wk * A[i + j + m] % MOD;
A[i + j] = (p + q) % MOD;
A[i + j + m] = (p - q + MOD) % MOD;
}
}
}
}
void DFT(Vec &A) {
NTT(A, false);
}
void IDFT(Vec &A) {
NTT(A, true);
int t = inv(A.size());
for (auto &x : A) x = 1LL * x * t % MOD;
}
}
using namespace FFT;
namespace Poly {
Vec operator + (Vec A, Vec B);
Vec operator + (Vec A, int v);
Vec operator + (int v, Vec A);
Vec operator - (Vec A, Vec B);
Vec operator - (Vec A, int v);
Vec operator - (int v, Vec A);
Vec operator - (Vec A);
Vec operator * (Vec A, Vec B);
Vec operator * (Vec A, int v);
Vec operator * (int v, Vec A);
Vec operator / (Vec A, Vec B);
Vec operator / (Vec A, int v);
Vec operator % (Vec A, Vec B);
Vec operator ~ (Vec A);
Vec operator ^ (Vec A, int k);
Vec operator << (Vec A, int x);
Vec operator >> (Vec A, int x);
Vec fix(Vec A, int n);
Vec der(Vec A, bool mod);
Vec inte(Vec A, bool mod);
Vec sqrt(Vec A);
Vec root(Vec A, int k);
Vec ln(Vec A);
Vec exp(Vec A);
Vec sin(Vec A);
Vec cos(Vec A);
Vec tan(Vec A);
Vec asin(Vec A);
Vec acos(Vec A);
Vec atan(Vec A);
void print(Vec A);
Vec operator + (Vec A, Vec B) {
int n = std::max(A.size(), B.size());
A.resize(n), B.resize(n);
for (int i = 0; i < n; i++) add(A[i], B[i]);
return A;
}
Vec operator + (Vec A, int v) {
add(A[0], v);
return A;
}
Vec operator + (int v, Vec A) {
add(A[0], v);
return A;
}
Vec operator - (Vec A, Vec B) {
int n = std::max(A.size(), B.size());
A.resize(n), B.resize(n);
for (int i = 0; i < n; i++) sub(A[i], B[i]);
return A;
}
Vec operator - (Vec A, int v) {
sub(A[0], v);
return A;
}
Vec operator - (int v, Vec A) {
A = -A, add(A[0], v);
return A;
}
Vec operator - (Vec A) {
for (auto &x : A) x = sub(-x);
return A;
}
Vec operator * (Vec A, Vec B) {
int n = A.size() + B.size() - 1, N = extend(n);
A.resize(N), DFT(A);
B.resize(N), DFT(B);
for (int i = 0; i < N; i++) A[i] = 1LL * A[i] * B[i] % MOD;
IDFT(A), A.resize(n);
return A;
}
Vec operator * (Vec A, int v) {
for (auto &x : A) x = 1LL * x * v % MOD;
return A;
}
Vec operator * (int v, Vec A) {
for (auto &x : A) x = 1LL * x * v % MOD;
return A;
}
Vec operator / (Vec A, Vec B) {
int n = A.size() - B.size() + 1;
if (n <= 0) return Vec(1, 0);
std::reverse(A.begin(), A.end());
std::reverse(B.begin(), B.end());
A.resize(n), B.resize(n);
A = fix(A * ~B, n);
std::reverse(A.begin(), A.end());
return A;
}
Vec operator / (Vec A, int v) {
return A * inv(v);
}
Vec operator % (Vec A, Vec B) {
int n = B.size() - 1;
return fix(A - A / B * B, n);
}
Vec operator ~ (Vec A) {
int n = A.size(), N = extend(n);
A.resize(N);
Vec I(N, 0);
I[0] = inv(A[0]);
for (int l = 2; l <= N; l <<= 1) {
Vec P(l), Q(l);
std::copy(A.begin(), A.begin() + l, P.begin());
std::copy(I.begin(), I.begin() + l, Q.begin());
int L = l << 1;
P.resize(L), DFT(P);
Q.resize(L), DFT(Q);
for (int i = 0; i < L; i++) {
P[i] = 1LL * Q[i] * (2 - 1LL * P[i] * Q[i] % MOD + MOD) % MOD;
}
IDFT(P), P.resize(l);
std::copy(P.begin(), P.begin() + l, I.begin());
}
I.resize(n);
return I;
}
Vec operator ^ (Vec A, int k) {
int n = A.size(), x = 0;
for (; x < n && A[x] == 0; x++);
if (1LL * x * k >= n) return Vec(n, 0);
A = A >> x;
int v = A[0];
return (exp(ln(A) * k) * pow(v, k)) << (x * k);
}
Vec operator << (Vec A, int x) {
int n = A.size();
Vec B(n, 0);
for (int i = 0; i < n - x; i++) B[i + x] = A[i];
return B;
}
Vec operator >> (Vec A, int x) {
int n = A.size();
Vec B(n, 0);
for (int i = 0; i < n - x; i++) B[i] = A[i + x];
return B;
}
Vec fix(Vec A, int n) {
A.resize(n);
return A;
}
Vec der(Vec A, bool mod = true) {
int n = A.size();
if (n == 1) return Vec(1, 0);
Vec D(n - 1, 0);
for (int i = 1; i < n; i++) D[i - 1] = 1LL * i * A[i] % MOD;
if (mod) D.resize(n);
return D;
}
Vec inte(Vec A, bool mod = true) {
int n = A.size();
Vec I(n + 1, 0);
for (int i = 1; i <= n; i++) I[i] = 1LL * inv(i) * A[i - 1] % MOD;
if (mod) I.resize(n);
return I;
}
Vec sqrt(Vec A) {
int n = A.size(), N = extend(n);
A.resize(N);
Vec R(N, 0);
R[0] = degree(A[0], 2, MOD);
int i2 = inv(2);
for (int l = 2; l <= N; l <<= 1) {
Vec P(l), Q(l);
std::copy(A.begin(), A.begin() + l, P.begin());
std::copy(R.begin(), R.begin() + l, Q.begin());
Vec I = ~Q;
int L = l << 1;
P.resize(L), DFT(P);
Q.resize(L), DFT(Q);
I.resize(L), DFT(I);
for (int i = 0; i < L; i++) {
Q[i] = 1LL * Q[i] * Q[i] % MOD;
P[i] = 1LL * (P[i] + Q[i]) * i2 % MOD * I[i] % MOD;
}
IDFT(P), P.resize(l);
std::copy(P.begin(), P.begin() + l, R.begin());
}
R.resize(n);
return R;
}
Vec root(Vec A, int k) {
return k == 1 ? A : k == 2 ? sqrt(A) : exp(ln(A) / k);
}
Vec ln(Vec A) {
assert(A[0] == 1);
int n = A.size();
return inte(fix(der(A) * ~A, n));
}
Vec exp(Vec A) {
assert(A[0] == 0);
int n = A.size(), N = extend(n);
A.resize(N);
Vec E(N, 0);
E[0] = 1;
for (int l = 2; l <= N; l <<= 1) {
Vec P = (-ln(fix(E, l)) + fix(A, l) + 1) * fix(E, l);
std::copy(P.begin(), P.begin() + l, E.begin());
}
E.resize(n);
return E;
}
Vec sin(Vec A) {
assert(A[0] == 0);
int i = degree(MOD - 1, 2, MOD);
Vec E = exp(i * A);
return (E - ~E) / (2LL * i % MOD);
}
Vec cos(Vec A) {
assert(A[0] == 0);
int i = degree(MOD - 1, 2, MOD);
Vec E = exp(i * A);
return (E + ~E) / 2;
}
Vec tan(Vec A) {
assert(A[0] == 0);
int n = A.size();
return fix(sin(A) * ~cos(A), n);
}
Vec asin(Vec A) {
assert(A[0] == 0);
int n = A.size();
return inte(fix(der(A) * ~sqrt(1 - fix(A * A, n)), n));
}
Vec acos(Vec A) {
assert(A[0] == 0);
return -asin(A);
}
Vec atan(Vec A) {
assert(A[0] == 0);
int n = A.size();
return inte(fix(der(A) * ~(1 + fix(A * A , n)), n));
}
void print(Vec A) {
int n = A.size();
for (int i = 0; i < n; i++) printf("%d%c", A[i], " \n"[i == n - 1]);
}
}
using namespace Poly;
int main() {
int n, m;
scanf("%d%d", &n, &m);
Vec cnt(m + 1, 0), A(m + 1, 0), inv(m + 1, 0);
for (int i = 1; i <= n; i++) {
int x;
scanf("%d", &x);
cnt[x]++;
}
inv[1] = 1;
for (int i = 2; i <= m; i++) {
inv[i] = 1LL * (MOD - MOD / i) * inv[MOD % i] % MOD;
}
for (int i = 1; i <= m; i++) {
if (cnt[i] == 0) continue;
for (int j = i; j <= m; j += i) {
add(A[j], 1LL * cnt[i] * inv[j / i] % MOD);
}
}
A = exp(A);
for (int i = 1; i <= m; i++) printf("%d\n", A[i]);
return 0;
}