数学 (OI)

数学 (OI)

前言

说实话这首歌有点反。

Jar of Love

Another sunrise, another sunset
Soon it'll all be yesterday
Another good day, another bad day,
What did you do today?
Why do we choose to chase what we'll lose?
What you want isn't what you have.
What you have may not be yours, to keep.
If I could find love, at a stop, in a park with open arms,
I would save all my love, in a jar,
made of sparks, sealed in my beating heart,
Could it be yours to keep, the Jar of Love.
Another left turn, another head turns
Could he be someone I deserve?
Another right turn, another lesson learned
Never leave an open flame to burn
Why do we choose to chase what we'll lose?
What you want isn't what you have.
What you have may not be yours, to keep.
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep, the Jar of Love.
Could you be my love
Could you be my love
Could you be my love
Could you be my love
Could you be her love
Could you be his love
Could you be my love
Could I be you love
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
If I could find love,
at a stop, in a park with open arms,
I would save all my love,
in a jar,made of sparks,
sealed in my beating heart,
Could it be yours to keep
the Jar of Love.
Could it be yours to keep
the Jar of Love.
Could it be yours to keep
the Jar of Love.

可能是最后一篇学习笔记(自古不退役不学多项式)。

FFT 学习笔记

傅里叶

很久之前就学了，现在好像忘了，连着 ntt,fwt 重新学一下。顺便打个学习笔记。

把 \(n\) 当做 \(2^k\)

用复数做多项式的根，具体的，\(x = w^1_n\) 一次单位根, 即 \((\cos \frac{2 \pi}{n}, \sin \frac{2 \pi}{n})\)。

理由的话可以证一下:

设 \(A(x) = \sum\limits_{i = 0}^{n - 1}a_ix^i\) 多项式，\(x\) 做 \(0 \dots n - 1\) 次单位根带出来的离散傅里叶变换为 \((y_0, \dots, y_{n - 1})\), 接着设 \(B(x) = \sum\limits_{i = 0}^{n - 1}y_ix^i\), 将 \(x^{-1}\) 带入作单位根得出新的离散傅里叶变换 \((z_0, \dots, z_{n - 1})\)。

那按题意模拟，

\[\begin{aligned} z_k &= \sum_{i = 0}^{n - 1}y_i \left(\omega_{n}^{k}\right)^{-i} \\ &= \sum_{i = 0}^{n - 1}\left(\sum_{j = 0}^{n - 1}a_j \left(\omega_{n}^{j}\right)^j\right) \left(\omega_{n}^{k}\right)^{-i} \\ &= \sum_{j = 0}^{n - 1} a_j\left(\sum_{i = 0}^{n - 1}(\omega_{n}^{j - k})^i\right) \\ &= na_k \end{aligned}\]

所以这个玩意儿就是做两遍傅里叶变换的事儿。

快速

然后就是优化复杂度的问题了。傅里叶叔叔作为一个连电脑都没见过的人，他是闲得蛋疼才会优化复杂度，于是他蛋完好。

后人用分治优化了 fft 复杂度。

具体原理如下：

设 \(A(x) = a_0x^0 + \dotsb + a_{n - 1}x^{n - 1}\)， 按奇偶分个类：

\[A(x) = (a_0 + a_2x^2 + \dotsb + a_{n - 2}x^{n - 2}) + (a_1x + \dotsb + a_{n - 1}x^{n - 1}) \]

整俩多项式：

\[A_1(x) = a_0 + a_2x + \dotsb + a_{n - 2}x^{\frac{n}{2} - 1} \]

\[A_2(x) = a_1 + a_3x + \dotsb + a_{n - 1}x^{\frac{n}{2} - 1} \]

显显显 \(A(x) = A_1(x^2) + xA_2(x^2)\)。

假设 \(k < \frac{n}{2}\)，将 \(x = \omega_{n}^{k}\) 带入，

\[\begin{aligned} A(\omega_{n}^k) &= A_1(\omega_n^{2k}) + \omega_n^kA_2(\omega_n^{2k}) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) + \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{aligned} \]

对于 \(A(\omega_{n}^{k + \frac{n}{2}})\)：

\[\begin{aligned} A(\omega_{n}^{k + \frac{n}{2}}) &= A_1(\omega_n^{2k + n}) + \omega_n^{k + \frac{n}{2}}A_2(\omega_n^{2k + n}) \\ &= A_1(\omega_{\frac{n}{2}}^{k}) - \omega_n^kA_2(\omega_{\frac{n}{2}}^{k}) \end{aligned}\]

所以只要对于 \(A_1A_2\) 在 \((\omega_{\frac{n}{2}}^{0} \dotsb \omega_{\frac{n}{2}}^{\frac{n}{2} - 1})\) 就可以求出 \(A\) 在 \((\omega_{n}^{0} \dotsb \omega_{n}^{n - 1})\) 离散傅里叶变换。

\(n - 1\) 直接 return。

实现

其实还有递归版，就是常数忒他妈大，写了也没啥意义，就不写了，这里是倍增法，比较快。

然后就是各位 大牛 们纷纷拽高级名词，什么蝴蝶变换，根本看不懂。

然后看代码就是用个变量数组存起来，减少常数。

\(QwQ\)， 有点唐。

考虑一下非递归实现就是要从下往上走，那么我们就需要最终位置。

简单手摸一下得出， \(x\) 的最后位置是其二进制翻转的值。

一点一点向上还原即可。

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 21) + 100; 
class _complex {
private:
double x, y; 
public:
_complex(double X = 0, double Y = 0): x(X), y(Y) {}
inline double Real() {return x; }
inline double Imag() {return y; }
inline _complex operator + (_complex a) {return _complex(x + a.x, y + a.y); }
inline _complex operator - (_complex a) {return _complex(x - a.x, y - a.y); }
inline _complex operator * (_complex a) {return _complex(x * a.x - y * a.y, a.x * y + x * a.y); }
} a[N], b[N]; 
int n, m, num = 1; 
const double Pi = acos(-1.0); 
int l, r[N]; 

void FFT(_complex *A, int type) {
	for (int i = 0; i < num; ++ i) {
		if (i < r[i]) swap(A[i], A[r[i]]); 
	}
	for (int mid = 1; mid < num; mid <<= 1) {
		_complex Wn = _complex(cos(Pi / mid), sin(Pi / mid) * type); 
		for (int R = mid << 1, j = 0; j < num; j += R) {
			_complex w = _complex(1, 0); 
			for (int k = 0; k < mid; ++ k, w = w * Wn) {
				_complex x = A[j + k], y = w * A[j + k + mid]; 
				A[j + k] = x + y, A[j + mid + k] = x - y; 
			}
		}
	}
}

signed main() {
	freopen("1.in", "r", stdin); 
	freopen("1.out", "w", stdout); 
	ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); 
	cin >> n >> m; 
	for (int i = 0, x; i <= n; ++ i) 
		cin >> x, a[i] = _complex(x, 0); 
	for (int i = 0, x; i <= m; ++ i) 
		cin >> x, b[i] = _complex(x, 0); 
	while (num <= n + m) num <<= 1, ++ l; 
	for (int i = 0; i < num; ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1)); 
	FFT(a, 1), FFT(b, 1); 
	for (int i = 0; i < num; ++ i) a[i] = a[i] * b[i]; 
	FFT(a, -1); 
	for (int i = 0; i <= n + m; ++ i) cout << (int)(a[i].Real() / num + 0.5) << ' '; 
}

NTT 学习笔记

长相和 fft 基本一致，利用原根的性质，获得了顶替 \(\omega\) 的资格

直接贴码：

CODE

const int G = 3; 
const int Gi = 332748118; // 3的逆元
inline void NTT(ll *A, int type, int length) {
	for (int i = 0; i < length; ++ i) 
		if (i < pos[i]) swap(A[i], A[pos[i]]); 
	for (int mid = 1; mid < length; mid <<= 1) {
		ll Wn = Quick_Pow(type == 1 ? G : Gi, (mod - 1) / (mid << 1)); 
		for (int R = mid << 1, j = 0; j < length; j += R) {
			ll w = 1; 
			for (int k = 0; k < mid; ++ k, w = Wn * w % mod) {
				ll x = A[j + k], y = w * A[j + mid + k] % mod; 
				A[j + k] = (x + y) % mod; A[j + mid + k] = (x + mod - y) % mod; 
			}
		}
	}
}

多项式求逆

字面意思，求 \(G(x)\)，使得

\[F(x)G(x) \equiv 1 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

设：

\[G'(x)F(x) \equiv 1 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^\frac{n}{2}\right) \]

把常数项移过去再平方得：

\[G'(x)^2F(x)^2 - 2G'(x)F(x) \equiv -1 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

同乘 \(-1\)：

\[2G'(x)F(x) - G'(x)^2F(x)^2 \equiv F(x)G(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

同消去 \(F(x)\)：

\[G(x) \equiv 2G'(x) - G'(x)^2F(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

然后做完了，复杂度估计为 \(O(n\log n)\)

贴码：

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 23) + 100; 
const int mod = 998244353; 
const int G = 3; 
const int Gi = 332748118; 

int n, r[N]; ll f[N], g[N], a[N], b[N]; 

inline ll Quick_Pow(ll a, ll b) {
	ll ans = 1; 
	while (b) {
		if (b & 1) ans = (ans * a) % mod; 
		b >>= 1, a = (a * a) % mod; 
	}
	return ans; 
}
inline void NTT(ll *A, int type, int length) {
	for (int i = 0; i < length; ++ i) 
		if (i < r[i]) swap(A[i], A[r[i]]); 
	for (int mid = 1; mid < length; mid <<= 1) {
		ll Wn = Quick_Pow(type == 1 ? G : Gi, (mod - 1) / (mid << 1)); 
		for (int R = mid << 1, j = 0; j < length; j += R) {
			ll w = 1; 
			for (int k = 0; k < mid; ++ k, w = Wn * w % mod) {
				ll x = A[j + k], y = w * A[j + mid + k] % mod; 
				A[j + k] = (x + y) % mod; A[j + mid + k] = (x + mod - y) % mod; 
			}
		}
	}
}

void solve(int lg) {
	if (lg == 0) return g[0] = Quick_Pow(f[0], mod - 2), void(); 
	solve(lg - 1); 
	static int length; length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(a, 0, sizeof(ll) * (length << 1)); memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	for (int i = 0; i < length; ++ i) a[i] = f[i], b[i] = g[i]; 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = ((2ll + mod - a[i] * b[i] % mod) % mod) * b[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) g[i] = a[i] * nueyong % mod; 
}

signed main() {
	freopen("1.in", "r", stdin); 
	freopen("1.out", "w", stdout); 
	ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); 
	cin >> n; int tmp = 1, lg = 0; 
	for (int i = 0; i < n; ++ i) cin >> f[i]; 
	while (tmp < n) tmp <<= 1, ++ lg; 
	solve(lg); 
	for (int i = 0; i < n; ++ i) cout << (g[i] + mod) % mod << ' '; 
}

多项式开根

不太像字面意思，求：

\[G(x)^2 \equiv F(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

感觉和多项式求逆思想差不多，都是倍增。

设：

\[H(x)^2 \equiv F(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^\frac{n}{2}\right) \]

所以可得：

\[H(x) \equiv G(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^\frac{n}{2}\right) \]

\[H(x) - G(x) \equiv 0 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^\frac{n}{2}\right) \]

\[H(x)^2 - 2H(x)G(x) + G(x)^2 \equiv 0 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

\[F(x) + H(x)^2 \equiv 2H(x)G(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

然后就可得

\[G(x) \equiv \frac{F(x) + H^2(x)}{2H(x)} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

具体解决办法好像应该 NTT 和多项式求逆一块用。

码：

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 19) + 100; 
const int mod = 998244353; 
const int G = 3; 
const int Gi = 332748118; 

int n, r[N]; 

inline ll Quick_Pow(ll a, ll b) {
	ll ans = 1; 
	while (b) {
		if (b & 1) ans = (ans * a) % mod; 
		b >>= 1, a = (a * a) % mod; 
	}
	return ans; 
}
inline void NTT(ll *A, int type, int length) {
	for (int i = 0; i < length; ++ i) 
		if (i < r[i]) swap(A[i], A[r[i]]); 
	for (int mid = 1; mid < length; mid <<= 1) {
		ll Wn = Quick_Pow(type == 1 ? G : Gi, (mod - 1) / (mid << 1)); 
		for (int R = mid << 1, j = 0; j < length; j += R) {
			ll w = 1; 
			for (int k = 0; k < mid; ++ k, w = w * Wn % mod) {
				ll x = A[j + k], y = w * A[j + k + mid] % mod; 
				A[j + k] = (x + y) % mod, A[j + k + mid] = (x - y + mod) % mod; 
			}
		}
	}
}

ll a[N], b[N], inver[N], beinver[N]; ll f[N], g[N]; 

void Inv(int lg) {
	if (lg == 0) return inver[0] = Quick_Pow(beinver[0], mod - 2), void(); 
	Inv(lg - 1); 
	int length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(a, 0, sizeof(ll) * (length << 1)), memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	for (int i = 0; i < length; ++ i) a[i] = beinver[i], b[i] = inver[i]; 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = (2ll - a[i] * b[i] % mod) % mod * b[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) inver[i] = a[i] * nueyong % mod; 
}

void Polysqrt(int lg) {
	if (lg == 0) return g[0] = 1, void(); 
	Polysqrt(lg - 1); 
	int length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(inver, 0, sizeof(ll) * (length << 1)), memset(beinver, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < length; ++ i) beinver[i] = 2ll * g[i] % mod; 
	Inv(lg); 
	for (int i = 0; i < length; ++ i) 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	memset(a, 0, sizeof(ll) * (length << 1)), memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < length; ++ i) a[i] = f[i], b[i] = g[i]; 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1), NTT(inver, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = (a[i] + b[i] * b[i] % mod) % mod * inver[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) g[i] = a[i] * nueyong % mod; 
}

signed main() {
	freopen("1.in", "r", stdin); 
	freopen("1.out", "w", stdout); 
	ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); 
	cin >> n; 
	for (int i = 0; i < n; ++ i) cin >> f[i]; int tmp = 1, lg = 0; 
	while (tmp < n) tmp <<= 1, ++ lg; 
	Polysqrt(lg); 
	for (int i = 0; i < n; ++ i) cout << (g[i] + mod) % mod << ' '; 
}

多项式 \(\ln\)

字面意思，求：

\[G(x) \equiv \ln F(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

设 \(f(x) = \ln x\)，则：

\[G(x) \equiv f(F(x)) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

两边同时求导得：

\[G'(x) \equiv f'(F(x))F'(x) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

因为 \(\ln(x)' = \frac{1}{x}\)，所以：

\[G'(x) \equiv \frac{F'(x)}{F(x)} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

对 \(F\) 求导得到 \(F'\)，多项式求逆后再乘积，得到 \(G'\)，积个分完事了。

简单一说就是：

\[\int x^a = \frac{x^{a + 1}}{a + 1} \]

贴码：

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 19) + 100; 
const int mod = 998244353; 
const int G = 3; 
const int Gi = 332748118; 

int n, r[N]; 
ll a[N], b[N], inver[N], beinver[N]; ll f[N], g[N], c[N]; 

inline ll Quick_Pow(ll a, ll b) {
	ll ans = 1; 
	while (b) {
		if (b & 1) ans = (ans * a) % mod; 
		b >>= 1, a = (a * a) % mod; 
	}
	return ans; 
}
inline void PolyInt(ll* A) {
	for (int i = 0; i < n; ++ i) c[i + 1] = Quick_Pow(i + 1, mod - 2) * A[i] % mod; c[0] = 0; 
}
inline void Polydiff(ll *A) {
	for (int i = 1; i < n; ++ i) c[i - 1] = A[i] * i % mod; c[n - 1] = 0; 
}
inline void NTT(ll *A, int type, int length) {
	for (int i = 0; i < length; ++ i) 
		if (i < r[i]) swap(A[i], A[r[i]]); 
	for (int mid = 1; mid < length; mid <<= 1) {
		ll Wn = Quick_Pow(type == 1? G: Gi, (mod - 1) / (mid << 1)); 
		for (int R = mid << 1, j = 0; j < length; j += R) {
			ll w = 1; 
			for (int k = 0; k < mid; ++ k, w = w * Wn % mod) {
				ll x = A[j + k], y = w * A[j + k + mid] % mod; 
				A[j + k] = (x + y) % mod, A[j + k + mid] = (x + mod - y) % mod; 
			}
		}
	}
}

void Inv(int lg) {
	if (lg == 0) return inver[0] = Quick_Pow(beinver[0], mod - 2), void(); 
	Inv(lg - 1); 
	int length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(a, 0, sizeof(ll) * (length << 1)), memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	for (int i = 0; i < length; ++ i) a[i] = beinver[i], b[i] = inver[i]; 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = (2ll - a[i] * b[i] % mod + mod) % mod * b[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) inver[i] = a[i] * nueyong % mod; 
}

inline void PolyLn(int lg) {
	int length = 1 << lg, nueyong = Quick_Pow(length << 1, mod - 2); 
	memcpy(beinver, f, sizeof(f)); Inv(lg); Polydiff(f); memcpy(a, c, sizeof(a)); 
	memset(b, 0, sizeof(b)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	NTT(a, 1, length << 1), NTT(inver, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = a[i] * inver[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) b[i] = a[i] * nueyong % mod; 
	PolyInt(b); memcpy(g, c, sizeof(g)); 
}

signed main() {
	freopen("1.in", "r", stdin); 
	freopen("1.out", "w", stdout); 
	ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); 
	cin >> n; int tmp = 1, lg = 0; 
	for (int i = 0; i < n; ++ i) cin >> f[i]; 
	while (tmp < n) tmp <<= 1, ++ lg; 
	PolyLn(lg); 
	for (int i = 0; i < n; ++ i) cout << g[i] << ' '; 
}

多项式EXP

字面意思，求：(g)

\[G(x) \equiv e^{F(x)} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

牛顿迭代法

变换一下得到：

\[\ln G(x) - F(x) \equiv 0 \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

然后直接把上面的左边当一个多项式，设为 \(P\) ，套牛顿迭代，把 \(F(x)\) 当常量使唤。得到：

\[G(x) \equiv G_0(x) - \frac{P(G_0(x))}{P'(G_0(x))} \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

直接得到：

\[G(x) \equiv G_0(x) \left(1 - \ln G_0(x) + F(x)\right) \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \left(\text{mod} \ x^n\right) \]

复杂度 \(O(n \log n)\)。贴码：

CODE

#include <bits/stdc++.h>
typedef long long ll; 
using namespace std; 
const int N = (1 << 19) + 100; 
const int mod = 998244353; 
const int G = 3; 
const int Gi = 332748118; 

int n, r[N]; 
ll a[N], b[N], inver[N], beinver[N]; ll f[N], g[N], c[N], d[N]; 

inline ll Quick_Pow(ll a, ll b) {
	ll ans = 1; 
	while (b) {
		if (b & 1) ans = (ans * a) % mod; 
		b >>= 1, a = (a * a) % mod; 
	}
	return ans; 
}
inline void PolyInt(ll* A, int length) {
	memset(c, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < length; ++ i) c[i + 1] = Quick_Pow(i + 1, mod - 2) * A[i] % mod; c[0] = 0; 
}
inline void Polydiff(ll *A, int length) {
	memset(c, 0, sizeof(ll) * (length << 1)); 
	for (int i = 1; i < length; ++ i) c[i - 1] = A[i] * i % mod; c[length - 1] = 0; 
}
inline void NTT(ll *A, int type, int length) {
	for (int i = 0; i < length; ++ i) 
		if (i < r[i]) swap(A[i], A[r[i]]); 
	for (int mid = 1; mid < length; mid <<= 1) {
		ll Wn = Quick_Pow(type == 1? G: Gi, (mod - 1) / (mid << 1)); 
		for (int R = mid << 1, j = 0; j < length; j += R) {
			ll w = 1; 
			for (int k = 0; k < mid; ++ k, w = w * Wn % mod) {
				ll x = A[j + k], y = w * A[j + k + mid] % mod; 
				A[j + k] = (x + y) % mod, A[j + k + mid] = (x + mod - y) % mod; 
			}
		}
	}
}

void Inv(int lg) {
	if (lg == 0) return inver[0] = Quick_Pow(beinver[0], mod - 2), void(); 
	Inv(lg - 1); 
	int length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(a, 0, sizeof(ll) * (length << 1)), memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	for (int i = 0; i < length; ++ i) a[i] = beinver[i], b[i] = inver[i]; 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = (2ll - a[i] * b[i] % mod + mod) % mod * b[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) inver[i] = a[i] * nueyong % mod; 
}

inline void PolyLn(int lg, ll *A, ll *B) {
	int length = 1 << lg; ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memcpy(beinver, A, sizeof(ll) * (length << 1)); memset(inver, 0, sizeof(ll) * (length << 1)); 
	Inv(lg); Polydiff(A, length); 
	memcpy(a, c, sizeof(ll) * (length << 1)); memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	NTT(a, 1, length << 1), NTT(inver, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = a[i] * inver[i] % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) b[i] = a[i] * nueyong % mod; 
	PolyInt(b, length); memcpy(B, c, sizeof(ll) * (length << 1)); 
}

void PolyExp(int lg) {
	if (lg == 0) return g[0] = 1, void(); 
	PolyExp(lg - 1); 
	int length = (1 << lg); ll nueyong = Quick_Pow(length << 1, mod - 2); 
	memset(d, 0, sizeof(ll) * (length << 1)); PolyLn(lg, g, d); 
	memset(a, 0, sizeof(ll) * (length << 1)), memset(b, 0, sizeof(ll) * (length << 1)); 
	for (int i = 0; i < length; ++ i) a[i] = g[i], b[i] = f[i]; 
	for (int i = 0; i < (length << 1); ++ i) r[i] = (r[i >> 1] >> 1) | ((i & 1) << lg); 
	NTT(a, 1, length << 1), NTT(b, 1, length << 1), NTT(d, 1, length << 1); 
	for (int i = 0; i < (length << 1); ++ i) a[i] = a[i] * ((1ll + mod - d[i] + b[i]) % mod) % mod; 
	NTT(a, -1, length << 1); 
	for (int i = 0; i < length; ++ i) g[i] = a[i] * nueyong % mod; 
}

signed main() {
	freopen("1.in", "r", stdin); 
	freopen("1.out", "w", stdout); 
	ios::sync_with_stdio(0), cin.tie(0), cout.tie(0); 
	cin >> n; int tmp = 1, lg = 0; 
	for (int i = 0; i < n; ++ i) cin >> f[i]; 
	while (tmp < n) tmp <<= 1, ++ lg; 
	PolyExp(lg); 
	for (int i = 0; i < n; ++ i) cout << g[i] << ' '; 
}

本来是有例题的，刚打完电脑死机了，没了。

生成函数

好吧，有点唐，具体是啥，就是个数列。

OGF

具体讲一下怎么转封闭形式，这个还是非常有用的。

1. \[f(x) = \sum_{i = 0}x^{2i} \]

显显显答案是 $$\cfrac{1}{1 - x^2}$$

为啥就不证了。

2. \[f(x) = \sum_{i = 0}(i + 1)x^i \]

求导题。

\[\begin{aligned} f(x) &= \sum_{i = 0}(i + 1)x^i \\ &= \sum_{i = 0}ix^{i - 1} \\ &= \left(\sum_{i = 0}x^i\right)^\prime \\ &= \left(\cfrac{1}{1 - x}\right)^\prime \\ &= \cfrac{1}{(1 - x)^2} \end{aligned}\]

3. \[f(x) = \sum_{i = 0}C_{m + i}^{i}x^i \]

这个太唐了，只能当结论记住。

\[f(x) = \frac{1}{(1 - x)^{m + 1}} \]

考虑使用归纳。 \(n = 0\) 时显然。

\[\begin{aligned} \cfrac{1}{(1 - x)^{m + 1}} &= \frac{1}{(1 - x)} \times \cfrac{1}{(1 - x)^m} \\ &= \sum_{i = 0} x^i \times \sum_{i = 0}C_{m + i - 1}^{i}x^i \\ &= \sum_{i = 0} x^i\sum_{j = 0}^i C_{m - 1 + j}^{j} \\ &= \sum_{i = 0} C_{m + j}^j x^i \end{aligned}\]