多项式全家桶

一、前置芝士

1. 基本概念

  • 多项式: 有限项相加的求和式 anxn,记作 f(x)=anxn
  • 多项式的度(次数):对于一个多项式 f(x),称其最高次项的次数为该多项式的度(degree),也称次数,记作 degf
  • 级数:将数列的项依次用加号连接起来的函数
  • 幂级数:每项均为非负整数次幂函数乘以常数系数的级数称为幂级数
更多&拓展

对于多项式,我们也可以将其直接定义为一个系数序列,表示为:

f(x)=f0+f1x+f2x2++fnxn

此处我们认为,x 只是一个形式符号,即一个对于系数位置的标识符

若支持无穷项,则称 f(x)形式幂级数

f(x)=f0+f1x+f2x2+


Q:多项式和级数的区别?

A:多项式是有限的,级数是无限的

e.g. i=0anxn 是幂级数但不是多项式

2. 多项式的表示方法

① 系数表示法

f(x)=a0+a1x+a2x2++an1xn1=i=0naixi

存储方式

存储时只存储每一位的系数即可

f(x)={a0,a1,a2,,an1}

② 点值表示法

即把多项式放入平面直角坐标系中,得到一个函数图像

存储方式

通常地,带入 nx,得到 ny

f(x)={(x0,f(x0)),(x1,f(x1)),(x2,f(x2)),,(xn1,f(xn1))}


计算多项式乘法时:点值表示法可以 O(n) 求,即每一位对应乘起来即可

f(x)g(x)={(x0,f(x0)g(x0)),(x1,f(x1)g(x1)),(x2,f(x2)g(x2)),,(xn1,f(xn1)g(xn1))}

而这两种表示法可以相互转化,但是暴力转化是 O(n2) 的,那么我们就要使用一个神奇的 FFT 优化多项式乘法,并且它是 O(nlogn)

3. 复指数函数 & 复三角函数

  • 复指数函数:对于复数 z=x+i×y,定义 ez 为其复指数函数。特别地,我们通常将 ez 记作 expz
    在代数上 expz=ex(cosy+isiny)
证明

前置知识:欧拉定理

expz=ex+i×y=exei×y=ex(cosy+isiny)

  • 复三角函数

cosz=exp(iz)+exp(iz)2

sinz=exp(iz)exp(iz)2i

4. 复数的三种形式

详细见 复数 - ricky_lin

① 代数形式

z=x+yi

适用范围:计算复数的加减乘除

② 三角形式

z=r(cosθ+isinθ)

③ 指数形式

z=rexp(iθ)

适用范围:这两种形式用于计算复数的乘除两个运算以及后面的运算较为方便

三角形式 知识层面要求更低,但 指数形式 会更加地方便

5. 单位根

本小节需要掌握的前置知识点

① 定义和约定

我们定义 xn=1 在复数意义下的解是 n 次复根

n 个解都为 n单位根

根据学过的知识 n 次单位根将单位圆分成 n 等份

图片来源-bilibili 3blue1brown 【官方双语】那么……什么是卷积? 截图(截图和草图)

ωn=exp2πin=cos2πn+sin2πni

那么 xn=1 解集表示为 {ωnk|k=0,1,,n1}.

ωnk=exp2kπin=cos(2kπn)+sin(2kπn)i

一般地,我们默认 n 次单位根 ωn 表示:从 1 开始逆时针方向的第一个解,即上方的 ωn

② 性质和证明

  • 性质1:ωnn=1
证明 由定义 xn=1 可证
  • 性质2:ωnk=ω2n2k
证明 ω2n2k=exp2×2k×πi2n=exp2kπin=ωnk
  • 性质3:ω2nk+n=ω2nk
证明

我们先需要知道一件事:

exp(πi)=eπi=1

然后:

ω2nn+k=exp2(n+k)πi2n=exp(2nπi2n+2kπi2n)=exp(πi)+exp(2kπi2n)=ω2nk

二、多项式乘法

上面我们讲到多项式若用点值表示,那么乘法是 O(n) 的,但是显然地,我们常用系数进行表示

然而,对于系数转点值,需要代入 n 个点,显然复杂度是 O(n2)

然后我们考虑,找一些 x 使得 xm=1,这样的话就会计算 m 次乘方之后进入循环,不需要做全部的 m 次运算了

但是,显然地,这样的 x 在实数域内只有可能是 ±1,所以我们要找到一些特殊性质的数。

1. FFT——快速傅里叶变换

① 递归

第一种方法是将 x 拓展到复数域

f(x)=(a0+a2x2++an2xn2)+x(a1+a3x2++an1xn2)

f1(x)=a0+a2x++an2xn21f2(x)=a1+a3x++an1xn21(保证 x 为偶数)

可以得到:f(x)=f1(x2)+xf2(x2)

代入 x=ωnk 得:

f(ωnk)=f1(ωn2k)+ωnkf2(ωn2k)=f1(ωn2k)+ωnkf2(ωn2k)

代入 ωnk+n2=ωnk 得:

f(ωnk+n2)=f1(ωn2k+n)+ωnk+n2f2(ωn2k+n)=f1(ωn2k)ωnkf2(ωn2k)=f1(ωn2k)ωnkf2(ωn2k)

然后便可以递归求解,复杂度 O(nlogn)

重要提示

考虑上面的式子的适用条件是每一次多项式的长度都是 2 倍数,那么原先的多项式的长度就必须是 2 的次幂,不满足的在后面补 0 即可

② 倍增

但是考虑递归常数过大,我们考虑倍增怎么做

递归的时候需要把系数按照 {a0,a2,a4,a6,},{a1,a3,a5,a7,}来分组,并交换位置,我们将其称为 位逆序变换

例子

8 项多项式为例:

状态 序列
初始序列 {a0,a1,a2,a3,a4,a5,a6,a7}
第一次二分 {a0,a2,a4,a6},{a1,a3,a5,a7}
第二次二分 {a0,a4},{a2,a6},{a1,a5},{a3,a7}
第三次二分 {a0},{a4},{a2},{a6},{a1},{a5},{a3},{a7}

根据上面的例子,我们可以找到一个规律:

初始序列的下标在二进制下 翻转便可以得到其在最终序列中的下标

例子

x1 为例,1 在二进制下是 001(1) 最终序列中是 100(4),同样的 x4 也跑到了下标为 1 的位置

位逆序变换实现

位逆序置换可以 O(n) 从小到大递推实现:

我们设 R(x) 表示长度为 k 的二进制数 x 翻转后的数(高位补 0)。我们要求的是 R(0),R(1),,R(n1)

首先 R(0)=0

我们从小到大求 R(x)。因此在求 R(x) 时,R(x2) 的值是已知的。因此我们把 x 右移一位(除以 2),然后翻转,再右移一位,就得到了 x 除了(二进制)个位之外
其它位的翻转结果,把最后一位 0/1 拼到最高位,就得到了 R(x)

数学语言:

x=((x >> 1) << 1) | (x & 1)   

R(x)=(R(x >> 1) >> 1) | ((x & 1) << len)   

③ 蝶形运算优化

我们从上面知道,我们需要通过下面两个式子进行推导:

f(ωnk)=f1(ωn2k)+ωnkf2(ωn2k)

f(ωnk+n2)=f1(ωn2k)ωnkf2(ωn2k)

使用 位逆序置换 后,对于给定的 n,k

  • f1(ωn2k) 的值存储在下标为 k 的位置上,f2(ωn2k) 的值存储在下标为 n2+k 的位置上
  • f(ωnk) 的值存储在数组下标为 kf(ωnk+n2) 的值将存在数组下标为 k+n2 的位置

所以说我们每次变换的时候不需要再新开一个数组,只需要在原数组上进行覆写,即可。

这样我们就完成了 DFT——离散傅里叶变换

④ IDFT——快速傅里叶逆变换

在前面,我们将两个多项式用 O(nlogn) 的复杂度,变成了点值表示法,并用 O(n) 的时间,求出了乘积。

现在,我们需要将点值表示法在 O(nlogn) 的时间内,变为系数表示法。

方法即为:将 ωnk共轭复数 ωn0 带回到 DFT 中,再对其所有系数都除以 len(多项式长度)

证明

我们知道 ωnk共轭复数ωnk

同时我们也知道了一个函数 f(x)ωn0,ωn1,,ωnn1 的点值 g(x),并设为 y0,y1,,yn1,我们需要求出原函数的系数 a0,a1,,an1

g(x)=i=0n1yixi,yi=j=0n1aj(ωni)j

g(x)=i=0n1j=0n1aj(ωni)jxi

代入 x=ωnk 得到:

g(ωnk)=i=0n1j=0n1aj(ωni)j(ωnk)i=j=0n1aji=0n1(ωni)jk=j=0n1aji=0n1(ωnjk)i

S(x)=i=0n1(ωnx)i,由等比数列求和得:

S(x)=(ωnx)n(ωnx)0ωnx1=(ωnn)x1ωnx1=0ωnx1

由上式得到,当 ωnx=1x=0 时,S(0)=i=0n11i=n,否则 S(x)=0

代回原式:

g(ωnk)=j=0n1aji=0n1(ωnjk)i=j=0n1aji=0n1S(jk)i

g(ωnk)=nak

code
#include<cstdio>
#include<vector>
#include<cmath>
using namespace std;
typedef long long ll;
const double PI = acos(-1);
struct Complex{
double x,y;
Complex(double x = 0,double y = 0): x(x),y(y){};
Complex operator + (const Complex &B){return Complex(x+B.x,y+B.y);}
Complex operator - (const Complex &B){return Complex(x-B.x,y-B.y);}
Complex operator * (const Complex &B){return Complex(x*B.x-y*B.y,x*B.y+y*B.x);}
};
typedef vector<Complex> Poly;
Poly A,B;
int len = 1,llog;
void FFT(vector<int> &rev,Poly &v,int inv){
for(int i = 0; i < len; ++i)
if(i < rev[i]) swap(v[i],v[rev[i]]);
for(int k = 1; k < len; k <<= 1){
Complex omega(cos(PI/k), inv * sin(PI/k));
for(int i = 0; i < len; i += k * 2){
Complex w(1,0);
for(int j = 0; j < k; ++j, w = w * omega){
Complex s = v[i+j], t = v[i+j+k] * w;
v[i+j] = s + t; v[i+j+k] = s - t;
}
}
}
// for(int i = 0; i < len; ++i) printf("%lf ",v[i].x);puts("");
if(inv == -1) for(int i = 0; i < len; ++i) v[i].x /= len;
}
int n,m;
int main() {
scanf("%d%d",&n,&m);
while(len <= n + m) len <<= 1,++llog;
A.resize(len);B.resize(len);
for(int i = 0; i <= n; ++i) scanf("%lf",&A[i].x);
for(int i = 0; i <= m; ++i) scanf("%lf",&B[i].x);
vector<int> rev(len);
for(int i = 0; i < len; ++i) rev[i] = (rev[i>>1] >> 1) | ((i&1) << (llog-1));
FFT(rev, A, 1), FFT(rev, B, 1);
for(int i = 0; i < len; ++i) A[i] = A[i] * B[i];
FFT(rev, A, -1);
for(int i = 0; i <= n + m; ++i) printf("%lld ",(ll)(A[i].x + 1e-3));
return 0;
}

⑤ 从线性代数角度理解FFT

DFT 本身是个线性变换,可以理解为将目标多项式当作向量,左乘一个矩阵得到变
换后的向量,以模拟把单位复根代入多项式的过程:

[y0y1y2y3yn1]=[ωn0ωn0ωn0ωn0ωn0ωn0ωn1ωn2ωn3ωnn1ωn0ωn2ωn4ωn6ωn2(n1)ωn0ωn3ωn6ωn9ωn3(n1)ωn0ωnn1ωn2(n1)ωn3(n1)ωn(n1)2][a0a1a2a3an1]

A     =                                           B                             ×    C

现在我们已经得到最左边的矩阵 A 了,中间的 x 值在目标多项式的点值表示中也是一一对应的

所以,根据矩阵的基础知识,我们只要在式子两边左乘中间那个大矩阵 B 的逆矩阵就行了。

由于这个矩阵的元素非常特殊,它的逆矩阵也有特殊的性质,就是每一项取倒数,再除以变换的长度 n,就能得到它的逆矩阵。

注:这里的长度 n 也指一个 2 的幂,不足则补 0

2. NTT——快速数论变换

由于 FFTω 的时候精度极有可能爆炸,导致一些神奇问题,所以我们需要一些精度更高的算法,而且我们更多时候是考虑的取模意义下的运算。

于是我们便有了 FFT

① 补充亿点芝士

  • a,p 互素,且 p>1,对于 an1(modp) 最小n,我们称之为 ap,记作 δp(a),例如 δ7(2)=3

性质:a,a2,,aδp(a)p 两两不同余

证明

反证法:

 ij  i,jδp(a)  aiaj(modp),则有 a|ij|1(modp)

显然地,0<|ij|<δp(a),这和阶的定义矛盾,故原命题成立。

原根

pZ+,aZ,若 δp(a)=φ(p),则称 a 为模 p 的一个原根,记为 g

δ7(3)=6=φ(7),因此 3 是模 7 的一个原根。

换一个角度理解

就是说 aφ(p)1modpiZ,1i<φ(p) 满足 ai1modp,则 ap 的原根

可以得到一个结论

ωn=gp1n

原根判定定理 & 原根计算

m3gcd(g,m)=1,则 g 是模 m 的原根的充要条件是:

对于 φ(m) 的每个素因数 p,都有原根。

没什么快速求法,只能暴力枚举判断。若素数 p 有原根,其最小原根 gp=O(p0.25+ϵ)

其中 ϵ>0,这保证了我们暴力找一个数的最小原根,复杂度是可以接受的。

通常模数常见的有 998244353,1004535809,469762049。这几个的原根都是 3

我们用原根可以代替 FFT 中的 ω,因为它满足虚数 ω 的所有性质(模意义下)

证明
  • 性质1:ωn0=ωnn=1

显然地 ωn0=g0=1ωnn=1

根据费马小定理 gp11(modp),得证 ωnngp1(modp),得证。

  • 性质2:ωnk=ω2n2k

由定义得:ω2n2k=(gp12n)2k=gk(p1)n=(gp1n)k=ωnk

  • 性质3:ωnk+n2=ωnk

(ωnn2)2=gp11modp

又因为 pprime,所以 g0p1modp 两两不同

所以 gp121modp

即证:gk+p12=gk

在实现 NTT 的时候,把 DFT 中的 ωnk=cos(kn2π)+sin(kn2π) 替换为 gk(p1)nmodpIDFT 就是将 gg1 即可。


例题:P3723 [AH2017/HNOI2017] 礼物

给定两个长为 n 的数列 a, b,a, b 均可循环移动,你需要选择一个整数 c,最小化

i=1n(aibi+c)2

n5×104,ai,bi100

solution 把柿子拆一拆: i=1n(aibi+c)2=i=1nai2+bi2+c22aibi+2c(aibi)=i=1n(ai2+bi2)+2ci=1n(aibi)+nc22i=1naibi 前三项都是关于一个 c 的二次函数,那么前三项通过 c 最小化即可,现在我们需要做的,便是最大化i=1naibi

我们先倍长 b,断环为链:

i=1naibi+k

因为此处 a,b 之差不变且为 k,所以我们使用减法卷积

a 数组逆序,得到:

Ank+1=i=1nani+1bi+k

然后直接 NTT 即可


三、多项式基本运算

1. 前置知识

本节只列出定义和常见公式,详细请自行查询相关资料学习

① 求导

(u(x)±v(x))=u(x)±v(x)

(u(x)v(x))=u(x)v(x)+u(x)v(x)

(u(x)v(x))=u(x)v(x)u(x)v(x)v2(x)

(u(v(x)))=v(x)u(v(x))

(xa)=axa1,(lnx)=1x,(ex)=ex

② 泰勒展开

对于一个函数 f(x) 以及一个点 x0,我们在 x0 处对函数 f 进行一个拟合,设拟合函数为 T,那么泰勒展开的一般形式如下:

T(x)=f(x0)+f(x0)1!(xx0)+f(x0)2!(xx0)2++fn(x0)n!(xx0)n=i=0fi(x0)i!(xx0)i

2. 多项式牛顿迭代

牛顿迭代解决了以下 problem:

给定多项式 g(x),已知有 f(x)

g(f(x))0  (modxn)

求出模 xn 意义下的 f(x)

solution 考虑倍增

首先当 n=1 是时 [x0]g(f(x))=0 的解需要单独求出。

假设现在已经得到了模 xn2 意义下的解 f0(x),要求模 xn 意义下的解 f(x)

g(f(x))f0(x) 处进行泰勒展开,有:

i=0+g(i)(f0(x))i!(f(x)f0(x))i0  (modxn)

因为 f(x)f0(x) 的最低非零项次数最低为 n2,所以:

i2,(f(x)f0(x))i0   (modxn)

那么:

i=0+g(i)(f0(x))i!(f(x)f0(x))ig(f0(x))(f(x)f0(x))g(f0(x))+g(f0(x))(f(x)f0(x))0  (modxn)

f(x)f0(x)g(f0(x))g(f0(x))  (modxn)

3. 多项式求逆

设给定函数为 h(x),求 f(x) 使得 f(x)h(x)1  (modxn)

solution 由上式可得: g(f(x))=1f(x)h(x)0  (modxn) 使用牛顿迭代得到: f(x)f0(x)1f0(x)h(x)1f02(x)(modxn)2f0(x)f02(x)h(x)(modxn) 时间复杂度: T(n)=T(n2)+O(nlogn)=O(nlogn)

4. 多项式求ln、exp

① 定义

我们先了解 ef(x)/exp(f(x)) 在多项式中的含义是什么

我们先泰勒展开,可以得到:

ex=i=0+xii!

exp(f(x))=ef(x)=i=0+f(x)ii!

exp 代表任意多个相应组合对象的组合。(f(x)i 即为将 i 个对象组合;i! 消除了各组 之间的差异,避免重复计数)

ln 便是 exp 的逆运算。

② 多项式求ln

设给定函数为 h(x),求 f(x) 使得 f(x)=ln h(x)  (modxn)

f(x)ln(h(x))  (modxn)

两边同时求导,得:

f(x)ln(h(x))h(x)  (modxn)

f(x)h(x)h(x)

最后使用多项式求逆,并且积分积回去即可:

f(x)h(x)h(x)

③ 多项式exp

设给定函数为 h(x),求 f(x) 使得 f(x)eh(x)  (modxn)

由上式可得:

g(f(x))=lnf(x)h(x)  (modxn)

使用牛顿迭代得:

f(x)f0(x)lnf0(x)h(x)1fo(x)(modxn)f0(x)(1lnf0(x)+h(x))(modxn)

时间复杂度:

T(n)=T(n2)+O(nlogn)=O(nlogn)

5. 多项式开方

设给定函数为 h(x),求 f(x) 使得 f2(x)h(x)  (modxn)

由上式可得:

g(f(x))=f2(x)h(x)0  (modxn)

使用牛顿迭代得:

f(x)f0(x)f02(x)h(x)2f0(x)(modxn)f02(x)+h(x)2f0(x)

时间复杂度:

T(n)=T(n2)+O(nlogn)=O(nlogn)

边界条件需要求一个数在模意义下的二次剩余,需使用 Cipolla 算法

6. 多项式除法&取模

给定多项式 f(x),g(x),求 g(x)f(x) 的商 Q(x) 和余数 R(x),即

f(x)=g(x)Q(x)+R(x)

发现若能消除 R(x) 的影响则可直接多项式求逆解决。

考虑构造变换:

fR(x)=xdegff(1x)

即,fR(x) 表示反转 f(x) 的系数,degf 表示多项式 f(x)

n=degf

posted @   ricky_lin  阅读(270)  评论(0编辑  收藏  举报
点击右上角即可分享
微信分享提示
评论
收藏
关注
推荐
深色
回顶
收起
  1. 1 有我 周深
有我 - 周深
00:00 / 00:00
An audio error has occurred.

作词 : 唐恬/闫光宇

作曲 : 钱雷

编曲 : 赵兆/付虹宇

制作人 : 赵兆

出品 : 共青团中央宣传部

版权 : 中国青少年新媒体协会

制作单位 : 能量悦动音乐

发行单位 : 银河方舟StarNation

出品人 : 郭峰

总监制 : 汤杰

总策划 : 钟亚楠

总统筹 : 金慧子

音乐监制 : 李天鹏/李三木

制作执行 : 张不贰/高聪怡

项目宣发 : 肖健/张国党/孙小千/戴胤/孙雯璟

音乐推广 : 代诗琪/杜思潮/马越/程铁峰/傅之豪

钢琴 : 赵兆

吉他 : 伍凌枫

贝斯 : 韩阳

鼓 : 武勇恒

合唱设计 : 赵兆

合唱 : 凡尔赛合唱团

人声录音 : 耿潇微

人声录音室 : 55TEC Studio Beijing

配唱 : 徐威@52Hz Studio (Shanghai)

混音 : 李游(小骷髅)@55TEC Studio Beijing

海报 : 格子

特别鸣谢 : 周深工作室

世界问 你是谁 来自哪 请回答

爱什么 梦什么 去何方 请回答

答案有 一百年的时光

我来自 硝烟中 课桌旁 的太阳

我来自 硝烟中 课桌旁 的太阳

他和她 宣的誓 迎的仗

来自那 燃烧的 和我一样 的年华

来自世间 一对平凡的夫妻 身旁

来自世间 一对平凡的夫妻 身旁

来自昨天 谁以青春赴万丈 理想

我是寸土 不让的 家乡啊

我是绝不 低头的 倔强啊

接过万千热血 的初衷

当有对答世界 的音量

要怎么形容明天 像我一样

要怎么形容明天 像我一样

承风骨亦有锋芒 有梦则刚

去何方 去最高 的想象

前往皓月星辰 初心不忘

那未来如何登场 有我担当

那未来如何登场 有我担当

定是你只能叫好 那种辉光

护身旁 战远方 有我啊

我的名字就是 站立的地方

Wu~

我的样子 就是 明天的模样

我是朝阳 落在乡间听书声 朗朗

我是朝阳 落在乡间听书声 朗朗

我是屏障 为谁挡一程厄运 的墙

我要一生 清澈地 爱着啊

我要长歌 领着风 踏着浪

朝着星辰大海 的方向

当有对答世界 的音量

要怎么形容明天 像我一样

要怎么形容明天 像我一样

承风骨亦有锋芒 有梦则刚

去远方 去最高 的想象

前往皓月星辰 初心不忘

那未来如何登场 有我担当

那未来如何登场 有我担当

定是你只能叫好 那种辉光

护身旁 战远方 有我啊

一生骄傲为我 站立的地方

Wu~

我的样子 就是 中国的模样

Wu~~~ Wu~~~

当炬火 去化作那道光

“谨以此歌献给一代代不负时代重托的中国青年”

“谨以此歌献给一代代不负时代重托的中国青年”