多项式全家桶-优快云博客

记号与约定

用 $F(x)=\sum\limits_{i=0} f_i x^i$ 表示一个多项式，用 $x^n]F(x)=f_n$ 表示 $F (x)$ 第 $n$ 项系数。

设 $n$ 为最大的非零的 $f_n$ ，那么 $F (x)$ 的次数为 $n$ ，记作 $\deg(F(x))=n$ ，注意 $n$ 次多项式有 $n + 1$ 项。

多项式加法： $F(x)+G(x)=\sum\limits_{i=0}(f_i+g_i)x^i$
多项式乘法： $F(x)\times G(x)=\sum\limits_{i=0}(\sum\limits_{j=0}^if_jg_{i-j})x^i$

一、多项式乘法

一切的基础

我们希望快速计算两个多项式 $F$ 和 $G$ 的卷积

考虑对于固定的 $x_0$ ， $F(x_0)G(x_0)=(F\times G)(x_0)$ 。

因为 $n + 1$ 个值可以确定一个 $n$ 次多项式，那么我们构造 $x_0\sim x_n$ ，分别求出 $F(x_0),F(x_1),...,F(x_n)$ 和 $G(x_0),G(x_1),...,G(x_n)$ ，这也称为点值表示法。

那么只需要求出所有 $(F\times G)(x_i)=F(x_i)G(x_i)$ ，然后就能插值还原。

复杂度瓶颈在于求点值表示和还原多项式。

方便起见我们认为项数为 $2^k$ 即 $f$ 长度是 $2$ 的幂，这个可以通过在 $F$ 和 $G$ 的高项补 $0$ 系数实现，下面令 $n=2^k$ 。

1.1 快速傅里叶变换（FFT）

快速傅里叶变换

用来处理一般的卷积

考虑构造 $n$ 次单位根 $x_i=w_n^{i}$ 作为点值，因为单位根有很好的性质：

$w_{2n}^{2i}=w_n^i$
$w_{n}^{i+n/2}=-w_n^i$

从复数的几何意义角度不难理解

那么 $\forall k=0,1\dots \frac{n}{2}-1$ ：

$F(w_n^k)=\sum\limits_{i=0}^{n-1}f_i(w_n^k)^i$

$=\sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i}(w_n^k)^{2i} + \sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i+1}(w_n^k)^{2i+1}$

$=\sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i}(w_n^{2k})^{i} + w_n^k\sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i+1}(w_n^{2k})^{i}$

$=\sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i}(w_{n/2}^{k})^{i} + w_n^k\sum\limits_{i=0}^{\frac{n}{2}-1}f_{2i+1}(w_{n/2}^{k})^{i}$

令 $F_0,F_1$ 分别为 $F$ 中只保留偶数、奇数项系数的多项式：

$F(w_n^k)=F_0(w_{n/2}^k)+w_n^kF_1(w_{n/2}^k)$

同理可以得到：

$F(w_n^{k+n/2})=F_0(w_{n/2}^k)-w_n^kF_1(w_{n/2}^k)$

我们只需要递归求出 $F_0(x_{0\sim n/2-1})$ 和 $F_{1}(x_{0\sim n/2-1})$ 即可。

这样复杂度是 $T(n)=T(n/2)+O(n)=O(n\log n)$

快速傅里叶逆变换

求出了所有点值后，我们还需要还原系数

快速傅里叶变换的过程可以看成是：

$\begin{bmatrix} 1 & 1 & \dots & 1&1\\ 1 & w_n & \dots &w_n^{n-2} & w_n^{n-1}\\ \dots \\ 1 & w_n^{n-1}& \dots &w_n^{(n-1)(n-2)} & w_n^{(n-1)(n-1)} \\ \end{bmatrix} \begin{bmatrix} f_0 \\ f_1 \\ \dots \\ f_{n-1} \end{bmatrix}= \begin{bmatrix} y_0 \\ y_1 \\ \dots \\ y_{n-1} \end{bmatrix}$

那么我们想要还原系数只需要乘上左边矩阵的逆即可

这个矩阵被称为 “范德蒙德矩阵”，它的逆很特殊，为：

$\frac{1}{n}\begin{bmatrix} 1 & 1 & \dots & 1&1\\ 1 & 1/w_n & \dots & 1/w_n^{n-2} & 1/w_n^{n-1}\\ \dots \\ 1 & 1/w_n^{n-1}& \dots & 1/w_n^{(n-1)(n-2)} & 1/w_n^{(n-1)(n-1)} \\ \end{bmatrix}$

即每项分别取倒数，然后再乘上 $1/ n$ ，那么只需要将上面过程中的所有 $w_n$ 替换成 $1/w_n$ 即可，具体实现我们可以只写一个 FFT 函数，然后传一个 opt 代表是 $w_n$ 还是 $1/w_n$

蝴蝶变换

递归版常数较大，通常我们会使用循环写法来代替

注意到我们实际上的运算过程是从下至上依次合并的，那么只需要求出最底层的各项 $f'_0,f'_1\dots f'_{n-1}$ 分别对应原来的哪一项

注意到 $f_i'=f_{rev(i)}$ ，其中 $re v (i)$ 代表 $i$ 的二进制高低位倒置后的数这步操作叫做蝴蝶变换

变换后按照分治的过程用循环实现即可

#include<bits/stdc++.h>
using namespace std ;

typedef long long LL ;
const int N = 4e6+10 ;
const double Pi = acos(-1.0) ;

int n , m , r[N] ;
struct complx
{
	double a , b ;
	complx (double xx=0, double yy=0 ) {
		a = xx , b = yy ; 
	}
	friend complx operator + ( const complx x , const complx y ) {
		return {x.a+y.a,x.b+y.b} ;
	}
	friend complx operator - ( const complx x , const complx y ) {
		return {x.a-y.a,x.b-y.b} ;
	}
	friend complx operator * ( const complx x , const complx y ) {
		return {x.a*y.a-x.b*y.b,x.a*y.b+x.b*y.a} ;
	}
}F[N] , G[N] ;
void FFT( complx *a , int l , int opt ) // [0,2^l)
{
	for(int i = 0 ; i < (1<<l) ; i ++ ) {
		if( i<r[i] ) swap(a[i],a[r[i]]) ;
	}
	for(int i = 1 ; i <= l ; i ++ ) {
		int p = (1<<(i-1)) ;
		complx w1(cos(Pi/p),opt*sin(Pi/p)) ;
		for(int j = 0 ; j < (1<<l) ; j += (1<<i) ) {
			complx w(1,0) ;
			for(int k = 0 ; k < p ; k ++ , w=w*w1 ) {
				complx A = a[j+k] , B = w*a[j+k+p] ;
				a[j+k] = A+B ;
				a[j+k+p] = A-B ;
			}
		}
	}
}

int main()
{
	scanf("%d%d" , &n , &m ) ;
	int x ;
	for(int i = 0 ; i <= n ; i ++ ) {
		scanf("%d" , &x ) ;
		F[i].a = x ;
	}
	for(int i = 0 ; i <= m ; i ++ ) {
		scanf("%d" , &x ) ;
		G[i].a = x ;
	}
	int len = 1 , l = 0 ;
	while( len<n+m+1 ) len<<=1 , l ++ ;
	for(int i = 0 ; i < len ; i ++ ) r[i] = (r[i>>1]>>1)|((i&1)<<(l-1)) ;  
	FFT(F,l,1) ; FFT(G,l,1) ;
	for(int i = 0 ; i < len ; i ++ ) F[i] = F[i]*G[i] ;
	FFT(F,l,-1) ;
	for(int i = 0 ; i <= n+m ; i ++ ) {
		printf("%d " , int(F[i].a/len+0.5) ) ;
	}
	return 0 ;
}

1.2 快速数论变换 (NTT)

FFT 在精度和常数等方面有很大的问题，并不十分常用

而当我们把数域限定为 $Z_p$ 即 $\{0,1,2\dots p-1\}$ （ $p$ 为素数）时，就有一个更常用的等效算法：快速数论变换

首先我们需要找到一种能够代替单位根 $w_n$ 的东西

考虑对于一个非零数 $g$ ， $g,g^2,g^3\dots$ 在 $\bmod {p}$ 意义下一定存在循环，如果循环节为 $p - 1$ ，我们就称 $g$ 为 $p$ 的一个原根，所有素数均有原根。

已知 $w_n^n=1$ ， $g^{p-1}=1$ ，那么我们就考虑用 $g^{(p-1)/n}$ 来代替 $n$ 次单位根 $w_n^i$ ，很容易说明这个东西有我们需要的单位根的所有性质。

注意这里 $n$ 是我们补齐后的 $2$ 的整次幂，需要保证 $n ∣ p - 1$ （也就是 $p - 1$ 的分解中 $2$ 的指数比较大），这样的模数 $p$ 也被称为 NTT 模数，最常见的就是 $998244353=1+2^{23}\times119$ ，其原根 $g = 3$

$w_n^1=g^{(mod-1)/n}$

怎么找原根呢？见 oi-wiki

原根判定定理

设 $m\geq 3$ ，则 $g$ 是模 $m$ 意义下的原根的充要条件是，对于 $\varphi(m)$ 的每个素因子 $p$ ，都有 $g^{\frac{\varphi(m)}{p}}\ne1\pmod m$

这样我们就可以暴力枚举 + check 找到某个模数的原根了

写法和 FFT 几乎没有区别

#include<bits/stdc++.h>
using namespace std ;

typedef long long LL ;
const int N = 4e6+100 , mod = 998244353 , G = 3 , invG = (mod+1)/G ; 
LL ksm( LL a , LL b )
{
	LL res = 1 , t = a % mod ;
	while( b ) {
		if( b&1 ) res = res * t % mod ;
		b = b >> 1 ;
		t = t * t % mod ;
	}
	return res ;
}

int n , m , r[N] ;
LL F[N] , H[N] ;
void NTT( LL *a , int l , bool fg )
{
	for(int i = 0 ; i < (1<<l) ; i ++ ) {
		if( i < r[i] ) swap(a[i],a[r[i]]) ;
	}
	for(int i = 1 ; i <= l ; i ++ ) {
		int p = (1<<(i-1)) ;
		LL w1 = ksm( fg?invG:G , (mod-1)/(1<<i) ) ;
		for(int j = 0 ; j < (1<<l) ; j += (1<<i) ) {
			LL w = 1 ;
			for(int k = 0 ; k < p ; k ++ , w=w*w1%mod ) {
				LL B = a[j+k+p]*w%mod ;
				a[j+k+p] = ( a[j+k]-B ) % mod ;
				a[j+k] = ( a[j+k] + B ) % mod ; 
			}
		}
	}
}

int main()
{
	scanf("%d%d" , &n , &m ) ;
	for(int i = 0 ; i <= n ; i ++ ) {
		scanf("%lld" , &F[i] ) ;
	}
	for(int i = 0 ; i <= m ; i ++ ) {
		scanf("%lld" , &H[i] ) ;
	}
	int len = 1 , l = 0 ;
	while( len < n+m+1 ) len<<=1 , l ++ ;
	for(int i = 0 ; i < len ; i ++ ) r[i] = (r[i>>1]>>1)|((i&1)<<(l-1)) ;
	NTT(F,l,0) ; NTT(H,l,0) ;
	for(int i = 0 ; i < len ; i ++ ) F[i] = F[i]*H[i]%mod ;
	NTT(F,l,1) ;
	LL iv = ksm(len,mod-2) ;
	for(int i = 0 ; i <= n+m ; i ++ ) {
		printf("%lld " , (F[i]*iv%mod+mod)%mod ) ;
	}
	return 0 ;
}