Power Projection Of Set Power Series
(set/power-projection-of-set-power-series.hpp)
- View this file on GitHub
- Last update: 2025-10-31 21:40:36+09:00
- Include:
#include "set/power-projection-of-set-power-series.hpp"
$a$ を集合冪級数とする. $\sum_{s}w_s[x^s]a^i$ を $i=0,\dots,M-1$ について列挙する.
$O(N^22^N+NM)$ 時間.
- [https://maspypy.com/%e9%9b%86%e5%90%88%e3%81%b9%e3%81%8d%e7%b4%9a%e6%95%b0%e9%96%a2%e9%80%a3-3-%e5%a4%9a%e9%a0%85%e5%bc%8f%e3%81%a8%e3%81%ae%e5%90%88%e6%88%90%e3%81%ae%e8%bb%a2%e7%bd%ae]
導入
$a$ を fix したとき,合成 $f(a)$ の計算は $\sum_{i}f_i[x^s]a^i$ の $s$ についての列挙.
この転置問題は $\sum_{s}w_s[x^s]a^i$ の $i$ についての列挙.
多項式 $f$ との合成は以下の要領で計算していた.
- $g(x)=f(x+c)$ とする
- $g$ を EGF にする
- $g$ と $a-c$ を合成する
$g(x)=f(x+c)$ のとき $g_i=\sum_{j\geq i}\binom{j}{i}c^{j-i}f_j$ である.この転置は難しくない.
また $g$ を EGF にするには $x^k$ の係数を $k!$ 倍すればよいので,特に転置は同じものである.
$[x^0]a=0$ のときの EGF との合成の転置を考える. subset convolution $a\star b$ は $b$ を固定すれば $a$ について線形写像. その転置アルゴリズムは $a\star b=\operatorname{rev}(\operatorname{rev}(a)\star b)$ である.
Depends on
Verified with
Code
#pragma once
#include "set/subset-convolution.hpp"
namespace SetPowerSeries {
template <class mint, int sz = 21>
vector<mint> power_projection(const vector<mint>& a, const vector<mint>& w, int m) {
static SubsetConvolution<mint, sz> sc;
int l = __builtin_ctz(a.size());
assert(a.size() == (1 << l));
mint c = a[0];
vector<vector<mint>> g(l + 1);
g[0] = w;
for (int i = 1; i <= l; i++) g[i] = vector<mint>(1 << (l - i), 0);
for (int k = l - 1; k >= 0; k--) {
auto p = sc.lift(vector<mint>(a.begin() + (1 << k), a.begin() + (2 << k)));
sc.ranked_zeta(p);
for (int i = 0; i < l - k; i++) {
auto h = vector<mint>(g[i].begin() + (1 << k), g[i].begin() + (2 << k));
reverse(h.begin(), h.end());
auto q = sc.lift(h);
sc.ranked_zeta(q);
sc.ranked_mul(q, p);
sc.ranked_mobius(q);
h = sc.unlift(q);
reverse(h.begin(), h.end());
for (int j = 0; j < h.size(); j++) g[i + 1][j] += h[j];
}
}
vector<mint> f1(l + 1, 1);
for (int i = 1; i <= l; i++) f1[i] = f1[i - 1] * i;
for (int i = 0; i <= l; i++) f1[i] *= g[i][0];
vector<mint> f(m, 0), binom(l + 1, 0);
binom[0] = 1;
mint r0 = 1;
for (int i = 0; i < m; i++) {
if (i > l) r0 *= c;
mint r = r0;
for (int j = min(i, l); j >= 0; j--, r *= c) f[i] += f1[j] * binom[j] * r;
for (int j = l; j > 0; j--) binom[j] += binom[j - 1];
}
return f;
}
}; // namespace SetPowerSeries
/**
* @brief Power Projection Of Set Power Series
* @docs docs/set/power-projection-of-set-power-series.md
*/#line 2 "set/subset-convolution.hpp"
template <class mint, int n_>
struct SubsetConvolution {
static constexpr int n = n_;
using poly = array<mint, n_ + 1>;
vector<int> pc;
SubsetConvolution() {
pc.assign(1 << n, 0);
for (int i = 1; i < pc.size(); i++) pc[i] = pc[i >> 1] + (i & 1);
}
void poly_add(poly& p, const poly& q, int d) {
for (int i = 0; i < d; i++) p[i] += q[i];
}
void poly_sub(poly& p, const poly& q, int d) {
for (int i = d; i <= n; i++) p[i] -= q[i];
}
void poly_mul(poly& p, const poly& q) {
poly r{};
for (int i = 0; i <= n; i++)
for (int j = 0; j <= n - i; j++)
r[i + j] += p[i] * q[j];
swap(p, r);
}
vector<poly> lift(const vector<mint>& a) {
int n = a.size();
assert(n == (n & -n));
vector<poly> b(n);
for (int i = 0; i < n; i++) {
b[i].fill(0);
b[i][pc[i]] = a[i];
}
return b;
}
vector<mint> unlift(const vector<poly>& b) {
int n = b.size();
assert(n == (n & -n));
vector<mint> a(n);
for (int i = 0; i < n; i++) a[i] = b[i][pc[i]];
return a;
}
void ranked_zeta(vector<poly>& a) {
int n = a.size();
for (int i = 1; i < n; i <<= 1)
for (int j = 0; j < n; j += i * 2)
for (int k = 0; k < i; k++)
poly_add(a[i + j + k], a[j + k], pc[i + j + k]);
}
void ranked_mobius(vector<poly>& a) {
int n = a.size();
for (int i = 1; i < n; i <<= 1)
for (int j = 0; j < n; j += i * 2)
for (int k = 0; k < i; k++)
poly_sub(a[i + j + k], a[j + k], pc[i + j + k]);
}
void ranked_mul(vector<poly>& a, const vector<poly>& b) {
for (int i = 0; i < a.size(); i++) poly_mul(a[i], b[i]);
}
vector<mint> multiply(const vector<mint>& a, const vector<mint>& b) {
auto p = lift(a);
auto q = lift(b);
ranked_zeta(p);
ranked_zeta(q);
ranked_mul(p, q);
ranked_mobius(p);
return unlift(p);
}
};
/**
* @brief Subset Convolution
* @docs docs/set/subset-convolution.md
*/
#line 3 "set/power-projection-of-set-power-series.hpp"
namespace SetPowerSeries {
template <class mint, int sz = 21>
vector<mint> power_projection(const vector<mint>& a, const vector<mint>& w, int m) {
static SubsetConvolution<mint, sz> sc;
int l = __builtin_ctz(a.size());
assert(a.size() == (1 << l));
mint c = a[0];
vector<vector<mint>> g(l + 1);
g[0] = w;
for (int i = 1; i <= l; i++) g[i] = vector<mint>(1 << (l - i), 0);
for (int k = l - 1; k >= 0; k--) {
auto p = sc.lift(vector<mint>(a.begin() + (1 << k), a.begin() + (2 << k)));
sc.ranked_zeta(p);
for (int i = 0; i < l - k; i++) {
auto h = vector<mint>(g[i].begin() + (1 << k), g[i].begin() + (2 << k));
reverse(h.begin(), h.end());
auto q = sc.lift(h);
sc.ranked_zeta(q);
sc.ranked_mul(q, p);
sc.ranked_mobius(q);
h = sc.unlift(q);
reverse(h.begin(), h.end());
for (int j = 0; j < h.size(); j++) g[i + 1][j] += h[j];
}
}
vector<mint> f1(l + 1, 1);
for (int i = 1; i <= l; i++) f1[i] = f1[i - 1] * i;
for (int i = 0; i <= l; i++) f1[i] *= g[i][0];
vector<mint> f(m, 0), binom(l + 1, 0);
binom[0] = 1;
mint r0 = 1;
for (int i = 0; i < m; i++) {
if (i > l) r0 *= c;
mint r = r0;
for (int j = min(i, l); j >= 0; j--, r *= c) f[i] += f1[j] * binom[j] * r;
for (int j = l; j > 0; j--) binom[j] += binom[j - 1];
}
return f;
}
}; // namespace SetPowerSeries
/**
* @brief Power Projection Of Set Power Series
* @docs docs/set/power-projection-of-set-power-series.md
*/