2. 程式人生 > 其它 >多項式複合和拉格朗日反演 學習筆記

多項式複合和拉格朗日反演 學習筆記

阿新 • • 發佈：2022-05-21

20220520

一些小東西

卡特蘭數和去掉根的二叉樹對應。考慮卡特蘭數的遞推公式 \(f(n) = \sum_{i = 0}^{n - 1} f(i) f(n - 1 - i)\)，並有 \(f(0) = f(1) = 1\)。

多項式複合

複合運算

冪級數的複合運算。

設 \(A(w)=\sum_{i\ge 0} a_i w^i\)，\(B(x)=\sum_{i \ge 1}b_i x ^ i\)。則 \(B(x)\) 與 \(A(w)\) 的複合就是

\[C(x) = (A\circ B)(x) = A(B(x))=\sum_{i \ge 0} a_i (B(x)) ^ i \]

可以將其整理為 \(c_i = \sum_{i \ge 0} c_i x ^ i\)

的形式。

由於假定 \(b_0 = 0\)，於是即使 \(A\) 有無限多個非 \(0\) 項， \(B(x)\) 和 \(A(w)\) 的複合可以定義。

多項式複合函式

給定 \(f(w), g(x)\)，如果直接用定義計算 \(f(g(x))\) 的前 \(n\) 項係數，時間複雜度為一個常數大的 \(O(n ^ 2\log n)\) 。

下面是一個複雜度為 \(O(N ^ 2)\) 的演算法，運用了大步小步的演算法思想。

令 \(d = \lceil \sqrt n\rceil\)，先求出來 \(g(x) ^ 2, g(x) ^ 3, \ldots , g(x) ^ d\)，複雜度是 \(O(d n\log n)\)

，同時計算 \(g(x) ^ {2d},g(x)^{3d},\ldots ,g(x) ^ {d(d - 1)}\)，複雜度為 \(O(d n\log n)\)。

接下來，根據

\[f(x) = \sum_{i \ge 0}f_i (g(x) ^ i) = \sum_{0 \le i < d}g(x) ^ {id}\sum_{0\le j <d} f_{id + j} g(x)^j \]

然後發現，如果我暴力計算 \(\sum_{0\le j <d} f_{id + j} g(x)^j\) 然後再暴力把前面的式子加起來，一共的多項式加法次數不會超過 \(O(N)\)，於是這部分複雜度為 \(O(N ^ 2)\)

。

如是說，我們在 \(O(N ^ 2 + n \sqrt n \log n)\) 的情況下計算了多項式複合逆。

題目：P5373 【模板】多項式複合函式 - 洛谷 | 電腦科學教育新生態 (luogu.com.cn)

// Siriqwq
#include <bits/stdc++.h>
using std::cin;
using std::cout;
using std::vector;
using std::copy;
using std::reverse;
using std::sort;
using std::get;
using std::unique;
using std::swap;
using std::array;
using std::cerr;
using std::function;
using std::map;
using std::set;
using std::pair;
using std::mt19937;
using std::make_pair;
using std::tuple;
using std::make_tuple;
using std::uniform_int_distribution;
using ll = long long;
namespace qwq {
	mt19937 eng;
	void init(int Seed) {return eng.seed(Seed);}
	int rnd(int l = 1, int r = 1000000000) {return uniform_int_distribution<int> (l, r)(eng);}
}
template<typename T>
inline void chkmin(T &x, T y) {if (x > y) x = y;}
template<typename T>
inline void chkmax(T &x, T y) {if (x < y) x = y;}
template<typename T>
inline T min(const T &x, const T &y) {return x < y ? x : y;}
template<typename T>
inline T max(const T &x, const T &y) {return x > y ? x : y;}
char buf[100000], *bufs, *buft;
#define gc() ((bufs == buft && (buft = (bufs = buf) + fread(buf, 1, 100000, stdin))), bufs == buft ? -1 : *bufs++);
template<typename T>
inline void read(T &x) {
	x = 0;
	bool f = 0;
	char ch = gc();
	while (!isdigit(ch)) f = ch == '-', ch = gc();
	while (isdigit(ch)) x = x * 10 + ch - '0', ch = gc();
	if (f) x = -x;
}
inline void reads(char *s) {
	char ch = gc();
	while (isspace(ch)) ch = gc();
	while (!isspace(ch)) *(s++) = ch;
	return;
}
template<typename T, typename ...Arg>
inline void read(T &x, Arg &... y) {
	read(x);
	read(y...);
}
#define O(x) cerr << #x << " : " << x << '\n'
const double Pi = acos(-1);
const int MAXN = 262144, MOD = 998244353, inv2 = (MOD + 1) / 2, I32_INF = 0x3f3f3f3f;
const long long I64_INF = 0x3f3f3f3f3f3f3f3f;
auto Ksm = [] (int x, int y) -> int {
	if (y < 0) {
		y %= MOD - 1;
		y += MOD - 1;
	}
	int ret = 1;
	for (; y; y /= 2, x = (long long) x * x % MOD) if (y & 1) ret = (long long) ret * x % MOD;
	return ret;
};
auto Mod = [] (int x) -> int {
	if (x >= MOD) return x - MOD;
	else if (x < 0) return x + MOD;
	else return x;
};
template<const int N_num, const int M_num>
struct Graph {
	int H[N_num];
	struct Edge {int to, lac;} e[M_num];
	inline void add_edge(int x, int y) {e[*H] = {y, H[x]};H[x] = (*H)++;}
	inline void init() {memset(H, -1, sizeof H);*H = 0;}
};
#define go(x, y) for (int i = x.H[y], v; (v = x.e[i].to) && ~i; i = x.e[i].lac)
inline int ls(int k) {return k << 1;}
inline int rs(int k) {return k << 1 | 1;}
int fac[MAXN], ifac[MAXN];
namespace POLY {
int SZ, R[MAXN], W[MAXN], INV[MAXN + 1];
void POLYINIT() {
	INV[1] = 1;
	for (int i = 2; i <= MAXN; ++i) INV[i] = (long long) (MOD - MOD / i) * INV[MOD % i] % MOD;
}
void INIT(int len) {
	if (SZ == len) return;
	SZ = len;
	for (int i = 1; i < len; ++i) R[i] = (R[i >> 1] >> 1) | (i & 1 ? (len >> 1) : 0);
	int wn = Ksm(3, (MOD - 1) / len);
	W[len >> 1] = 1;
	for (int i = (len >> 1) + 1; i < len; ++i) W[i] = (long long) W[i - 1] * wn % MOD;
	for (int i = (len >> 1) - 1; i > 0; --i) W[i] = W[i << 1];
}
unsigned long long c[MAXN];
void Ntt(vector<int>& F, int limit, int type) {
	copy(F.begin(), F.begin() + limit, c);
	for (int i = 1; i < limit; ++i) if (i < R[i]) swap(c[i], c[R[i]]);
	for (int o = 2, j = 1; o <= limit; o <<= 1, j <<= 1) {
		for (int i = 0; i < limit; i += o) {
			for (int k = 0; k < j; ++k) {
				unsigned long long OI = c[i + j + k] * W[k + j] % MOD;
				c[i + j + k] = c[i + k] + MOD - OI;
				c[i + k] += OI;
			}
		}
	}
	if (type == -1) {
		reverse(c + 1, c + limit);
		int inv = INV[limit];
		for (int i = 0; i < limit; ++i) c[i] = c[i] % MOD * inv % MOD;
	}
	for (int i = 0; i < limit; ++i) F[i] = c[i] % MOD;
}
struct Poly {
	vector<int> v;
	int& operator [] (const int &pos) { return v[pos]; }
	int len() { return v.size(); }
	void set(int l) { return v.resize(l); }
	void adjust() { while (v.size() > 1 && !v.back()) v.pop_back(); }
	void rev() { reverse(v.begin(), v.end()); }
	void Ntt(int L, int type) {
		int limit = 1 << L;
		INIT(limit);
		set(limit);
		POLY::Ntt(v, limit, type);
	}
	void Squ() {
		int L = ceil(log2(len())) + 1, limit = 1 << L;
		Ntt(L, 1);
		for (int i = 0; i < limit; ++i) v[i] = (long long) v[i] * v[i] % MOD;
		Ntt(L, -1);
		adjust();
	}
	void operator += (Poly &x) {
		if (len() < x.len()) set(x.len());
		for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] + x[i]);
		adjust();
	}
	void operator -= (Poly &x) {
		if (len() < x.len()) set(x.len());
		for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] - x[i]);
		adjust();
	}
	Poly operator * (Poly &x) {
		Poly ret, tmp0 = *this, tmp1 = x;
		int L = ceil(log2(tmp0.len() + tmp1.len() - 1)), n = 1 << L;
		tmp0.Ntt(L, 1);
		tmp1.Ntt(L, 1);
		ret.set(n);
		for (int i = 0; i < n; ++i) ret[i] = (long long) tmp0[i] * tmp1[i] % MOD;
		ret.Ntt(L, -1);
		ret.adjust();
		return ret;
	}
	Poly operator - (Poly &x) {
		Poly ret;
		ret.set(max(len(), x.len()));
		for (int i = 0; i < len(); ++i) ret[i] = v[i];
		for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] - x[i]);
		return ret;
	}
	Poly operator * (int x) {
		Poly ret = *this;
		for (auto &i: ret.v) i = (ll) i * x % MOD;
		return ret;
	}
	Poly operator + (Poly &x) {
		Poly ret;
		ret.set(max(len(), x.len()));
		for (int i = 0; i < len(); ++i) ret[i] = v[i];
		for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] + x[i]);
		return ret;
	}
	void operator *= (Poly &x) {
		Poly tmp = x;
		int L = ceil(log2(len() + x.len() - 1)), n = 1 << L;
		Ntt(L, 1);
		tmp.Ntt(L, 1);
		for (int i = 0; i < n; ++i) v[i] = (long long) v[i] * tmp[i] % MOD;
		Ntt(L, -1);
		adjust();
	}
};
}
using namespace POLY;
int N, M;
Poly f, g, b1[150], b2[150];
int main() {
	// std::ios::sync_with_stdio(0);
	// cout << std::fixed << std::setprecision(8);
	// cin.tie(0);
	// cout.tie(0);
	freopen("1.in", "r", stdin);
	// freopen("1.out", "w", stdout);
	read(N);
	read(M);
	++N;
	++M;
	qwq::init(20050112);
	int lim = ceil(sqrt(N));
	POLYINIT();
	f.set(N);
	g.set(M);
	for (int i = 0; i < N; ++i) read(f[i]);
	for (int j = 0; j < M; ++j) read(g[j]);
	b1[0].set(N);
	b1[0][0] = 1;
	for (int i = 1; i < lim; ++i) b1[i] = b1[i - 1] * g, b1[i].set(N);
	Poly tmp = g * b1[lim - 1];
	tmp.set(N);
	b2[0].set(N);
	b2[0][0] = 1;
	for (int i = 1; i < lim; ++i) b2[i] = b2[i - 1] * tmp, b2[i].set(N);
	Poly ans;
	ans.set(N);
	for (int i = 0; i < lim; ++i) {
		Poly ret;
		ret.set(N);
		for (int j = 0; j < lim; ++j) {
			if (i * lim + j >= N) break;
			Poly tmp = b1[j] * f[i * lim + j];
			ret += tmp;
			// ret.set(N);
		}
		ret *= b2[i];
		ret.set(N);
		ans += ret;
	}
	ans.set(N);
	for (int i = 0; i < N; ++i) printf("%d ", ans[i]);
	puts("");
	// cout << (-3 / 2);
	cerr << ((double) clock() / CLOCKS_PER_SEC) << '\n';
	return (0-0);
}

多項式複合逆

給定 \(f(x)\)，滿足 \(f_0 = 0,f_1 \not = 0\)，求 \(g(x)\)，滿足 \(f(g(x))=x\)。

先來一個公式，就是下面的反演公式

\[[x^n] g(x) = \frac{1}{n}[w ^ {n - 1}]\left(\frac{w}{f(w)}\right)^n \]

然後，我們還是類似於上面的套路，利用大步小步思想，令 \(d=\lceil \sqrt n\rceil\)

\[\begin{aligned} g(x)&=\sum_{i = 1}^{n}\left(\frac{1}{i}\left(\frac{w}{f(w)}\right)^i[w ^ { ^{i - 1}}]\right) x ^ i \\ &=\sum_{i = 0}^{d - 1}\sum_{j = 1}^{d}[id + j \le n] \left(\left(\frac{1}{id + j} [w^{id + j - 1}]\left(\frac{w}{f(w)}\right)^{id} \left(\frac{w}{f(w)}\right)^{j}\right)\right)x ^ {id + j} \end{aligned} \]

通過預處理 \(\left(\frac{w}{f(w)}\right)\) 的若干次冪可以做到 \(O(N^2+N\sqrt N \log N)\) 的複雜度。

// Siriqwq
#include <bits/stdc++.h>
using std::cin;
using std::cout;
using std::vector;
using std::copy;
using std::reverse;
using std::sort;
using std::get;
using std::unique;
using std::swap;
using std::array;
using std::cerr;
using std::function;
using std::map;
using std::set;
using std::pair;
using std::mt19937;
using std::make_pair;
using std::tuple;
using std::make_tuple;
using std::uniform_int_distribution;
using ll = long long;
namespace qwq {
	mt19937 eng;
	void init(int Seed) {return eng.seed(Seed);}
	int rnd(int l = 1, int r = 1000000000) {return uniform_int_distribution<int> (l, r)(eng);}
}
template<typename T>
inline void chkmin(T &x, T y) {if (x > y) x = y;}
template<typename T>
inline void chkmax(T &x, T y) {if (x < y) x = y;}
template<typename T>
inline T min(const T &x, const T &y) {return x < y ? x : y;}
template<typename T>
inline T max(const T &x, const T &y) {return x > y ? x : y;}
char buf[100000], *bufs, *buft;
#define gc() ((bufs == buft && (buft = (bufs = buf) + fread(buf, 1, 100000, stdin))), bufs == buft ? -1 : *bufs++);
template<typename T>
inline void read(T &x) {
	x = 0;
	bool f = 0;
	char ch = gc();
	while (!isdigit(ch)) f = ch == '-', ch = gc();
	while (isdigit(ch)) x = x * 10 + ch - '0', ch = gc();
	if (f) x = -x;
}
inline void reads(char *s) {
	char ch = gc();
	while (isspace(ch)) ch = gc();
	while (!isspace(ch)) *(s++) = ch;
	return;
}
template<typename T, typename ...Arg>
inline void read(T &x, Arg &... y) {
	read(x);
	read(y...);
}
#define O(x) cerr << #x << " : " << x << '\n'
const double Pi = acos(-1);
const int MAXN = 262144, MOD = 998244353, inv2 = (MOD + 1) / 2, I32_INF = 0x3f3f3f3f;
const long long I64_INF = 0x3f3f3f3f3f3f3f3f;
auto Ksm = [] (int x, int y) -> int {
	if (y < 0) {
		y %= MOD - 1;
		y += MOD - 1;
	}
	int ret = 1;
	for (; y; y /= 2, x = (long long) x * x % MOD) if (y & 1) ret = (long long) ret * x % MOD;
	return ret;
};
auto Mod = [] (int x) -> int {
	if (x >= MOD) return x - MOD;
	else if (x < 0) return x + MOD;
	else return x;
};
template<const int N_num, const int M_num>
struct Graph {
	int H[N_num];
	struct Edge {int to, lac;} e[M_num];
	inline void add_edge(int x, int y) {e[*H] = {y, H[x]};H[x] = (*H)++;}
	inline void init() {memset(H, -1, sizeof H);*H = 0;}
};
#define go(x, y) for (int i = x.H[y], v; (v = x.e[i].to) && ~i; i = x.e[i].lac)
inline int ls(int k) {return k << 1;}
inline int rs(int k) {return k << 1 | 1;}
int fac[MAXN], ifac[MAXN];
namespace POLY {
int SZ, R[MAXN], W[MAXN], INV[MAXN + 1];
void POLYINIT() {
	INV[1] = 1;
	for (int i = 2; i <= MAXN; ++i) INV[i] = (long long) (MOD - MOD / i) * INV[MOD % i] % MOD;
}
void INIT(int len) {
	if (SZ == len) return;
	SZ = len;
	for (int i = 1; i < len; ++i) R[i] = (R[i >> 1] >> 1) | (i & 1 ? (len >> 1) : 0);
	int wn = Ksm(3, (MOD - 1) / len);
	W[len >> 1] = 1;
	for (int i = (len >> 1) + 1; i < len; ++i) W[i] = (long long) W[i - 1] * wn % MOD;
	for (int i = (len >> 1) - 1; i > 0; --i) W[i] = W[i << 1];
}
unsigned long long c[MAXN];
void Ntt(vector<int>& F, int limit, int type) {
	copy(F.begin(), F.begin() + limit, c);
	for (int i = 1; i < limit; ++i) if (i < R[i]) swap(c[i], c[R[i]]);
	for (int o = 2, j = 1; o <= limit; o <<= 1, j <<= 1) {
		for (int i = 0; i < limit; i += o) {
			for (int k = 0; k < j; ++k) {
				unsigned long long OI = c[i + j + k] * W[k + j] % MOD;
				c[i + j + k] = c[i + k] + MOD - OI;
				c[i + k] += OI;
			}
		}
	}
	if (type == -1) {
		reverse(c + 1, c + limit);
		int inv = INV[limit];
		for (int i = 0; i < limit; ++i) c[i] = c[i] % MOD * inv % MOD;
	}
	for (int i = 0; i < limit; ++i) F[i] = c[i] % MOD;
}
struct Poly {
	vector<int> v;
	int& operator [] (const int &pos) { return v[pos]; }
	int len() { return v.size(); }
	void set(int l) { return v.resize(l); }
	void adjust() { while (v.size() > 1 && !v.back()) v.pop_back(); }
	void rev() { reverse(v.begin(), v.end()); }
	void Ntt(int L, int type) {
		int limit = 1 << L;
		INIT(limit);
		set(limit);
		POLY::Ntt(v, limit, type);
	}
	void Squ() {
		int L = ceil(log2(len())) + 1, limit = 1 << L;
		Ntt(L, 1);
		for (int i = 0; i < limit; ++i) v[i] = (long long) v[i] * v[i] % MOD;
		Ntt(L, -1);
		adjust();
	}
	void operator += (Poly &x) {
		if (len() < x.len()) set(x.len());
		for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] + x[i]);
		adjust();
	}
	void operator -= (Poly &x) {
		if (len() < x.len()) set(x.len());
		for (int i = 0; i < x.len(); ++i) v[i] = Mod(v[i] - x[i]);
		adjust();
	}
	Poly operator * (Poly &x) {
		Poly ret, tmp0 = *this, tmp1 = x;
		int L = ceil(log2(tmp0.len() + tmp1.len() - 1)), n = 1 << L;
		tmp0.Ntt(L, 1);
		tmp1.Ntt(L, 1);
		ret.set(n);
		for (int i = 0; i < n; ++i) ret[i] = (long long) tmp0[i] * tmp1[i] % MOD;
		ret.Ntt(L, -1);
		ret.adjust();
		return ret;
	}
	Poly operator - (Poly &x) {
		Poly ret;
		ret.set(max(len(), x.len()));
		for (int i = 0; i < len(); ++i) ret[i] = v[i];
		for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] - x[i]);
		return ret;
	}
	Poly operator * (int x) {
		Poly ret = *this;
		for (auto &i: ret.v) i = (ll) i * x % MOD;
		return ret;
	}
	Poly operator + (Poly &x) {
		Poly ret;
		ret.set(max(len(), x.len()));
		for (int i = 0; i < len(); ++i) ret[i] = v[i];
		for (int i = 0; i < x.len(); ++i) ret[i] = Mod(ret[i] + x[i]);
		return ret;
	}
	void operator *= (Poly &x) {
		Poly tmp = x;
		int L = ceil(log2(len() + x.len() - 1)), n = 1 << L;
		Ntt(L, 1);
		tmp.Ntt(L, 1);
		for (int i = 0; i < n; ++i) v[i] = (long long) v[i] * tmp[i] % MOD;
		Ntt(L, -1);
		adjust();
	}
	Poly GetInv(int deg = -1) {
		if (deg == 1) return {{Ksm(v[0], MOD - 2)}};
		Poly ret = GetInv((deg + 1) / 2), tmp;
		int L = ceil(log2(deg)) + 1, n = 1 << L, mx = min(len(), deg);
		tmp.set(deg);
		for (int i = 0; i < mx; ++i) tmp[i] = v[i];
		tmp.Ntt(L, 1);
		ret.Ntt(L, 1);
		for (int i = 0; i < n; ++i) ret[i] = (2 - (long long) tmp[i] * ret[i] % MOD + MOD) * ret[i] % MOD;
		ret.Ntt(L, -1);
		ret.set(deg);
		return ret;
	}
};
}
using namespace POLY;
Poly f, g, b1[132], b2[132];
int N;
int main() {
// 	freopen("1.in", "r", stdin);
	read(N);
	POLYINIT();
	f.set(N);
	read(f[0]);
	--N;
	for (int i = 0; i < N; ++i) read(f[i]);
	int d = ceil(sqrt(N));
	g = f.GetInv(N);
	b2[0].set(N);
	b2[0][0] = 1;
	for (int i = 1; i <= d; ++i) b2[i] = b2[i - 1] * g, b2[i].set(N);
	b1[0].set(N);
	b1[0][0] = 1;
	for (int i = 1; i <= d; ++i) b1[i] = b1[i - 1] * b2[d], b1[i].set(N);
	printf("0 ");
	for (int i = 0; i < d; ++i) {
		for (int j = 1; j <= d; ++j) {
			int nw = i * d + j - 1;
			if (nw >= N) {
				puts("");
				cerr << ((double) clock() / CLOCKS_PER_SEC) << '\n';
				return 0;
			}
			int ans = 0;
			for (int k = 0; k <= nw; ++k) ans = (ans + (ll) b1[i][k] * b2[j][nw - k]) % MOD;
			ans = (ll) ans * INV[nw + 1] % MOD;
			printf("%d ", ans);
		}
	}
	puts("");
	cerr << ((double) clock() / CLOCKS_PER_SEC) << '\n';
	return 0;
}

右複合一些特定冪級數

  • 複合 \(cx\)，直接將 \(f_n = f_n c ^ n\) 即可。
  • 複合 \(x ^ a\)，將 \(f_i\) 的係數挪到 \(f_{ia}\) 即可
  • 複合 \(x+c\)，\(F(x + c) = \sum_{n\ge 0}f_n (x + c)^n = \sum_{n\ge 0}f_n\sum_{i = 0}^ n \binom{n}{i}x ^ i c ^ {n - i}= \sum_{i \ge 0} x ^ i \sum_{n \ge i} f_n c ^ {n - i}\) 。然後是一個差卷積的形式，可以通過翻轉 \(c\) 轉化為和卷積去做。
  • 複合 \(\sqrt {x + 1}\) ，不會。

CF438E

對於一棵有根無標號二叉樹，對其每個結點賦予一個正整數權值，且在 \(\{c_1,c_2,\dots,c_n\}\) 中。
給定 \(m\)，問權值為 \(0,1,\dots,m\) 的二叉樹數量，模 998244353998244353。

\(n\) 個節點無標號二叉樹生成函式是 \(C(x)\)，其中它滿足 $C(x) = x (C(x))^2+1 $

考慮對於一個節點定義 \(W(x) = \sum_{i = 1} ^ nx ^ {c_i}\)，那麼，我們就知道答案 \(F = C\circ W\)。那麼應該有

\[F(W(x)) = W(x) F(W(x))^ 2+1 \]

把他寫成

\[F =W F ^ 2 + 1 \] \[F= \frac{1±\sqrt {1 - 4W} }{2W} \]

考慮在化一下，變成

\[F = \frac{2}{1± \sqrt{1 - 4 W}} \]

考慮如果取減號，下面式子不能求逆。所以是加號。

然後就可以做了，使用多項式開根和多項式求逆。

\(\text{Lagrange}\) 反演

\(\text{Lagrange}\) 反演 基礎形式

若兩個形式冪級數滿足 \(f (g(x)) = x\) ，那麼稱 \(F\) 是 \(G\) 的複合逆，不難發現兩個冪級數互為複合逆

\[[x ^ n] g(x) =\frac{1}{n}[w ^ {-1}]\left(\frac{1}{f(w)}\right)^n \]

由於 \(\frac{1}{f(w)}\) 不能求逆，所以更常見的方式是

\[[x ^ n] g(x) =\frac{1}{n}[w ^ {n-1}]\left(\frac{w}{f(w)}\right)^n \]

\(\text{Lagrange}\) 反演 擴充套件形式

下面記 \(f(x)\) 的複合逆為 \(\hat{f}(x)\)，有

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

特殊情況，令 \(h(x) = x ^ k\)

那麼有

\[[x ^ n] g^k(x) = \frac{k}{n}[u ^ {n - k}]\left(\frac{u}{f(u)}\right) ^ n \]

若 \(f(g(x)) = h(x)\)。應該有 \(g(x) = h(\hat{f}(x))\)。

這時候居然有

\[[x ^ n]g(x) =\frac{1}{n} [x ^ {n - 1}]H'(x) \left (\frac{x}{f(x)}\right) ^ n \]

和之前的式子一樣。

\(\text{Lagrange}\) 反演 另類形式

設 \(g(x)\) 是 \(f(x)\) 的複合逆，則

\[[x ^ n]g ^ k(x) = [x ^{-k - 1}] \frac{F'(x)}{F ^ {n + 1}(x)} \]

改寫一下因該有

\[[x ^ n] g ^ k(x) = [x ^ {n - k}] F'(x)\left (\frac{x}{F(x)}\right) ^ {n + 1} \]

這可以求解不同的 \(k\) 的問題。

ABC222H

題目連結 H - Beautiful Binary Tree (atcoder.jp)

中文翻譯一下：給你一個正整數 \(n\)，滿足一下條件的為 \(N\) 階美麗二叉樹。

  • 每個頂點上寫著 \(0/1\)。
  • 每個葉子節點上寫著 \(1\)。
  • 可以進行以下操作 \(N - 1\) 次，使得根節點寫上 \(N\)，其他節點寫上 \(0\)。操作為：選擇兩個節點 \(u,v\)，其中 \(v\) 一定是 \(u\) 的兒子或者是孫子。將 \(a_u + a_v \to a_u,a_v \to 0\)。

首先根節點一定是 \(1\)，並且不應該有相鄰的 \(0\)。

怎麼辦？考慮一個樸素的 \(dp\) 吧。令 \(f_n\) 表示原問題的答案，然後 \(g_n\) 是強制令根為 \(0\)，但是每個兒子都是美麗二叉樹的方案。

那麼應該有

\[g_n = 2 f_n + \sum_{0 < j < n} f_j f_{n - j} \]

然後可以得到

\[\begin{aligned} f_1 &= 1\\ f_n &= 2(f_{n - 1} + g_{n - 1}) + \sum_{0 < j < i - 1} (f_j + g_j) (f _ {i - 1 - j} + g_{i - 1 - j}) \end{aligned} \]

可以使用 分治NTT來計算，複雜度為 \(O(N \log^2 N)\)。

但是題目要求 \(n\le 10^7\)，時間絕對超時。考慮生成函式

\[g(x) = 2f(x) + f ^ 2(x) \]

然後

\[f(x) = x + 2xf(x) + 2x g(x) + x (f(x) + g(x)) ^ 2 \]

即

\[f(x) = x(1 + 3f(x) + f^2(x))^2 \]

考慮拉格朗日反演，\(f(x) = x(1 + 3f(x) +f ^ 2(x))\)，他的複合逆就是 \(g(x) = \frac{x}{(1 + 3x + x ^ 2)^2}\)。

所以說 \([x ^ N]F(x) = \frac{1}{N}[x ^ {N - 1}](1 + 3x + x^2)^{2N}\) 這玩意怎麼算？

二項式定理展開答案就是 \(\frac{\sum_{j = 0} ^ {2N}\binom{2n}{j} \binom{2n - j}{n - 1 - 2j}3 ^ {n - 1 - 2j}}{n}\) 。

fuss-catalan 數

求有 \(N\) 個內點的 \(K\) 叉樹數量，滿足兒子內有順序，滿足非葉子節點都有 \(K\) 個兒子。

一個等價的問題，從 \((0,0)\) 走到 \((kn,0)\)，然後每次可以走 \((+1,+1)\) 或者是 \((+1,-(k - 1))\)，並且 \(y\) 座標不能在 \(0\) 以下的路徑方案數量。

用生成函式來表達答案 \(f(x) = x (f(x) + 1) ^ k\)，求出來 \(f_n\) 即可。。

得到 \(x = \frac{f(x)}{(f(x) + 1)^ k}\)。所以 $f(x) $ 的複合逆是 \(\frac{x}{(x + 1) ^ k}\)。根據拉格朗日反演。

\[f_n = [x ^ n]f(x) = \frac{1}{n} [x ^ {n- 1}]\left(\frac{x}{g(x)}\right)=\frac{1}{n}[x ^ {n - 1}](x + 1)^{kn}=\frac{\binom{kn}{n - 1}}{n} \]

所以解決了。

相關推薦

多項式複合和拉格朗日反演 學習筆記

20220520 一些小東西 卡特蘭數和去掉根的二叉樹對應。考慮卡特蘭數的遞推公式 \\(f(n) = \\sum_{i = 0}^{n - 1} f(i) f(n - 1 - i)\\)，並有 \\(f(0) = f(1) = 1\\)。

拉格朗日反演學習筆記

拉格朗日反演 用於在\\(O(n\\log n)\\)的時間內求\\([x^n]G(x)\\)，其中\\(G(x)\\)滿足\\(F(G(x))=x\\)，\\(F(x)\\)已知，且\\([x^0]F(x)=[x^0]G(x)=0,[x^1]F(x) \\neq 0,[x^1]G(x) \\neq 0\\)。

BZOJ3684 大朋友與多叉樹(拉格朗日反演)

BZOJ3684 大朋友與多叉樹(拉格朗日反演) 題目大意 我們的大朋友很喜歡電腦科學，而且尤其喜歡多叉樹。對於一棵帶有正整數點權的有根多叉樹，如果它滿足這樣的性質，我們的大朋友就會將其稱作神犇的：點權為1的結點是

拉格朗日插值學習筆記

拉格朗日插值學習筆記 簡介 在數值分析中，拉格朗日插值法是以法國十八世紀數學家約瑟夫·拉格朗日命名的一種多項式插值方法。許多實際問題中都用函式來表示某種內在聯絡或規律，而不少函式都只能通過實驗和觀測來了

qbxt DAY 6 T3 柯西不等式和拉格朗日不等式

qbxt DAY 6 T3 其實主要是對拉格朗日不等式的瞭解和認識。 \\((\\sum a_i^2)(\\sum b_i^2)=(\\sum a_ib_i)^2+\\sum\\limits_{1\\leq i< j\\leq n}(a_ib_j-a_jb_i)^2\\)

拉格朗日插值 學習

技術標籤：數學總結 前置知識 1 + 1 == 2 內容 拉格朗日插值是用來解決對於給定的n+1個點值，確定一個函式f()，使其能夠經過這n+1個點，且f()是n次多項式，同時，對於一個給定的點x，它可以快速求出f(x)的值，即

自然冪數和（拉格朗日插值）

求 $\\sum_{i=1}^{n} i^k$. 這是一個關於 $n$ 的 $k+1$ 次多項式. 所以可以取 $k+2$ 個點帶進去，然後用拉格朗日插值法來求值.

多項式入門——拉格朗日插值

多項式入門——拉格朗日插值 插值用來求解這樣一類問題：給定 \\(n\\) 點 \\((x_i,y_i)\\) 求過這些點的多項式。

【知識點】拉格朗日乘數法

簡介： 在滿足約束條件$\\varphi(x_1 ,x_2 ,\\cdots ,x_n )=0$時求$f(x_1 ,x_2 ,\\cdots ,x_n )$的極值。

LuoguP4463 [集訓隊互測2012] calc DP+拉格朗日插值

這道題的題意不太明確. 應該是兩個序列 $a,b$ 不同，當且僅當存在位置 $i$ 使得 $a[i]$ 不等於 $b[i]$.

樣條之拉格朗日Lagrange(一元全區間)插值函式[轉]

樣條之拉格朗日Lagrange(一元全區間)插值函式 這是使用拉格朗日插值函式生成的樣條曲線。在數值分析中，拉格朗日插值法是以法國十八世紀數學家約瑟夫·拉格朗日命名的一種多項式插值方法。許多實際問題中都用

拉格朗日插值法學習筆記

適用範圍：給出n個點$(x_i,y_i)$，過這n個點能夠確定一個最高n-1次的多項式$f(x)$，求$f(k)$

【模板】 拉格朗日插值

由n + 1 個點可以確定一個n次多項式。 拉格朗日插值法可以在 n^2 複雜度確定一個多項式並進行求解。

CodeForces - 622F The Sum of the k-th Powers 拉格朗日插值模板題

定理： {an} 是一個p階等差數列的充要條件是數列的通項 an 為n的一個p次多項式。

cf 995F Cowmpany Cowmpensation - 拉格朗日插值 + 樹形dp

傳送門 第一次獨立完成div1的F題，

LOJ #165. 拉格朗日插值

重心拉格朗日插值模板題。 用上面那篇文章的公式計算。插入時更新每個 \\(w_i\\)，詢問時分別 \\(O(n)\\) 計算 \\(l(k)\\) 和 \\(\\sum\\limits_{i=1}^n \\frac{w_i}{k-x_i}\\) 即可。注意使用上述公式時 \\(k-x_i\\

[OI筆記]利用拉格朗日乘數法求函式的最值

\\(about\\) 為什麼寫這篇\\(Blog\\)呢\\(...\\) 拉格朗日乘數法在今天訓練的一道題上用到了\\(,\\)當場\\(wyj/pcf/csl\\)都正確的推出了式子\\(.\\)

拉格朗日插值法

\\(n\\)階的多項式 \\(f(x)\\) 可以由 \\(n+1\\) 個點確認。 若現有 \\(n+1\\) 個點\\((x_i,y_i) \\ \\ , \\ i \\in [0,n]\\) 在 \\(f(x)\\) 上

【洛谷7116】微信步數（拉格朗日插值）

點此看題面 有一個\\(k\\)維場地，第\\(i\\)維範圍為\\([1,w_i]\\)。 有\\(n\\)個操作，每個操作可以表示為將某一維的座標\\(±1\\)。

【學習筆記】拉格朗日插值

學習多項式的第一步。 參考資料： attack的luogu部落格 oi wiki拉格朗日插值 Apocryphal的luogu部落格