多项式学习笔记（二）：多项式初等函数

牛顿迭代

高中数学中，我们学过牛顿法求函数零点，其实，对于多项式的计算，也有牛顿迭代法。

回忆一下牛顿法的过程，对于二次函数 \(f(x) = x^2 - 10\)，我们首先猜测一个零点 \(x_0\)，求出 \((x_0, f(x_0))\) 处的切线。

比如我们取 \(x_0 = 3\)，由于 \(f'(x) = 2x\)，我们可以很轻松（如果你会求导）的写出 \((3, -1)\) 处切线方程 \(y = 6x - 19\)，得到一个新的近似零点 \(x_1 = \displaystyle\frac{19}{6}\)。

再带一次：\((\displaystyle\frac{19}{6}, \frac{1}{36})\) 处的切线方程 \(y = \displaystyle\frac{19}{3}x - \frac{721}{36}\)，得到新的近似零点 \(x_2 = \displaystyle\frac{721}{228}\approx 3.16228\)，\(3.16228^2 = 10.0000147984\)，已经相当接近 \(10\) 了，可见其逼近速度之快。

现在考虑将 \(x\) 变为一个多项式 \(F(x)\)，那么此时就可以快速求出这个多项式的系数。

假设 \(G(F(x)) = 0\)，如何在知道 \(G\) 的情况下求出 \(F(x)\)（听说 \(G\) 叫做泛函数，蒟蒻瑟瑟发抖）。

考虑类似数学归纳法的方法。

首先，常数项是好求的，相当于是求 \(G(x) = 0\) 时 \(x\) 的取值，记为 \(F_0(x)\)，也就相当于 \(F(x) \equiv F_0(x) \pmod x^1\)。

假设我们已经知道了 \(F(x)\) 在 \(\pmod{x^n}\) 时的函数 \(F_n(x)\)，将 \(G(x)\) 在 \(F_n(x)\) 处进行泰勒展开，由于 \(G(F_n(x)) = 0\)，因此 \(\displaystyle\sum_{i \geq 0} \frac{G^{(i)}F_n(x)}{i!}[F(x) - F_n(x)]^i = G(F(x)) \equiv 0 \pmod{x^n}\)。

现在我们将模数扩大到 \(x^{2n}\)，右边式子中的 \(0\) 显然不变，而左边式子中 \(F(x) - F_n(x)\) 的最低项至少是 \(x^n\)，一旦进行平方或更高次方，最低项至少是 \(x^{2n}\)，会被取模掉，最终该式子仅剩下 \(G(F_n(x)) + G'(F_n(x))(F(x) - F_n(x))\)，整理可得 \(F(x) \equiv F_n(x) - \displaystyle\frac{G[F_n(x)]}{G'[F_n(x)]} \pmod{x^{2n}}\)，这时，我们就找到了一个将 \(F_n(x)\) 扩大到 \(f_{2n}(x)\) 的方法，通过不断的迭代，最终可以求出前任意项的系数。

举个例子，现在已知 \(G(F(x)) = F^2(x) - F(x) - 2x = 0\)。

我们先求常数项，由于常数项可以用 \(F(0)\) 表示，那么原式就等于 \(F^2(0) - F(0) = 0\)，解出 \(F(0) = 1\) 或 \(0\)，\(0\) 没有太多含义，拿 \(1\) 举例，记为 \(F_0(x)\)，套用牛顿迭代的式子，\(F_1(x) = 1 - \displaystyle\frac{1 - 1 - 2x}{2 - 1} = 1 + 2x\)，代入原式可得：\(1 + 4x + 4x^2 - 1 - 2x - 2x \equiv 0 \pmod{x^2}\)，这是满足条件的，依次推下去就可以得到想要的多项式。

快速多项式计算

加法

\(H(x) = F(x) + G(x)\)，求 \(H(x)\)。

解：都已经 \(O(n)\) 了，还怎么优化！！！

乘法

\(H(x) \equiv F(x)G(x) \pmod{x^n}\)，求 \(H(x)\)。

解：参考多项式学习笔记（一）（2024.7.6）

完整代码：洛谷 P3803 【模板】多项式乘法（FFT）

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N];
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int f[N << 1], g[N << 1], n, m;
signed main(){
	scanf("%lld%lld", &n, &m);
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	for(int i = 0; i <= m; i++)
		scanf("%lld", &g[i]);
	int lim = 1;
	while(lim <= n + m)
		lim <<= 1;
	NTT(f, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i <= lim; i++)
		f[i] = f[i] * g[i] % MOD;
	NTT(f, lim, -1); 
	for(int i = 0; i <= n + m; i++)
		printf("%lld ", f[i]);
	return 0;
}

求逆

\(F(x)G(x) \equiv 1 \pmod{x^n}\)。

求 \(G(x)\)。

解：从这里开始就要用到牛顿迭代了。

先考虑初始情况，常数项 \(G_0(x) = \displaystyle\frac{1}{F(0)}\) 是好求的，套入牛顿迭代的式子可以得到 \(G_{k + 1}(x) \equiv G_k(x) - \displaystyle\frac{\frac{1}{G_k(x)} - F(x)}{-\frac{1}{G_k^2(x)}} \pmod{x^{2n}}\)。

整理可得：\(G_{k + 1}(x) \equiv 2G_k(x) + G_k^2(x)F(x) \pmod{x^{2n}}\)，现在就可以通过不断多项式乘法来求出 \(G(x)\) 了。

完整代码：洛谷 P4238 【模板】多项式乘法逆

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N];
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void solve(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	solve(f, g, n >> 1);
	int lim = 1;
	while(lim <= n * 2 + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f[N << 1], g[N << 1], n;
signed main(){
	scanf("%lld", &n);
	n--;
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	solve(f, g, n);
	for(int i = 0; i <= n; i++)
		printf("%lld ", g[i]);
	return 0;
}

除法/取模

\(F(x) \equiv Q(x)G(x) + R(x) \pmod{x^n}\)。

求 \(Q(x),R(x)\)。

解：设 \(n\) 是 \(F(x)\) 的次数，\(m\) 是 \(G(x)\) 的次数。

设 \(rF(x)\) 表示 \(x^nF(\displaystyle\frac 1x)\)，意义是将 \(F(x)\) 的系数翻转，\(rG(x), rQ(x)\) 同理。

由于 \(F(x) \equiv Q(x)G(x) + R(x) \pmod{x^n}\)，因此 \(F(\displaystyle\frac 1x) \equiv Q(\frac 1x)G(\frac 1x) + R(\frac 1x) \pmod{x^n}\)。

两边同乘 \(x^n\)，由于 \(Q(x)\) 是一个 \(n - m\) 项多项式，\(R(x)\) 是一个 \(m - 1\) 次多项式，那么将 \(x^n\) 分配一下就可以得到 \(x^nF(\displaystyle\frac 1x) \equiv x^{n - m}Q(\frac 1x)x^mG(\frac 1x) + x^nR(\frac 1x) \pmod{x^n}\)。

整理一下可得 \(rF(x) \equiv rQ(x)rG(x) \pmod{x^{n - m + 1}}\)。

于是 \(rQ(x) \equiv rF(x)rG^{-1}(x)\)，通过多项式求逆和多项式乘法，就可以求出 \(rQ(x)\)，进一步求出 \(Q(x)\)，根据被除数等于商乘除数加余数，可以求出 \(R(x)\)。

完整代码：洛谷 P4512 【模板】多项式除法

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N];
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void inv(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	inv(f, g, n >> 1);
	int lim = 1;
	while(lim <= n * 2 + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
void div(int *f, int *g, int len1, int len2){
	int lim = 1;
	while(lim <= len1 + len2)
		lim <<= 1;
	NTT(f, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i <= lim; i++)
		f[i] = f[i] * g[i] % MOD;
	NTT(f, lim, -1);
}
int f[N << 1], g[N << 1], r[N << 1], q[N << 1], rf[N << 1], rg[N << 1], rr[N << 1], n, m;
signed main(){
	scanf("%lld%lld", &n, &m);
	for(int i = 0; i <= n; i++){
		scanf("%lld", &f[i]);
		rf[n - i] = f[i];
	}
	for(int i = 0; i <= m; i++){
		scanf("%lld", &g[i]);
		rg[m - i] = g[i];
	}
	for(int i = n - m + 1; i <= m; i++)
		rg[i] = 0;
	inv(rg, rr, n - m);
	div(rf, rr, n, n - m);
	for(int i = 0; i <= n - m; i++)
		printf("%lld ", q[i] = rf[n - m - i]);
	printf("\n");
	div(g, q, m, n - m);
	for(int i = 0; i < m; i++)
		printf("%lld ", r[i] = (f[i] - g[i] + MOD) % MOD);
	return 0;
}

求导/积分

\(G(x) = F'(x)\)

\(H(x) = \displaystyle\int F(x) \mathrm d x\)

求 \(G(x),H(x)\)

解：都已经 \(O(n)\) 了，还怎么优化！！！

完整代码：洛谷 U455122 【模板】多项式求导/积分

点击查看代码

开根

\(G^2(x) \equiv F(x) \pmod{x^n}\)。

求 \(G(x)\)。

解：先考虑初始情况，常数项 \(G_0(x) = \sqrt{F(0)}\) 是好求的，套入牛顿迭代的式子可以得到 \(G_{k + 1}(x) \equiv G_k(x) - \displaystyle\frac{G_k^2(x) - F(x)}{2G_k(x)} \pmod{x^{2n}}\)。

整理可得：\(G_{k + 1}(x) \equiv \displaystyle\frac{G_k(x)}{2} + \frac{F(x)}{2G_k(x)} \pmod{x^{2n}}\)，运用多项式除法即可求出答案。

完整代码：P5205 【模板】多项式开根

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 3e5 + 9, MOD = 998244353;
int rev[N], in[N];
void init(){
	in[1] = 1;
	for(int i = 2; i <= N; i++)
		in[i] = (MOD - MOD / i) * in[MOD % i] % MOD;
}
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void inv(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	inv(f, g, n >> 1);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	for(int i = n + 1; i < lim; i++)
		tmp[i] = 0;
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f2[N << 1];
void ln(int *f, int *g, int n){
	for(int i = 0; i < n * 4 + 3; i++)
		g[i] = 0;
	inv(f, g, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i < n; i++)
		f2[i] = f[i + 1] * (i + 1) % MOD;
	for(int i = n; i < lim; i++)
		f2[i] = 0;
	NTT(f2, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * f2[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n; i > 0; i--)
		g[i] = g[i - 1] * in[i] % MOD;
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
	g[0] = 0;
}
int lng[N << 1];
void exp(int *f, int *g, int n){
	if(!n){
		g[0] = 1;
		return;
	}
	exp(f, g, n >> 1);
	ln(g, lng, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++){
		if(f[i] >= lng[i])
			lng[i] = f[i] - lng[i];
		else
			lng[i] = f[i] - lng[i] + MOD;
	}
	for(int i = n + 1; i < lim; i++)
		lng[i] = g[i] = 0;
	lng[0]++;
	NTT(lng, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * lng[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f[N << 1], g[N << 1], res[N << 1], n, k;
char s[1009];
signed main(){
	init();
	scanf("%lld", &n);
	n--;
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	ln(f, g, n);
	int inv = qpow(2, MOD - 2);
	for(int i = 0; i <= n; i++)
		g[i] = g[i] * inv % MOD;
	exp(g, res, n);
	for(int i = 0; i <= n; i++)
		printf("%lld ", res[i]);
	return 0;
}

对数函数（\(\ln\)）

\(G(x) \equiv \ln F(x) \pmod{x^n}\)，保证 \(f_0 = 1\)。

求 \(G(x)\)。

解：将等式两边求导可得 \(G'(x) = F'(x)F^{-1}(x)\)，运用多项式乘法、求导、积分、求逆即可得出答案。

完整代码：洛谷 P4725 【模板】多项式对数函数（多项式 ln）

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N], in[N];
void init(){
	in[1] = 1;
	for(int i = 2; i <= N; i++)
		in[i] = (MOD - MOD / i) * in[MOD % i] % MOD;
}
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void inv(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	inv(f, g, n >> 1);
	int lim = 1;
	while(lim <= n * 2 + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f2[N << 1];
void ln(int *f, int *g, int n){
	inv(f, g, n);
	int lim = 1;
	while(lim <= n * 2 + 1)
		lim <<= 1;
	for(int i = 0; i < n; i++)
		f2[i] = f[i + 1] * (i + 1) % MOD;
	for(int i = n; i < lim; i++)
		f2[i] = 0;
	NTT(f2, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * f2[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n; i > 0; i--)
		g[i] = g[i - 1] * in[i] % MOD;
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
	g[0] = 0;
}
int f[N << 1], g[N << 1], n;
signed main(){
	init();
	scanf("%lld", &n);
	n--;
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	ln(f, g, n);
	for(int i = 0; i <= n; i++)
		printf("%lld ", g[i]);
	return 0;
}

指数函数（\(\exp\)）

\(G(x) \equiv e^{F(x)} \pmod{x^n}\)，保证 \(f_0 = 0\)。

求 \(G(x)\)。

解：等式两边同时取 \(\ln\) 可以得到 \(\ln G(x) - F(x) \equiv 0 \pmod{x^n}\)，初始值 \(\ln G(0) = 0\) 很容易求出，带入牛顿迭代的式子可得 \(G_{k + 1} \equiv G_k(x) + \displaystyle\frac{\ln G(x) - F(x)}{G^-1(x)} \pmod{x^{2n}}\)，整理可得 \(G_{k + 1}(x) \equiv G_k(x)(1 - \ln G_k(x) + F(x)) \pmod{x^{2n}}\)，套用多项式对数函数、乘法即可求出答案。

完整代码：洛谷 P4726 【模板】多项式指数函数（多项式 exp）

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N], in[N];
void init(){
	in[1] = 1;
	for(int i = 2; i <= N; i++)
		in[i] = (MOD - MOD / i) * in[MOD % i] % MOD;
}
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void inv(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	inv(f, g, n >> 1);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	for(int i = n + 1; i < lim; i++)
		tmp[i] = 0;
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f2[N << 1];
void ln(int *f, int *g, int n){
	for(int i = 0; i < n * 4 + 3; i++)
		g[i] = 0;
	inv(f, g, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i < n; i++)
		f2[i] = f[i + 1] * (i + 1) % MOD;
	for(int i = n; i < lim; i++)
		f2[i] = 0;
	NTT(f2, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * f2[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n; i > 0; i--)
		g[i] = g[i - 1] * in[i] % MOD;
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
	g[0] = 0;
}
int lng[N << 1];
void exp(int *f, int *g, int n){
	if(!n){
		g[0] = 1;
		return;
	}
	exp(f, g, n >> 1);
	ln(g, lng, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++){
		if(f[i] >= lng[i])
			lng[i] = f[i] - lng[i];
		else
			lng[i] = f[i] - lng[i] + MOD;
	}
	for(int i = n + 1; i < lim; i++)
		lng[i] = g[i] = 0;
	lng[0]++;
	NTT(lng, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * lng[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f[N << 1], g[N << 1], n;
signed main(){
	init();
	scanf("%lld", &n);
	n--;
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	exp(f, g, n);
	for(int i = 0; i <= n; i++)
		printf("%lld ", g[i]);
	return 0;
}

幂函数

\(G(x) \equiv F^k(x) \pmod{x^n}\)。

求 \(G(x)\)。

解：对等式两边取对数可得 \(\ln G(x) \equiv k \ln F(x) \pmod{x^n}\)，于是，求出 \(\ln F(x)\) 后，乘上 \(k\) 再 \(\exp\) 回去即可求出 \(G(x)\)。

完整代码：洛谷 P5245 【模板】多项式快速幂

点击查看代码

#include <bits/stdc++.h>
using namespace std;
#define int long long
const int N = 5e6 + 9, MOD = 998244353;
int rev[N], in[N];
void init(){
	in[1] = 1;
	for(int i = 2; i <= N; i++)
		in[i] = (MOD - MOD / i) * in[MOD % i] % MOD;
}
int qpow(int x, int y){
	int res = 1;
	while(y){
		if(y & 1)
			res = res * x % MOD;
		x = x * x % MOD;
		y >>= 1;
	}
	return res;
}
void NTT(int *f, int n, int opt){
	for(int i = 0; i < n; i++){
		rev[i] = rev[i >> 1] >> 1;
		if(i % 2 == 1)
			rev[i] |= n >> 1;
	}
	for(int i = 0; i < n; i++)
		if(i < rev[i])
			swap(f[i], f[rev[i]]);
	for(int h = 2; h <= n; h <<= 1){
		int mi = h >> 1;
		int gn = qpow(3, (MOD - 1) / h);
		for(int j = 0; j < n; j += h){
			int g = 1;
			for(int k = 0; k < mi; k++){
				int tmp = f[j + k + mi] * g % MOD;
				f[j + k + mi] = (f[j + k] - tmp + MOD) % MOD;
				f[j + k] = (f[j + k] + tmp) % MOD;
				g = g * gn % MOD;
			}
		}
	}
	if(opt == -1){
		reverse(f + 1, f + n);
	    int inv = qpow(n, MOD - 2);
	    for(int i = 0; i < n; i++)
			f[i] = f[i] * inv % MOD;
	}
}
int tmp[N << 1];
void inv(int *f, int *g, int n){
	if(!n){
		g[0] = qpow(f[0], MOD - 2);
		return;
	}
	inv(f, g, n >> 1);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++)
		tmp[i] = f[i];
	for(int i = n + 1; i < lim; i++)
		tmp[i] = 0;
	NTT(tmp, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = (2 - g[i] * tmp[i] % MOD + MOD) % MOD * g[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f2[N << 1];
void ln(int *f, int *g, int n){
	for(int i = 0; i < n * 4 + 3; i++)
		g[i] = 0;
	inv(f, g, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i < n; i++)
		f2[i] = f[i + 1] * (i + 1) % MOD;
	for(int i = n; i < lim; i++)
		f2[i] = 0;
	NTT(f2, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * f2[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n; i > 0; i--)
		g[i] = g[i - 1] * in[i] % MOD;
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
	g[0] = 0;
}
int lng[N << 1];
void exp(int *f, int *g, int n){
	if(!n){
		g[0] = 1;
		return;
	}
	exp(f, g, n >> 1);
	ln(g, lng, n);
	int lim = 1;
	while(lim < 2 * n + 1)
		lim <<= 1;
	for(int i = 0; i <= n; i++){
		if(f[i] >= lng[i])
			lng[i] = f[i] - lng[i];
		else
			lng[i] = f[i] - lng[i] + MOD;
	}
	for(int i = n + 1; i < lim; i++)
		lng[i] = g[i] = 0;
	lng[0]++;
	NTT(lng, lim, 1);
	NTT(g, lim, 1);
	for(int i = 0; i < lim; i++)
		g[i] = g[i] * lng[i] % MOD;
	NTT(g, lim, -1);
	for(int i = n + 1; i < lim; i++)
		g[i] = 0;
}
int f[N << 1], g[N << 1], res[N << 1], n, k;
char s[10000009];
signed main(){
	init();
	scanf("%lld%s", &n, s + 1);
	n--;
	int len = strlen(s + 1);
	for(int i = 1; i <= len; i++)
		k = ((k << 3) + (k << 1) + (s[i] - 48)) % MOD;
	for(int i = 0; i <= n; i++)
		scanf("%lld", &f[i]);
	ln(f, g, n);
	for(int i = 0; i <= n; i++)
		g[i] = g[i] * k % MOD;
	exp(g, res, n);
	for(int i = 0; i <= n; i++)
		printf("%lld ", res[i]);
	return 0;
}

三角函数

\(G(x) = \sin F(x)\)。

\(H(x) = \cos F(x)\)。

求 \(G(x), H(x)\)。

解：套用 \(\sin\) 与 \(\cos\) 与 \(e^x\) 的关系就可以用 \(\exp\) 得出答案。

完整代码：洛谷 P5264 多项式三角函数

点击查看代码

参考资料

Introductory Combinatorics(Fifth Edition) Richard A. Brualdi
利用生成函数求斐波那契数列通项公式 zwfymqz
生成函数封闭形式的正确理解 Sunlight-zero

posted @ 2025-01-20 16:52 Orange_new 阅读(36) 评论(0) 收藏举报

刷新页面返回顶部

JPGOJCZX

多项式学习笔记（二）：多项式初等函数

牛顿迭代

快速多项式计算

加法

乘法

求逆

除法/取模

求导/积分

开根

对数函数（\(\ln\)）

指数函数（\(\exp\)）

幂函数

三角函数

公告