多项式学习-完全版

持续更新中...

start in 14 July 2026

卷积变换

快速求 $$c_k = \sum _{i \oplus j = k} a_i b_j$$

关于卷积…… 如何通俗易懂地解释卷积? - 知乎 https://www.zhihu.com/question/22298352/answer/228543288

FFT -> NTT -> FMT/FWT

https://www.luogu.com.cn/training/1008#information

https://blog.messywind.top/posts/cp生成函数多项式题单/

快速傅里叶变换 FFT

题目大意:对多项式 \(F(x), G(X)\),以 \(O(n \log n)\) 的时间复杂度求得 \(F(x) G(x)\) 的乘积。

https://www.luogu.com.cn/article/jm6yxjpx
https://www.cnblogs.com/pks-t/p/9251147.html

普通的做法,直接相乘,时间复杂度是 \(O(n^2)\) 的,难以接受。

\(F(X), G(X)\) 用点值表示,这样就可以用 \(O(n)\) 来相乘了。而将系数表达式与点值表达式的互相转换,时间复杂度为 \(O(n \log n)\)

如下是非递归写法。

#include <stdio.h>
#include <algorithm>
#include <math.h>

const double PI = acos(-1.0);

struct CPX {
	double x, y; // x+yi 
	CPX (double x = 0.0, double y = 0.0) : x(x), y(y)  { }
	CPX operator + (const CPX &t) { return { x + t.x, y + t.y }; }
	CPX operator - (const CPX &t) { return { x - t.x, y - t.y }; }
	CPX operator * (const CPX &t) { return { x * t.x - y * t.y , x * t.y + y * t.x }; }
} a[3000003], b[3000003], c[3000003]; 	// note that N >= 2^upper(logn)

int n, m, len, rev[3000003];

inline void fft (CPX a[], double md) { // md: 1: FFT, -1: IFFT
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]); 
	for (int l = 1; l < len; l <<= 1) { // l: half of the len of the interval now
		CPX W = { cos(PI / l), md * sin(PI / l) }; 
		for (int i = 0; i < len; i += (l << 1)) {	// for every interval...
			CPX w(1.0, 0.0);
			for (int j = 0; j < l; ++j) {	// for every node in the interval now...
				CPX u = a[i + j];
				CPX v = w * a[i + j + l];
				a[i + j] = u + v;
				a[i + j + l] = u - v;
				w = w * W;
			}
		}
	} 
}

inline void kagari () {
	scanf("%d %d", &n, &m);
	for (int i = 0; i <= n; ++i) scanf("%lf", &a[i].x);
	for (int i = 0; i <= m; ++i) scanf("%lf", &b[i].x);
	
	len = 1;
	while (len <= n + m) len <<= 1;
	
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	
	fft(a, 1);
	fft(b, 1);
	
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i];
	
	fft(c, -1);
	
	for (int i = 0; i <= n + m; ++i) 
		printf("%d ", (int)(c[i].x / len + 0.5));
	puts("");
	return;
}
int main () {
	kagari();
	return 0;
}

此外还有一个有趣的优化(from https://www.luogu.com.cn/article/pvwh9hpp ):
\(G(x)\) 放入 \(F(x)\) 的虚部中,这时 \((F+Gi)^2=(F^2-G^2)+2FGi\),虚部除以 \(2\) 即为答案。如此,只需做两次傅里叶变换(先对 \(F+Gi\) 做 FFT,开方后,再对 \(FG\) 做 IFFT)而非三次,在常数级别大大优化了时间复杂度。

原根

https://blog.csdn.net/zhang20072844/article/details/11541133

假设一个数 \(G\) 对于质数 \(P\) 来说是原根,那么 \(G^i \mod P\)的结果两两不同,且有 \(1<G<P, 1<i<P\),那么 \(G\) 可以称为是 \(P\) 的一个原根。\(P\) 的原根 $g $ 的 \(1,2,\dots,P-1\) 次幂 \(\mod P\) 的结果一定互不相同。

要求最小原根,则可以从 \(2\) 开始枚举 \(G\),直到 \(G^{P-1} \equiv 1 \pmod P\)

p.s. 不是每个素数都存在原根。\(p=k\times2^a+1\)

\[p = 167772161 = 5 \times 2^{25}+1, g=3 \]

\[p = 469762049 = 7 \times 2^{26}+1, g=3 \]

\[p = 754974721 = 3^2 \times 5 \times 2^{24}+1, g=11 \]

\[p = 998244353 = 7 \times 17 \times 2^{23}+1, g=3 \]

\[p = 1004535809 = 479 \times 2^{21}+1, g=3 \]

快速数论变换 NTT

https://zhuanlan.zhihu.com/p/347726949

NTT与多项式全家桶

在 FFT 的基础上,用原根代替单位根,以解决取模与 double 的精度问题。

那么显然 NTT 就难以处理 double 类型的数据了。

#include <stdio.h>
#include <algorithm>
#define ll long long

const ll MODN = 998244353;
const ll G = 3; 			// primitive root of MODN
const ll GI = 332748118; 	// inv of MODN 

inline ll power (ll x, ll k) {
	ll res = 1;
	while (k) {
		if (k & 1) res = res * x % MODN;
		x = x * x % MODN;
		k >>= 1;
	}
	return res;
}

int n, m, len, rev[3000003];
ll a[3000003], b[3000003], c[3000003];

inline void ntt (ll a[], ll md) { // md: 1: NTT, -1: INTT
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]); 
	for (int l = 1; l < len; l <<= 1) { 
		ll W = power(md == 1 ? G : GI, (MODN - 1) / (l << 1));
		for (int i = 0; i < len; i += (l << 1)) {
			ll w = 1;
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j];
				ll v = w * a[i + j + l] % MODN;
				a[i + j] = (u + v) % MODN;
				a[i + j + l] = (u - v + MODN) % MODN;
				w = w * W % MODN;
			}
		}
	} 
}

inline void kagari () {
	scanf("%d %d", &n, &m);
	for (int i = 0; i <= n; ++i) scanf("%lld", &a[i]);
	for (int i = 0; i <= m; ++i) scanf("%lld", &b[i]);
	
	len = 1;
	while (len <= n + m) len <<= 1;
	
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	
	ntt(a, 1);
	ntt(b, 1);
	
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	
	ntt(c, -1);
	
	ll invlen = power(len, MODN - 2);
	for (int i = 0; i <= n + m; ++i) 
		printf("%lld ", c[i] * invlen % MODN);
	puts("");
	return;
}
int main () {
	kagari();
	return 0;
}

快速莫比乌斯变换 FMT | 快速沃尔什变换 FWT

https://www.luogu.com.cn/article/2pavj2pd
https://www.luogu.com.cn/article/idd90f6t

FWT 是用于解决对下标进行位运算卷积问题的方法,与 FFT 类似。即从求 \(C_i= \Sigma_{i=j+k} A_j B_k\) 变为求 \(C_i= \Sigma_{i=j\oplus k} A_j B_k\),其中 \(\oplus\) 是二元位运算 and, or, xor 中的一个。

以 P4717 【模板】快速莫比乌斯 / 沃尔什变换 (FMT / FWT) 为例:

题目大意:给定长度为 \(2^n\) 两个序列 \(A,B\),设 \(C_i=\sum_{j\oplus k = i}A_j \times B_k\),分别当 \(\oplus\) 是 or, and, xor 时求出 \(C\),并对 \(998244353\) 取模。

对 or, and, xor 三种操作,代码差别只在于对 \(u, v\) 的处理区别。实际应用上在套用模板时,也很少会改变模板中的代码。

#include <stdio.h>
#define ll long long
#define MODN 998244353 
#define inv2 499122177
int n, len;
ll a[2000003], b[2000003], c[2000003];

inline void fwtor (ll a[], ll md) {
	for (int l = 1; l < len; l <<= 1) 
		for (int i = 0; i < len; i += (l << 1)) 
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j];
				ll v = a[i + j + l];
				if (md == 1)
					a[i + j + l] = (u + v) % MODN;
				else
					a[i + j + l] = (v - u + MODN) % MODN;
			}
}

inline void fwtand (ll a[], ll md) {
	for (int l = 1; l < len; l <<= 1) 
		for (int i = 0; i < len; i += (l << 1)) 
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j];
				ll v = a[i + j + l];
				if (md == 1)
					a[i + j] = (u + v) % MODN;
				else
					a[i + j] = (u - v + MODN) % MODN;
			}
}

inline void fwtxor (ll a[], ll md) {
	for (int l = 1; l < len; l <<= 1) 
		for (int i = 0; i < len; i += (l << 1)) 
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j];
				ll v = a[i + j + l];
				if (md == 1)
					a[i + j] = (u + v) % MODN,
					a[i + j + l] = (u - v + MODN) % MODN;
				else
					a[i + j] = (u + v) % MODN * inv2 % MODN,
					a[i + j + l] = (u - v + MODN) % MODN * inv2 % MODN;
			}
}

inline void kagari () {
	scanf("%d", &n);
	len = (1 << n);
	for (int i = 0; i < len; ++i) scanf("%lld", &a[i]);
	for (int i = 0; i < len; ++i) scanf("%lld", &b[i]);
	
	fwtor(a, 1), fwtor(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	fwtor(a, -1), fwtor(b, -1), fwtor(c, -1);
	for (int i = 0; i < len; ++i) printf("%lld ", c[i]);
	puts("");
	
	fwtand(a, 1), fwtand(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	fwtand(a, -1), fwtand(b, -1), fwtand(c, -1);
	for (int i = 0; i < len; ++i) printf("%lld ", c[i]);
	puts("");
	
	fwtxor(a, 1), fwtxor(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	fwtxor(a, -1), fwtxor(b, -1), fwtxor(c, -1);
	for (int i = 0; i < len; ++i) printf("%lld ", c[i]);
	puts("");

	return;
}
int main () {
	kagari();
	return 0;
}

P1919 【模板】高精度乘法 / A*B Problem 升级版

题目大意:求 \(a \times b\),其中 \(a, b \leq 10^{10^6}\)

\(a, b\) 拆成数组,直接利用 NTT 即可。是 FFT 与 NTT 的经典简单应用。

const ll MODN = 998244353; 
const ll G = 3;
const ll GI = 332748118; 

int n, m, len, rev[8000003];
ll a[8000003], b[8000003], c[8000003];

inline void ntt (ll a[], ll md) {
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]);
	for (int l = 1; l < len; l <<= 1) {
		ll W = power(md == 1 ? G : GI, (MODN - 1) / (l << 1));
		for (int i = 0; i < len; i += (l << 1)) {
			ll w = 1;
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j], v = w * a[i + j + l] % MODN;
				a[i + j] = (u + v) % MODN;
				a[i + j + l] = (u - v + MODN) % MODN;
				w = w * W % MODN;
			}
		}
	} 
}

inline void nttmul (ll a[], ll b[], ll c[]) {
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	ntt(a, 1), ntt(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	ntt(c, -1);
	ll invlen = power(len, MODN - 2);
	for (int i = 0; i < len; ++i) c[i] = c[i] * invlen % MODN;
}

inline void kagari () {
	std:: string p, q; std:: cin >> p >> q; n = p.size(), m = q.size();
	std:: reverse(p.begin(), p.end());
	std:: reverse(q.begin(), q.end());
	for (int i = 0; i < n; ++i) a[i] = p[i] - '0';
	for (int i = 0; i < m; ++i) b[i] = q[i] - '0';
	len = 1; while (len <= n + m) len <<= 1;
	
	nttmul(a, b, c);
	
	std:: string ans = "";
	ll r = 0;
	for (int i = 0; i < len; ++i) {
		c[i] += r;
		r = c[i] / 10;
		c[i] %= 10;
		ans.push_back(c[i] + '0');
	}
	while (ans.size() && ans.back() == '0') ans.pop_back();
	std:: reverse(ans.begin(), ans.end());
	std:: cout << ans << '\n';
	return;
}

牛顿迭代法 | Newton's Method

在此介绍牛顿迭代法。

回到高等数学中,先考虑如何求一个非线性方程 \(f(x)=0\) 的解 \(x\)

我们任取一个点 \(x_0\),这时,\(x_1 = x_0 - \dfrac{f(x_0)}{f'(x_0)}\) 一定使得 \(f(x_1)\)\(f(x_0)\) 更趋于 \(0\)

以此类推,利用递推式 \(x_{n+1} = x_n - \dfrac{f(x_n)}{f'(x_n)}\),无论 \(x_0\) 如何取值,重复足够多次后都可以取得近似解。

再将其应用于多项式中

在多项式运算里,我们往往是在模 \(x^n\) 意义下,求解一个未知的多项式 \(B(x)\)。这可以改写为:

\[F(B(x)) \equiv 0 \pmod{x^n} \]

假设我们已知 \(\mod{x^{\left\lfloor \frac{n}{2} \right\rfloor}}\) 时的 \(B_0(x)\)

此时借用牛顿迭代公式,将 \(B(x)\) 当成变量 \(x_n\)

\[B(x) \equiv B_0(x) - \frac{F(B_0(x))}{F'(B_0(x))} \pmod{x^n} \]

上式证明略。

\(B(x) \equiv 0 \pmod x^0\) 一般是极易求得的。

如此,只需倍增枚举取模次数,就可以 \(O(n \log n)\) 递推解决问题了。

应用

然而在求多项式 \(\ln\) 时是不需要使用牛顿迭代法的。在求多项式逆元 \(inv\),计算 \(exp\)等情况时则会用到。

多项式求逆为例,我们要找 \(A(x)\) 的逆元 \(B(x)\),也就是满足 \(A \cdot B \equiv 1\),移项构造方程:

\[F(B) = \frac{1}{B} - A \equiv 0 \]

未知多项式 \(B\) 求导,把 \(A\) 看作常数:

  • 函数:\(F(B) = B^{-1} - A\)
  • 导数:\(F'(B) = -B^{-2} = -\frac{1}{B^2}\)

现在,把 \(F(B_0)\)\(F'(B_0)\) 塞进多项式牛顿迭代公式里:

\[B \equiv B_0 - \frac{\frac{1}{B_0} - A}{-\frac{1}{B_0^2}} \pmod{x^{2k}} \]

化简即得:

\[B \equiv B_0 \cdot (2 - A \cdot B_0) \pmod{x^{2k}} \]

注意!!!!!!!!!!!!!!!!!

时刻注意清零,不要将污染了的数据带入下一重计算。

同时注意数组大小,防止越界与 MLE

P4238 【模板】多项式乘法逆

题目大意:给定一个 \(n \leq 10^5\) 次多项式 \(F(x)\) ,求一个多项式 \(G(x)\), 使得 \(F(x) * G(x) \equiv 1 \pmod{x^n}\),系数对 \(998244353\) 取模。

通过数学关系推导,转化为使用 FFT/NTT 的方法。

已知:$$A B \equiv 1 \pmod{x^{\left\lfloor \frac{n}{2} \right\rfloor}}$$

而又有:$$A B' \equiv 1 \pmod{x^{\left\lfloor \frac{n}{2} \right\rfloor}}$$

相减得:$$B' - B \equiv 0 \pmod{x^{\left\lfloor \frac{n}{2} \right\rfloor}}$$

平方,得:$$B'^2 - 2BB' + B^2 \equiv 0 \pmod{x^n}$$

等式左右同乘 \(A\):$$A B'^2 - 2B' + B \equiv 0 \pmod{x^n}$$

移项得到递推式:$$B \equiv 2B' - A B'^2 \pmod{x^n}$$

也就是求 $$ f_n = 2f_{n-1}-A \cdot f_{n-1}^2$$

初始状态 \(f_0 = inv(A_0)\)

由主定理知,此方法时间复杂度为 \(O(n \log n)\)

注意在算 \(A B \equiv 1 \pmod{x^k}\)\(B\) 时,\(B\) 所有 \(>k\) 次的系数均为 \(0\),否则会发生错误。

inline void ntt (ll a[], ll md) {
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]);
	for (int l = 1; l < len; l <<= 1) {
		ll W = power(md == 1 ? G : GI, (MODN - 1) / (l << 1));
		for (int i = 0; i < len; i += (l << 1)) {
			ll w = 1;
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j], v = w * a[i + j + l] % MODN;
				a[i + j] = (u + v) % MODN, a[i + + j + l] = (u - v + MODN) % MODN;
				w = w * W % MODN;
			}
		}
	}
	if (md == -1) {
		ll invlen = power(len, MODN - 2);
		for (int i = 0; i < len; ++i) a[i] = a[i] * invlen % MODN;
	} 
}

inline void kagari () {
	scanf("%d", &n);
	for (int i = 0; i < n; ++i) scanf("%lld", &a[i]);
	
	b[0] = inv(a[0]);
	for (int k = 1; k < (n << 1); k <<= 1) {
		for (int i = 0; i <= (k >> 1); ++i) d[i] = b[i];
		len = 1;
		while (len <= (k << 1)) len <<= 1;
		for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
		
		for (int i = (k >> 1) + 1; i < len; ++i) b[i] = 0;
		for (int i = 0; i <= k; ++i) c[i] = a[i]; // note: regard ai = 0 when i > k
		for (int i = k + 1; i < len; ++i) c[i] = 0;
		ntt(c, 1), ntt(b, 1);
		for (int i = 0; i < len; ++i) b[i] = c[i] * b[i] % MODN * b[i] % MODN;
		ntt(b, -1);
		for (int i = 0; i < len; ++i) b[i] = (d[i] * 2ll % MODN - b[i] + MODN) % MODN;
	} 
	for (int i = 0; i < n; ++i) printf("%lld ", b[i]);
	puts("");
	return;
}

P4725 【模板】多项式对数函数(多项式 ln)

要求 $$ G = \ln F \pmod{998244353}$$

对两边求导,有

\[G' = \dfrac{F'}{F} \pmod{998244353} \]

也就是求

\[G = \int F' \cdot F^{-1} \pmod{998244353} \]

其中对 \(F\) 求导与对右边积分的时间复杂度都是 \(O(n)\),求 \(F\) 的逆与卷积都是 \(O(n \log n)\)

#include <stdio.h>
#include <algorithm>
#define ll long long

const ll MODN = 998244353;
const ll G = 3;
const ll GI = 332748118;
inline ll power (ll x, ll k) {
	ll res = 1;
	while (k) {
		if (k & 1) res = res * x % MODN;
		x = x * x % MODN;
		k >>= 1;
	}
	return res;
}
inline ll inv (const ll &x) {
	return power(x % MODN, MODN - 2);
}

int n, len, rev[800003];
ll a[800003], b[800003], c[800003], d[800003], e[800003], f[800003];

inline void ntt (ll a[], ll md) {
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]);
	for (int l = 1; l < len; l <<= 1) {
		ll W = power(md == 1 ? G : GI, (MODN - 1) / (l << 1));
		for (int i = 0; i < len; i += (l << 1)) {
			ll w = 1;
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j], v = w * a[i + j + l] % MODN;
				a[i + j] = (u + v) % MODN, a[i + + j + l] = (u - v + MODN) % MODN;
				w = w * W % MODN;
			}
		}
	}
	if (md == -1) {
		ll invlen = power(len, MODN - 2);
		for (int i = 0; i < len; ++i) a[i] = a[i] * invlen % MODN;
	} 
}

inline void inv (ll a[], ll b[]) {
	b[0] = inv(a[0]);
	for (int k = 1; k < (n << 1); k <<= 1) {
		for (int i = 0; i <= (k >> 1); ++i) d[i] = b[i];
			
		len = 1;
		while (len <= (k << 1)) len <<= 1;
		for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
		
		for (int i = (k >> 1) + 1; i < len; ++i) b[i] = 0;
		for (int i = 0; i <= k; ++i) c[i] = a[i]; // note: regard ai = 0 when i > k
		for (int i = k + 1; i < len; ++i) c[i] = 0;
		ntt(c, 1), ntt(b, 1);
		for (int i = 0; i < len; ++i) b[i] = c[i] * b[i] % MODN * b[i] % MODN;
		ntt(b, -1);
		for (int i = 0; i < len; ++i) b[i] = (d[i] * 2ll % MODN - b[i] + MODN) % MODN;
	} 
	for (int i = n; i < len; ++i) b[i] = 0; // 注意 i>n 时 bi 应为 0 
}

inline void qiudao (ll a[], ll b[]) {
	for (int i = 1; i < n; ++i) b[i - 1] = a[i] * (ll)i % MODN;
}

inline void jifen (ll a[], ll b[]) {
	for (int i = 0; i < n; ++i) b[i + 1] = a[i] * inv(i + 1) % MODN;
}

inline void nttmul (ll a[], ll b[], ll c[]) {
	len = 1;
	while (len <= (n << 1)) len <<= 1;
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	ntt(a, 1), ntt(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	ntt(c, -1);
}

inline void kagari () {
	scanf("%d", &n);
	for (int i = 0; i < n; ++i) scanf("%lld", &a[i]);
	
	qiudao(a, e);
	inv(a, b); 
	nttmul(e, b, b);
	jifen(b, f);
	
	for (int i = 0; i < n; ++i) printf("%lld ", f[i]);
	puts("");
	return;
}
int main () {
	kagari();
	return 0;
}

P4726 【模板】多项式指数函数(多项式 exp)

题目大意:给定一个 \(n-1\) 次多项式 \(A(x)\),其中常数项 \(a_0 = 0\),求多项式 \(B(x)\),使得:\(B(x) \equiv e^{A(x)} \pmod{x^n}\)

使用牛顿迭代方法。

要求 $$B(x) \equiv e^{A(x)} \pmod{x^n}$$

两边同时取对数,移项,则构造出一个关于 \(B(x)\) 的函数方程:

\[F(B(x)) = \ln B(x) - A(x) \equiv 0 \pmod{x^n} \]

由牛顿迭代公式:

\[B(x) \equiv B_0(x) - \frac{F(B_0(x))}{F'(B_0(x))} \equiv B_0(x) - \frac{\ln B_0(x) - A(x)}{\frac{1}{B_0(x)}} \equiv B_0(x) \cdot (1 - \ln B_0(x) + A(x)) \pmod{x^n} \]

于是递推即可。注意常数较大,有时需优化。

P5245 【模板】多项式快速幂

题目大意:给定一个 \(n-1\) 次多项式 \(A(x)\),其中 \(a_0=1\),求一个在 \(\bmod\ x^n\) 意义下的多项式 \(B(x)\),使得 \(B(x) \equiv (A(x))^k \ (\bmod\ x^n)\)。系数对 \(998244353\) 取模。

\(B(x) \equiv (A(x))^k \ (\bmod\ x^n)\) 两边取对数:

\[\ln B(x) \equiv k \ln A(x) \pmod {x^n} \]

那么 \(e^{k \ln A(x)}\) 即为所求。

注意 \(k\) 可能很大,直接将其对 \(998244353\) 即可。

P5205 【模板】多项式开根

题目大意:给定一个 \(n-1\) 次多项式 \(A(x)\),其中 \(a_0=1\),求一个在 \(\bmod x^n\) 意义下的多项式 \(B(x)\),使得 \(B^2(x) \equiv A(x) \pmod{x^n}\)。系数对 \(998244353\) 取模。

\(F(x) = B^2-A\),则有

\[B_0 \equiv B - \dfrac{B_0^2 - A}{2B_0} \equiv \frac{1}{2} (B_0-A B_0^{-1}) \pmod{x^k} \]

套上 NTT 与逆元即可。

多项式大模板

函数 空间大小
$$ntt(a, b, N)$$ $$4N$$
\(ntt(a, b, N)\) \(4N\)

例题

P3723 [AHOI2017/HNOI2017] 礼物

FFT/NTT 模板题。

将式子改写为:

\[\sum_{i=1}^n \left(x_{i} - y_{i} + c\right)^{2} \equiv \sum_{i=1}^{n} x_{i}^{2} + \sum_{i=1}^n y_{i}^{2} + 2c \sum_{i=1}^n (x_i-y_i) + n c^2 - 2 \sum_{i=1}^n x_iy_i \]

左边两项是定值,中间两项通过枚举 \(c \in [-m, m] \in [-100, 100]\) (表示第一个手环加上的值)即可得最小值。于是,除了最右边的项全都是定值。

所以我们只需要求 k 取何值时, \(\sum_{i=1}^n x_i y_{i+k}\) 最大就行了。

翻转 \(x\),式子变为 \(\sum_{i=1}^n x_{n-i+1} y_{i+k}\),这是熟悉的卷积的形式,使用 NTT 进行求解。

const ll MODN = 998244353;
const ll G = 3;
const ll GI = 332748118;

int n, len, rev[140003];
ll m, a[140003], b[140003], c[140003];

inline void ntt (ll a[], ll md) {
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]);
	for (int l = 1; l < len; l <<= 1) {
		ll W = power(md == 1 ? G : GI, (MODN - 1) / (l << 1));
		for (int i = 0; i < len; i += (l << 1)) {
			ll w = 1;
			for (int j = 0; j < l; ++j) {
				ll u = a[i + j], v = w * a[i + j + l] % MODN;
				a[i + j] = (u + v) % MODN;
				a[i + j + l] = (u - v + MODN) % MODN;
				w = w * W % MODN;
			}
		}
	}
}

inline void nttmul () {
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	
	ntt(a, 1), ntt(b, 1);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i] % MODN;
	ntt(c, -1);
	ll invlen = power(len, MODN - 2);
	for (int i = 0; i < len; ++i) c[i] = c[i] * invlen % MODN; 
}

inline void kagari () {
	scanf("%d %lld", &n, &m);
	for (int i = 1; i <= n; ++i) scanf("%lld", &a[i]);
	for (int i = 1; i <= n; ++i) scanf("%lld", &b[i]), b[i + n] = b[i];
	
	ll res = 0, ans = MODN, r = 0;
	for (int i = 1; i <= n; ++i) r += a[i] - b[i];
	for (ll i = -m; i <= m; ++i) ans = std:: min(ans, i * i * (ll)n + i * 2ll * r);
	for (int i = 1; i <= n; ++i) ans += a[i] * a[i] + b[i] * b[i];
	
	len = 1;
	while (len <= (n << 1)) len <<= 1; 
	for (int i = 1; i <= (n >> 1); ++i) std:: swap(a[i], a[n - i + 1]);
	
	nttmul();
	
	for (int i = n + 1; i <= (n << 1); ++i) res = std:: max(res, c[i]);
	ans = (ans - 2ll * res + MODN) % MODN;
	printf("%lld\n", ans);
	return;
}

P3338 [ZJOI2014] 力

较简单的题目。

题目大意:给出 \(n\) 个数 \(a_1,a_2, \dots a_n, n \leq 10^5\),对 \(1 \leq i \leq n\),求 \(E_i=\dfrac{F_i}{q_i}\),其中\(F_j~=~\sum_{i = 1}^{j - 1} \dfrac{a_i \times a_j}{(i - j)^2}~-~\sum_{i = j + 1}^{n} \dfrac{a_i \times a_j}{(i - j)^2}\)

即求

\[F_j=\sum_{i = 1}^{j - 1} \dfrac{a_j}{(i - j)^2}~-~\sum_{i = j + 1}^{n} \dfrac{a_j}{(i - j)^2} \]

联想到 FFT 的形式 \(C_k= \sum A_i \times B_{k-i}\),令 \(b_i=\dfrac{1}{i^2}\),则有

\[F_j= \sum_{i = 1}^{j - 1} a_j \times b_{i-j} - \sum_{i = j + 1}^{n} a_j \times b_{i-j} \]

注意到上式的左式虽易得,但右式并不好求,下标 \(i - j\) 有时为负,不妨扩大 \(b\) 的范围,将上式与 \(b_i\) 改写为:

\[F_j= \sum_{i = 1}^{j - 1} a_j \times b_{i-j+n} + \sum_{i = j + 1}^{n} a_j \times b_{i-j+n} \\ b_{n+i} = \begin{cases} \dfrac{1}{i^2}, & 0<i\leq n \\ 0, & i = 0 \\ - \dfrac{1}{i^2}, & -n \leq i < 0 \end{cases} \]

答案就是 \(F_n+1\)\(F_{2n}\)

#include <stdio.h>
#include <algorithm>
#include <math.h>

const double PI = acos(-1.0);

struct CPX {
	double x, y; // x+yi 
	CPX (double x = 0.0, double y = 0.0) : x(x), y(y)  { }
	CPX operator + (const CPX &t) { return { x + t.x, y + t.y }; }
	CPX operator - (const CPX &t) { return { x - t.x, y - t.y }; }
	CPX operator * (const CPX &t) { return { x * t.x - y * t.y , x * t.y + y * t.x }; }
} a[800003], b[800003], c[800003], d[800003], e[800003], f[800003]; 	// note that N >= 2^upper(logn)

int n, rev[800003];

inline void fft (CPX a[], double md, int len) { 
	for (int i = 0; i < len; ++i) if (i < rev[i]) std:: swap(a[i], a[rev[i]]); 
	for (int l = 1; l < len; l <<= 1) { 
		CPX W = { cos(PI / l), md * sin(PI / l) }; 
		for (int i = 0; i < len; i += (l << 1)) {	
			CPX w(1.0, 0.0);
			for (int j = 0; j < l; ++j) {
				CPX u = a[i + j];
				CPX v = w * a[i + j + l];
				a[i + j] = u + v;
				a[i + j + l] = u - v;
				w = w * W;
			}
		}
	} 
}

inline void fftmul (CPX a[], CPX b[], CPX c[], int N) {
	int len = 1; while (len <= (N << 1)) len <<= 1;
	for (int i = 0; i < len; ++i) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) ? (len >> 1) : 0);
	fft(a, 1, len);
	fft(b, 1, len);
	for (int i = 0; i < len; ++i) c[i] = a[i] * b[i];
	fft(c, -1, len);
	for (int i = 0; i < len; ++i) c[i].x /= len, c[i].y = 0;
} 

inline void kagari () {
	scanf("%d", &n);
	for (int i = 1; i <= n; ++i) scanf("%lf", &a[i].x);
	for (int i = 1; i <= n; ++i) b[n + i].x = 1.0 / i / i;
	for (int i = -1; i >= -n; --i) b[n + i].x = -1.0 / i / i;
	b[n].x = 0.0;
	
	fftmul(a, b, c, n << 1);
	
	for (int i = n + 1; i <= n * 2; ++i) printf("%f\n", c[i].x);	
	return;
}
int main () {
	kagari();
	return 0;
}

P5488【模板】差分与前缀和

题目大意:给定一个长为 \(n\) 的序列 \(a\),求出其 \(k\) 次前缀和与差分(即前缀和/差分 \(k\) 次后的结果)。系数对 \(1004535809\) 取模。

在此以生成函数的方式来说明此题:

  • 先考虑如何求前缀和 \(S(x)\)

令生成函数 \(A(x) = a_0+a_1x+a_2x^2+\dots+a_{n-1}x^{n-1}\),进行一次前缀和的结果为 \(s_0, s_1, \dots, s_{n-1}\)

对前缀和函数 \(S(x) = s_0+s_1x+s_2x^2+\dots+s_{n-1}x^{n-1}\)

两边乘以 \(x\),有 \(xS(x) = s_0x+s_1x^2+s_2x^3+\dots+s_{n-1}x^n\)

相减,有 \((1-x)S(x) = a_0+a_1x+a_2x^2+\dots+a_{n-1}x^{n-1} = A(X)\)

于是 \(S(X)=\dfrac{A(X)}{1-x}\)。则 \(k\) 次前缀和就是 \(S_k(x) = \dfrac{A(X)}{(1-x)^k}=A(x)\times(1-x)^{-k}\)

  • 现在考虑差分函数 \(C(x) = c_0+c_1x+c_2x^2+\dots+c_{n-1}x^{n-1}\),其中 \(c_i = a_i - a_{i-1}\)

思考如何凑出 \(a_i - a_{i-1}\),考虑到对 \(A(x)\) 进行操作

\[(1-x)A(x) = \sum_{0 \leq i < n} (a_i-a_{i-1})x^i = C(i) \]

那么要求 \(k\) 次差分,即求 \(C_k(x) = A(x)\times(1-x)^k\)

综上,我们套用 NTT,多项式求逆与快速幂即可。

事实上,前缀和的式子也可以是 \(S(x) = A(x) \times (1+x+\dots+x^{n-1})\),而且感觉这个式子是更容易想到的。

组合数学方式的进一步优化,可提供一种时间优化方式 by gemini

直接去算多项式幂 \((1-x)^K\)\((1-x)^{-K}\),常数极大,效率很低。但这道题里的多项式只有两项 \((1-x)\),我们可以直接利用二项式定理展开!

假设我们要求的多项式为 \(B(x) = \sum_{i=0}^{N-1} b_i x^i\)

  • 情况一:\(K\) 阶前缀和 —— 计算 \((1-x)^{-K}\)

根据负整数次幂的广义二项式定理展开:

\[(1-x)^{-K} = \sum_{i=0}^{\infty} \binom{-K}{i} (-x)^i = \sum_{i=0}^{\infty} \binom{K+i-1}{i} x^i \]

其第 \(i\) 项系数 \(b_i = \binom{K+i-1}{i}\)
观察相邻两项的关系:

\[\frac{b_i}{b_{i-1}} = \frac{\binom{K+i-1}{i}}{\binom{K+i-2}{i-1}} = \frac{K+i-1}{i} \]

已知 \(b_0 = 1\),我们可以用 \(O(N)\) 的线性递推直接求出整个 \(B(x)\) 多项式:

\[b_i = b_{i-1} \cdot \frac{K+i-1}{i} \pmod{P} \]

  • 情况二:\(K\) 阶差分 —— 计算 \((1-x)^K\)

根据二项式定理展开:

\[(1-x)^K = \sum_{i=0}^{\infty} \binom{K}{i} (-1)^i x^i \]

其第 \(i\) 项系数 \(b_i = (-1)^i \binom{K}{i}\)
观察相邻两项的关系:

\[\frac{b_i}{b_{i-1}} = \frac{(-1)^i \binom{K}{i}}{(-1)^{i-1} \binom{K}{i-1}} = -\frac{K-i+1}{i} = \frac{i-1-K}{i} \]

已知 \(b_0 = 1\),同样可以用 \(O(N)\) 线性递推:

\[b_i = b_{i-1} \cdot \frac{i-1-K}{i} \pmod{P} \]

inline void sum (const vell &a, vell &b, const std:: string &K, int N) { // K 次前缀 
	vell c(N << 3, 0), d(N << 3), e(N << 3);
	for (int i = 0; i < N; ++i) d[i] = 1;
//	or: c[0] = 1, c[1] = MODN - 1; inv(c, d, N);
	power(d, e, K, N);
	nttmul(a, e, b, N);
}

inline void fini (const vell &a, vell &b, const std:: string &K, int N) { // k 次差分 
	vell c(N << 3, 0), d(N << 3);
	c[0] = 1, c[1] = MODN-1; // 注意不是 -1 而是 MODN - 1!! 
	power(c, d, K, N);
	nttmul(a, d, b, N);
}

P5641 【CSGRound2】开拓者的卓识

题目大意:给定序列 \(a\),记一个序列 \([l,r]\)\(k\) 阶子段和为 \(sum_{k,l,r}\),有\(sum_{k,l,r}=\begin{cases}\sum\limits_{i=l}^{r}a_i&,k=1\\\sum\limits_{i=l}^{r}\sum\limits_{j=i}^{r}sum_{k-1,i,j}&,k\geq 2\end{cases}\)。对所有 \(1 \leq i \leq n\),求 \(sum_{k,1,i}\)。答案对 \(998244353\) 取模。

from https://www.luogu.com.cn/article/bev8cd54 :

考虑 \(a_i\) 对答案的贡献。

\(\mathit{sum}_{1,l,r}\) 就是从 \(a_l\) 加到 \(a_r\)

\(\mathit{sum}_{2,l,r}\) 根据 \(\mathit{sum}_{k,l,r}\) 的公式,可以看作是选择一个区间 \([l_1,r_1] \subseteq [l,r]\),再将 \(a_{l_1}\)\(a_{r_1}\) 的和加入答案中。

\(\mathit{sum}_{3,l,r}\) 是选择一个区间 \([l_1,r_1] \subseteq [l,r]\),再对 \([l_1,r_1]\) 求一遍 \(\mathit{sum}_{2,l_1,r_1}\),也就是再选择一个区间 \([l_2,r_2] \in [l_1,r_1]\)\(a_{l_2}\)\(a_{r_2}\) 的和加入答案中。

以此类推,\(\mathit{sum}_{k,l,r}\) 就是选择 \(k-1\) 个区间 \([l_{k-1},r_{k-1}] \subseteq [l_{k-2},r_{k-2}] \subseteq \dots \subseteq [l_1,r_1] \subseteq [l,r]\),并将 \([l_k,r_k]\) 的和加入答案中。

考虑转换贡献体,枚举 \(a_i\) 会对 \(\mathit{sum}_{k,1,r}\) 产生多大的贡献,也就是有多少组区间 \([l_i,r_i]\) 满足 \(i \in [l_{k-1},r_{k-1}] \subseteq [l_{k-2},r_{k-2}] \subseteq \dots \subseteq [l_1,r_1] \subseteq [l_0=1,r_0=r]\),由于是相邻区间之间是包含关系,故 \(i \ge l_{k-1} \ge l_{k-2} \ge \dots \ge l_1 \ge 1\),记 \(d_k = i - l_{k-1},d_{k-1} = l_{k-1} - l_{k-2},\dots,d_2 = l_2 - l_1,d_1 = l_1 - 1\),那么有 \(d_1 + d_2 + d_3 + \dots + d_k = i - 1\),其中 \(d_i \ge 0\),采用隔板法可求得所有 \((l_1,l_2,\dots,l_{k-1})\) 的个数为\(\binom{i + k - 2}{k - 1}\)。同理可得所有 \((r_1,r_2,\dots,r_{k-1})\) 的个数为\(\binom{r - i + k - 1}{k - 1}\)

\(a_i\)\(r\) 产生的贡献为

\[a_i \times \binom{i + k - 2}{k - 1} \times \binom{r - i + k - 1}{k - 1} \]

也就是要对每个 \(1 \leq r \leq n\)

\[\sum_{i=1}^r a_i \times \binom{i + k - 2}{k - 1} \times \binom{r - i + k - 1}{k - 1} \]

注意到其中类似卷积的形式,记

\[A_i = a_i \times \binom{i + k - 2}{k - 1}, B_i = \binom{i + k - 2}{k - 1} \]

那么

\[\mathit{sum}_{k,1,r}=C_{r+1} = \sum_{i=1}^r A_i\times B_{r-i+1} = \sum_{i=1}^r a_i \times \binom{i + k - 2}{k - 1} \times \binom{r - i + k - 1}{k - 1} \]

计算组合数 \(B_i\) ,有递推式

\[B_{i+1} = \binom{m+i}{m} = \dfrac{(m+i)!}{m!\ i !} = \dfrac{(m+i-1)!}{m!(i-1)!} \times \dfrac{m+i}{i} = \binom{m+i-1}{m} \times \dfrac{m+i}{i}= B_i\times \dfrac{m+i}{i} \]

其中 \(m = k-1\)\(B_1=\binom{k-1}{k-1}=1\)

搞清楚了会发现代码很短。

inline void kagari () {
	int n; ll k; scanf("%d %lld", &n, &k);
	vell a(n << 3), b(n << 3), c(n << 3);
	for (int i = 1; i <= n; ++i) scanf("%lld", &a[i]);
	
	b[1] = 1; 
	for (int i = 1; i < n; ++i) b[i + 1] = b[i] * (ll)(k - 1 + i) % MODN * inv(i) % MODN;
	for (int i = 1; i <= n; ++i) a[i] = a[i] * b[i] % MODN;
	
	nttmul(a, b, c, n + 1);
	
	for (int i = 1; i <= n; ++i) printf("%lld ", c[i + 1]);
	puts("");
	return;
}
posted @ 2026-07-14 16:00  dbg_8  阅读(6)  评论(0)    收藏  举报