这是一篇基于OI Wiki中的内容进行整理的算法竞赛中关于多项式的相关内容的第二篇笔记。
1.前置知识
1.1 阶
若正整数 a , p a,p a,p互质,且满足 p > 1 p>1 p>1,则对于 a n ≡ 1 ( m o d p ) a^n ≡ 1 \pmod p an≡1(modp)的 n n n,我们称这样的 n n n为 a a a模 p p p的阶,记作 δ p ( a ) \delta_p(a) δp(a)。
如: 2 2 2是 3 3 3模 8 8 8的阶,因为 3 3 3与 8 8 8互质,且 n = 2 n=2 n=2是最小的满足 3 n ≡ 1 ( m o d p ) 3^n ≡ 1\pmod p 3n≡1(modp)的值。
阶具有以下性质:
- 对于 i ∈ [ 0 , δ p ( a ) ) i∈[0,\delta_p(a)) i∈[0,δp(a)),所有的 a i ( m o d p ) a^i \pmod p ai(modp)结果互不相同。(可以用反证法证明)
- 若 a n ≡ 1 ( m o d p ) a^n≡1\pmod p an≡1(modp),则 δ p ( a ) ∣ n \delta_p(a) \vert n δp(a)∣n。( 3 2 ≡ 1 ( m o d 8 ) 3^2≡1\pmod 8 32≡1(mod8), 3 4 ≡ 1 ( m o d 8 ) 3^4≡1\pmod 8 34≡1(mod8))
- 若 p p p为质数,则 δ p ( g i ) = δ p ( g ) \delta_p(g^i)=\delta_p(g) δp(gi)=δp(g)的充要条件为 gcd ( δ p ( g ) , i ) = 1 \gcd(\delta_p(g),i)=1 gcd(δp(g),i)=1,其中 g g g为任意的正整数。
- δ p ( a b ) = δ p ( a ) gcd ( b , δ p ( a ) ) \delta_p(a^b)=\frac{\delta_p(a)}{\gcd(b,\delta_p(a))} δp(ab)=gcd(b,δp(a))δp(a),其中 p , a , b p,a,b p,a,b为任意正整数。
1.2 剩余类
剩余类
亦称同余类。设模为 n n n,则根据余数可将所有的整数分为 n n n类,把所有与整数 a a a模 n n n同余的整数构成的集合叫做模 n n n的一个剩余类,记作 [ a ] [a] [a]。
完全剩余系
从模 n n n的每一个剩余类中各取一个数,得到一个由 n n n个数组成的集合,叫做模 n n n的一个完全剩余系。
简化剩余系
是 m m m的完全剩余系中与 m m m互质的数构成的子集。
例如, m = 10 m=10 m=10的一个完全剩余系是 { 0 , 1 , 2 , 3 , 4 , 5 , 6 , 7 , 8 , 9 } \{0,1,2,3,4,5,6,7,8,9\} {0,1,2,3,4,5,6,7,8,9},这个完全剩余系中元素与 10 10 10互质的子集为 { 1 , 3 , 7 , 9 } \{1,3,7,9\} {1,3,7,9},这个子集就是 m = 10 m=10 m=10的简化剩余系。
1.3 原根
设 n n n是正整数, g g g是整数,如果 g g g模 n n n的阶 δ n ( g ) = φ ( n ) \delta_n(g)=\varphi(n) δn(g)=φ(n),且 gcd ( g , n ) = 1 \gcd(g,n)=1 gcd(g,n)=1,则称 g g g为 n n n的一个原根。
在抽象代数中,原根就是循环群的生成元。
原根具有以下性质:
若 gcd ( g , n ) = 1 \gcd(g,n)=1 gcd(g,n)=1且 n > 0 n>0 n>0,则 g g g为 n n n的一个原根的充要条件为: S = { g 1 , g 2 , g 3 , ⋯ , g φ ( n ) } S=\{g^1,g^2,g^3,\cdots,g^{\varphi(n)}\} S={g1,g2,g3,⋯,gφ(n)}为 n n n的一组简化剩余系。
原根判定定理
设 n ≥ 3 , gcd ( a , n ) = 1 n\geq 3,\gcd(a,n)=1 n≥3,gcd(a,n)=1,则 a a a是模 m m m的原根的充要条件是:对于 φ ( n ) \varphi(n) φ(n)的每个质因数 p p p,都有 a φ ( n ) p ≢ 1 ( m o d n ) a^{\frac{\varphi(n)}{p}}\not≡ 1\pmod n apφ(n)≡1(modn)。
原根个数
若一个数 m m m有原根,则它的原根的个数是 φ ( φ ( m ) ) \varphi(\varphi(m)) φ(φ(m))。
原根存在定理
一个数 m m m存在原根当且仅当 m = 2 , 4 , p α , 2 p α m=2,4,p^{\alpha},2p^{\alpha} m=2,4,pα,2pα,其中 p p p为奇素数, α ∈ N ∗ \alpha∈\mathbb{N}^* α∈N∗。
对于奇素数 p p p, p p p、 p α p^\alpha pα、 2 p α 2p^\alpha 2pα有原根 ( α ∈ N ∗ ) (\alpha∈\mathbb{N}^*) (α∈N∗)。
最小原根的数量级
若
m
m
m有原根,其最小原根是不多于
m
0.25
m^{0.25}
m0.25级别的。
暴力找一个数的最小原根的时间复杂度是可以接受的。
2.快速数论变换
快速傅里叶变换存在缺点,我们将多项式的系数,变成一个矩阵,乘上一个复数系数的矩阵加以处理,且每个复数系数的实部和虚部都是一个正弦/余弦函数,大部分的系数都是一个浮点数,不仅产生了大量的运算量,还产生了不小的运算误差。
数论变换解决的是多项式乘法带模数的情况,可以说有些受模数的限制,数也比较大。
数论变换(Number-theoretic transform, NTT)是通过将离散傅里叶变换化为 F / p \mathbb{F}/p F/p的形式,其中 p p p是一个整数模质数。这是一个有限域,只要 n n n可以除 p − 1 p-1 p−1,就存在本元 n n n次方根。具体来说,对于质数 p = q n + 1 ( n = 2 m ) p=qn+1(n=2^m) p=qn+1(n=2m),原根 g g g满足 g q n ≡ 1 ( m o d p ) g^{qn}≡1\pmod p gqn≡1(modp),将 g n = g q ( m o d p ) g_n=g^q\pmod p gn=gq(modp)看作 ω n \omega_n ωn的等价,则其满足相似的性质,比如 g n n ≡ 1 ( m o d p ) g_n^n≡1\pmod p gnn≡1(modp), g n n / 2 ≡ − 1 ( m o d p ) g_n^{n/2}≡-1\pmod p gnn/2≡−1(modp)。
区分FFT中的 n n n,我们将此处的 n n n称为 N N N。
涉及到数论变化,所以 N N N可以比FFT中的 n n n大,但是只需要将 q N n \frac{qN}{n} nqN看作这里的 q q q即可,就可以避免大小问题。
常见的有:
p
=
1004535809
=
479
×
2
2
1
+
1
,
g
=
3
p=1004535809=479\times 2^21+1,g=3
p=1004535809=479×221+1,g=3
p
=
998244353
=
7
×
17
×
2
23
+
1
,
g
=
3
p=998244353=7\times17\times2^{23}+1,g=3
p=998244353=7×17×223+1,g=3
NTT的实现方法很简单,利用到结论:
$
ω
n
≡
g
p
−
1
n
(
m
o
d
p
)
\omega_n≡g^{\frac{p-1}{n}}\pmod p
ωn≡gnp−1(modp)
将FFT中的 ω n \omega_n ωn全部替换掉即可。
p = 998244353 p=998244353 p=998244353时,原根 g = 3 g=3 g=3。
LuoguP3803 - 【模板】多项式乘法(FFT)
NTT代码
#include<bits/stdc++.h>
using namespace std;
typedef long long LL;
const double PI = acos(-1.0);
const LL mod = 998244353;
const LL g = 3, gi = 332748118;
//gi = power(3ll, mod - 2);
template<class T>T power(T a, LL b) {
T res = 1;
for (; b; b >>= 1) {
if (b & 1) res = (res * a) % mod;
a = a * a % mod;
}
return res;
}
LL rev[4000005];
void change(LL *y, int len) {
for (int i = 0; i < len; ++i) {
if (i < rev[i]) swap(y[i], y[rev[i]]);
}
}
void NTT(LL *y, int len, int on) {
change(y, len);
for (LL h = 2; h <= len; h <<= 1) {
LL wn = power((on == 1) ? g : gi, (mod - 1) / h);//dif
for (int j = 0; j < len; j += h) {
LL w = 1ll;
for (int k = j; k < j + h / 2; ++k) {
LL u = y[k] % mod;
LL t = w * y[k + h / 2] % mod;
y[k] = (u + t) % mod;
y[k + h / 2] = (u - t + mod) % mod;
w = w * wn % mod;
}
}
}
}
LL x1[3000005], x2[3000005];
int n, m;
void main2() {
cin >> n >> m;
int len = 1;
while (len <= (n + m)) len <<= 1;
for (int i = 0; i < len; ++i) {
rev[i] = rev[i >> 1] >> 1;
if (i & 1) rev[i] |= len >> 1;
}
for (int i = 0; i <= n; ++i) {
cin >> x1[i];
x1[i] = (x1[i] + mod) % mod;
}
for (int i = n + 1; i < len; ++i) {
x1[i] = 0;
}
for (int i = 0; i <= m; ++i) {
cin >> x2[i];
x2[i] = (x2[i] + mod) % mod;
}
for (int i = m + 1; i < len; ++i) {
x2[i] = 0;
}
NTT(x1, len, 1); NTT(x2, len, 1);
for (int i = 0; i < len; ++i) {
x1[i] = (x1[i] * x2[i]) % mod;
}
NTT(x1, len, -1);
LL inv = power((LL)len, mod - 2);
for (int i = 0; i <= n + m; ++i) {
cout << (x1[i] * inv) % mod << ' ';
}
}
int main() {
ios::sync_with_stdio(false);
cin.tie(0); cout.tie(0);
LL _ = 1;
// cin >> _;
while (_--) main2();
return 0;
}
对三种多项式乘法的实现方式进行一个比较。
FTT:2.13s
FTT(3to2):1.59s
NTT:1.80s