多项式大杂烩

· · 个人记录

前言

这个坑咕了好久了,一直没敢填。真实原因是我数学太菜

趁着新冠病毒爆发必须宅在家里这几天赶紧补一补qwq。

为了防止文章爆炸,同时方便大家阅读,我把代码都扔到云剪切板了,文章里只丢个链接。

文章里面除了多项式全家桶,还会补充一些需要用到的知识,比如'二次剩余'。

方便起见,本文里的 f[i] 表示 [x^i]f(x),即 f(x)i 次项系数。

逐渐更新的vector版多项式板子

FFT

任意模数版本

目录

1. 拉格朗日插值

1.1 公式

1.2 应用

1.3 code

2. 快速傅里叶变换 (FFT)

2.1 概况

2.2 单位根

2.3 DFT

2.4 IDFT

2.5 蝴蝶变换优化

2.6 code

3. 快速数论变换 (NTT)

3.1 原根

3.2 NTT与原根应用

3.3 code

4. 多项式求逆

4.1 推式子

4.2 code

5. 多项式(带余)除法

5.1 推式子

5.2 code

6. 多项式对数函数 (多项式Ln)

6.1 简单微积分

6.2 推式子

6.3 code

7.多项式指数函数 (多项式exp)

7.1 泰勒展开(泰勒多项式,泰勒公式)

7.2 牛顿迭代

7.3 推式子

7.4 code

8.多项式开根

8.1 二次剩余

8.2 推式子

8.3 code

9.多项式快速幂

9.1 推式子

9.2 code

10.多项式多点求值

10.1 推式子

10.2 code

11.多项式快速插值

11.1 洛必达法则

11.2 推式子

11.3 code

12.多项式三角函数

12.1 推式子

12.2 code

13.多项式反三角函数

13.1 推式子

13.2 code

14.MTT (FFT进阶)

14.1 FFT的精度优化

14.2 code

114514.TIPS

正文

一个 n多项式是形如这样的一个东西:

f(x)=\sum\limits_{i=0}^n{a_i x^i}\ \ \ (a_n\neq 0)

其中,n 是多项式的最高次幂,a_{i} 为第 i 项系数。

1.拉格朗日插值

(题目)

问题:给定 n\ (n\leq 2\times 10^3) 个点 (x_i,y_i),求经过这 n 的点的 n-1 次多项式。(或求在 k 这个点的取值)

首先珂以暴力高斯消元,不过这样是 \mathcal{O}(n^3),显然布星。

1.1 公式

背住这个公式就好了,证明..太复杂(珂以去问Karry神仙,至少你要知道什么是范德蒙德矩阵的行列式qwq)。

f(k)=\sum\limits_{i=1}^n{y_i \prod\limits_{j\neq i}{\frac{k-x_j}{x_i-x_j}}}

这样在 \mathcal{O}(n^2) 的时间复杂度内就可以解决掉了!

如果 x_{i}=i 的之后有更优秀的解决方法:

t=\prod\limits_{i=1}^n{(k-i)}

那么式子就变成了:

f(k)=\sum\limits_{i=1}^n{\left( {\frac{y_i\cdot t}{(i-1)!\ (n-i)!\ (k-i)}\cdot(-1)^{n-i}}\right)}

通过求前缀积和后缀积,这样 \mathcal{O}(n) 或者 \mathcal{O}(n\log n) 就可以解决。

1.2 应用:

求:\sum\limits_{i=1}^{n}{i^k} \quad (n\leq 10^{13},k\leq 10^5)

暴力 \mathcal{O}(n\log k) 布星,需要更优秀的方法。

Karry巨神说:珂以用伯努利数吖!!

↑这个毒瘤东西我们暂且不说↑

老祖宗告诉我们,它是一个关于 nk+1 次的多项式。(原理其实是将式子给差分 k+1 次)

举个大家幼儿园都会的例子:

k=1 的时候,f(n)=\frac{n(n+1)}{2},是个2次的;

k=2 的时候,f(n)=\frac{n(n+1)(2n+1)}{6},是个3次的。

很神奇吧!带入 x_i=1\dots k+2 求出对应的 y_{i},再套改采讲的的式子就 \mathcal{O}(k\log k) 解决掉啦!

1.3 code here (是那道模板题的)

2.快速傅里叶变换 (FFT)

(题目)

问题:给定两个多项式 f(x)g(x), 求 f(x)*g(x),即两个多项式的卷积。

什么是卷积?就是求一个多项式 p(x)=f(x)*g(x) 满足 p[k]=\sum\limits_{i+j=k}f[i]\times g[j]

说白了就是小学一年级学的多项式乘法嘛!

显然有个 \mathcal{O}(n^2) 的暴力求法,不过FFT在 \mathcal{O}(n\log n) 的时间复杂度内求出来。

2.1 概况

FFT的主要流程就是系数表达式->点值表达式->系数表达式。(第一个过程叫DFT,第二个过程叫IDFT)

点值表达式就是用 n+1 个点 (x_i,y_i) 来表示一个 n 次多项式。

为什么要转点值表达式?因为把两个点值表达式乘一起 \mathcal{O}(n) 就搞定了!

把对应的 x_iy_i 乘起来就变成了要求的多项式的点值表达式!

DFT和IDFT肿么搞?

2.2 单位根

直接带入 0,1,2,\dots n 显然不可取,这样带入复杂度就是 \mathcal{O}(n^2) 的了。

先来讲一讲单位根。首先你要知道什么是虚数,还要会它的加减乘除运算,不会的又转百度。

先问个问题:x^n=1 在复数范围内的 n 个解分别是什么?

答案:就是把复平面内的单位圆从 (1,0) 开始 n 等分,圆周上 n 个点的坐标就是那 n 个解。设这 n 个解分别是 \omega_n^0,\omega_n^1,\cdots ,\omega_n^{n-1}

显然,\omega_n^0=1\omega_n^1 可以用三角函数搞一下,就是 \cos\left( \frac{2\pi}{n}\right)+\sin\left( \frac{2\pi}{n}\right)i这个也很显然, 其实就是三角函数与单位圆的关系,不懂三角函数回小学补一补

他们有一些神奇的性质:

  1. \omega_n^k=\omega_{pn}^{pk}$ $\quad (p\in \mathbf{N^*})
  2. \omega_n^k=(\omega_{n}^{1})^k

2.3 DFT

下面开始重点!FFT实际上是一个分治的过程

我们把 f(x) 拆成两部分,分奇偶打乱。

即:f(x)=(f[0]+f[2]x^2+\dots +f[n-2]x^{n-2})+(f[1]x+f[3]x^3+\dots +f[n-1]x^{n-1})

假设我们已经通过分治求出了:

fl(x)=f[0]+f[2]x+\dots +f[n-2]x^{n/2-1} fr(x)=f[1]+f[3]x+\dots +f[n-1]x^{n/2-1}

那么:f(x)=fl(x^2)+x\cdot fr(x^2)

由于 \omega 好的性质,这样分治就珂以搞啦!!

结论就是:(式子根据那些性质很好推,我就不推了)

f(\omega_n^k)=fl(\omega_{n/2}^k)+\omega_n^k\cdot fr(\omega_{n/2}^k) f(\omega_n^{k+n/2})=fl(\omega_{n/2}^k)-\omega_n^k\cdot fr(\omega_{n/2}^k)

(这个很好记吧,第一个式子就是直接带入的结果,第二个和第一个只有一个正负号的区别)

这样就得到了一坨点值表达式啦!

2.4 IDFT

我们得到了一坨点值表达式:s[k]=\sum\limits_{i=0}^{n-1}{f[i]\cdot (\omega_n^k)^i}

那么最后要求的多项式长这样 n \times p[k]=\sum\limits_{i=0}^{n-1}{s[i]\cdot (\omega_n^{-k})^i}

注: \omega_n^k=\omega_n^{k\mod n},所以 \omega_n^{-k}=\omega_n^{n-k}

证明的话呢带入就好啦,我就不啰嗦了。

然后和DFT的流程一摸一样了,唯一的区别就是带入了逆单位根

码量就可以大大缩小了,只要传一个参数表示现在做的是DFT还是IDFT来决定带入的单位根就行。

2.5 蝴蝶变换优化

通过前面的内容已经珂以写成递归版的FFT啦!

但是常数极大,连模板题都过不掉,所以要写非递归版的。

BUT,要分奇偶打乱,这该肿么搞?直接把每一项的系数放到它最后应该在的位置上就好了。

怎么放?先看一张图。

最后一层是系数分奇偶打乱后最终的位置。(位置编号从0开始)

仔细研究发现,奇偶打乱后每个数的位置是其原数经过二进制反转后的位置

比如 6=(110)_2,打乱后的位置就是 (011)_2=3

这个通过 \mathcal{O}(n) 递推珂以得到。

tr_ii 经过奇偶打乱后的位置。

递推式:

tr[i]=((tr[i>>1]>>1)|((i&1)?n>>1:0));

希望大家自行理解。 提示:把 i 拆成 最低位和剩余部分,两部分分别反转,最低位反转是 ((i&1)?n>>1:0),剩余部分反转是 (tr[i>>1]>>1)

然后就珂以去实现非递归版(迭代版)的啦!!

还有一个细节,n 不是 2^k 不好搞,那就先把 n 增加到 2^k,多出的位置补0就行了。(同时别忘了把空间开大!)

2.6 code here

3.快速数论变换 (NTT)

FFT需要算一堆三角函数太慢,而且系数比较大的时候精度也会成问题导致WA。这时候NTT就出场了。

3.1 原根

定义:

g 是奇质数 p 的原根,则 g^1,g^2,\dots,g^{p-1} 在模 p 意义下互不相同。

判断一个数 x 是否是奇质数 p 的原根的方法:

先把 p-1 质因数分解,p-1=p_1^{q_1}\times p_2^{q_2}\times \dots\times p_k^{q_k}

若对于 \forall{i}\in[1,k],满足:x^{(p-1)/{p_i}}\not\equiv1\pmod{p},则 xp 的原根。

常用原根:

我太懒了,直接去这里康吧。。

一定要记住 998244353 的最小原根是 3 就行!

3.2 NTT与原根应用

知道了原根有神马用?

我们把 g^{(p-1)/t} 当作 \omega_{n}^1 (\ t 是形如 2^k 的一个数,奇质数 p 要满足形如 r\cdot 2^k+1\ ,比如 998244353=119\times 2^{23}+1\ )

惊奇地发现:w 满足的优秀性质原根都满足!!

这样就好啦!!把FFT原本的单位根换成原根,再也不会有精度问题了!! 其他部分就和FFT一毛一样了。

但是要注意,NTT只能处理 n\leq 2^k 的数据,当模数比较毒瘤,比如 10^9+7 ,此时 k=1,NTT不再能用。(不过不用怕,此后会介绍任意模数NTT,但更加复杂。)

3.3 code here

4.多项式求逆

(题目)

问题:给定 f(x),求 g(x) 满足:f(x)*g(x)\equiv1\pmod{x^n}

求逆的主要思想是倍增。(以后好多题目都要用到倍增思想)

4.1 推式子:

假设已经求出了一个 h(x),满足 h(x)\equiv f(x)^{-1}\pmod{x^{\frac{n}{2}}}

由于 g(x)\equiv f(x)^{-1}\pmod{x^n}

通过做差得到:g(x)-h(x)\equiv0\pmod{x^{\frac{n}{2}}}

然后在式子两边同时平方(为了把模数的次数凑成 x^n ):

g(x)^2-2*g(x)*h(x)+h(x)^2\equiv0\pmod{x^n}

然后机智得把两边同时乘上个 f(x),由于 f(x)*g(x)\equiv1\pmod{x^n}

所以珂以得到: g(x)-2*h(x)+f(x)*h(x)^2\equiv0\pmod{x^n}

即:

g(x)\equiv 2*h(x)-f(x)*h(x)^2\pmod{x^n}

套个NTT板子就把 g(x) 求出来啦!!

通过这个倍增就好了。

时间复杂度:

T(n)&=T(\frac{n}{2})+\mathcal{O}(n\log n) \\ &=\mathcal{O}(n\log n) \end{aligned}

不要误以为是两个 \log 了! (不明白?想一想线段树建树的复杂度)

4.2 code here

5.多项式(带余)除法

(题目)

问题:给定 f(x),g(x),求 q(x),r(x) 满足:f(x)=q(x)*g(x)+r(x)deg_q=deg_f-deg_g,deg_r<deg_g

注: deg_f 为多项式 f(x) 的次数

5.1 推式子:

做法玄学

先带入 \frac{1}{x} 得:f(\frac{1}{x})=q(\frac{1}{x})*g(\frac{1}{x})+r(\frac{1}{x})

两边同乘 x^n 得:x^nf(\frac{1}{x})=x^{n-m}q(\frac{1}{x})*x^mg(\frac{1}{x})+x^nr(\frac{1}{x})

定义 f^R(x)=x^{deg_f}f(\frac{1}{x}) (等同于把 f(x) 的系数反转后再带入 x )

则:f^R(x)=q^R(x)*g^R(x)+r^R(x)\cdot x^{n-m+1}

同时模 x^{n-m+1} 珂消去 r^R(x)\cdot x^{n-m+1},得:

f^R(x)\equiv q^R(x)*g^R(x)\pmod{x^{n-m+1}}

g^R(x) 移一下得:

q^R(x)\equiv \frac{f^R(x)}{g^{R}(x)}\pmod{x^{n-m+1}}

右边这坨把 g^R(x) 求个逆然后再乘上 f^R(x) 就行了!

由于 deg_{q^R}<n-m+1 ,所以 q^R(x)\mod x^{n-m+1}=q^R(x)

这样 q^R(x) 就珂以被求出来啦!只用把 q^R(x) 的系数反转一下就是 q(x),直接用STL库里的reverse函数就行了。

$$r(x)=f(x)-q(x)*g(x)$$ 很神奇吧! 时间复杂度和求逆是一样的:$\mathcal{O}(n\log n)

5.2 code here

代码写起来非常恶心。。要调对着标程调半年qwq。

6.多项式对数函数 (多项式Ln)

(题目)

问题:给定 f(x),求 g(x) 满足:g(x)\equiv \ln f(x)\pmod{x^n}

想一想 \ln 最初是怎么定义的? 对 \ln x 求导的结果是 \frac{1}{x}。也就是 \frac{\mathrm{d}\ln x}{\mathrm{d}x}=\frac{1}{x}

不知道微积分??

6.1 简单微积分

微分

对于一个函数 f(x),如果它在 x_0 处可导,则 f(x)x_0 处的导数为:

f'(x_0)=\lim\limits_{\Delta x\to 0}{\frac{f(x_0+\Delta x)-f(x_0)}{\Delta x}}

如果你不知道 \lim 是什么,你可以把 \lim\limits_{\Delta x\to 0} 理解为:\Delta x 无限接近 0 但永远也不等于 0

如果函数在连续区间上可导,则函数在这个区间上存在导函数,记作 f'(x)\frac{\mathrm{d}f(x)}{\mathrm{d}x}

不知道极限?那真没救了。

如果 f(x) 是多项式函数,即:

f(x)=\sum\limits_{i=0}^n a_i x^i

则:

f'(x)=\sum\limits_{i=1}^{n}i a_i x^{i-1} #### 积分 积分是微分的逆运算。也就是对一个函数先求个微分然后再积分起来还是原来的函数。 求法就是去求一段函数图像的代数面积。(几何意义) $\int f(x)\mathrm{d}x$ 表示对函数 $f(x)$ 求积分。 同样,对多项式函数积分也是有公式滴! $$\int f(x)\mathrm{d}x=\sum\limits_{i=1}^{n+1}{\frac{a_{i-1} x^{i}}{i}}$$ 由于常数项不确定,它还是个**不定积分**。 $\mathcal{O}(n)$ (需要预处理逆元)或 $\mathcal{O}(n\log n)$ 均可。 #### 链式法则 对一个复合函数该怎么求导? 其实很简单,“拆”一下,分别求个导就好了! $$f'(g(x))=\frac{\mathrm{d}f(g(x))}{\mathrm{d}x}=\frac{\mathrm{d}f(g(x))}{\mathrm{d}g(x)}\cdot \frac{\mathrm{d}g(x)}{\mathrm{d}x}$$ 不需要我解释吧,约分一下就好了。。很显然。 ### 6.2 推式子 知道了简单微积分,解决这道题就小菜一碟啦! 原式:$g(x)\equiv \ln f(x)\pmod{x^n}

两边同时求导(用链式法则):g'(x)\equiv (\ln' f(x))\cdot f'(x)\pmod{x^n}

根据 \frac{\mathrm{d}\ln x}{\mathrm{d}x}=\frac{1}{x},式子就变成了:

g'(x)\equiv \frac{f'(x)}{f(x)}\pmod{x^n} g(x)\equiv \int\frac{f'(x)}{f(x)}\pmod{x^n}

然后求个逆,求个导,求个乘法,求个积分,开始开心地套板子啦^_^

时间复杂度仍然是:\mathcal{O}(n\log n)

6.3 code here

这道题的代码还是非常好写的23333...

7.多项式指数函数 (多项式exp)

(题目)

问题:给定 f(x),求 g(x) 满足:g(x)\equiv \mathrm{e}^{ f(x)}\pmod{x^n}

多项式$\exp$ 需要更多~~恶臭的~~前置知识qwq。 ### 7.1 泰勒展开(泰勒多项式,泰勒公式) 先给出泰勒展开的~~恶心的~~式子: $$f(x)=f(x_0)+\sum\limits_{i=1}^{n}{\frac{f^{(i)}(x_0)}{i!}(x-x_0)^i+R_n(x)}$$ 这个式子需要保证 $f(x)$ 在 $x=x_0$ 的邻域内 $n+1$ 阶可导。其中,$f^{(i)}(x)$ 表示 $f(x)$ 的 $n$ 阶导数,$R_n(x)$ 是佩亚诺(Peano)余项,就是 $(x-x_0)^n$ 的高阶无穷小。(知道就好啦,不用太深究) 珂以形象的理解为通过求导将构造一个函数是它与 $f(x)$ 无限地接近。 这个玩意和生成函数有密切联系,之后再说。 麦克劳林公式就是泰勒公式当 $x_0=0$ 的时候的特殊情况。同样需要保证$f(x)$ 在 $x=0$ 的邻域内 $n+1$ 阶可导。 $$f(x)=f(0)+\sum\limits_{i=1}^{n}{\frac{f^{(i)}(0)}{i!}x^i+R_n(x)}$$ ### 7.2 牛顿迭代 问题:已知 $f(x)$,求 $g(x)$,满足:$f(g(x))\equiv0\pmod{x^n}

牛顿迭代的主要思想也是倍增。同时要用到上面的泰勒展开,先记住喽!

沿用多项式求逆的思路,假设我们已经求出来了 h(x) 满足:f(h(x))\equiv0\pmod{x^{\frac{n}{2}}}

然后把 f(g(x))h(x) 处泰勒展开, (smg),得到如下的一坨式子:

f(g(x))=f(h(x))+f'(h(x))(g(x)-h(x))+\frac{f''(h(x))}{2!}(g(x)-h(x))^2+\cdots

由于 g(x)h(x) 的前 \frac{n}{2} 项是相同的,所以: g(x)-h(x)\equiv0\pmod{x^{\frac{n}{2}}}

那么在刚才泰勒展开的那一坨式子中,在 \mod{x^n} 的意义下,由于 (g(x)-h(x))^i \ \ (i\ge2) 的存在,从从第三项开始的值都为0。

原本的式子就被简化成了:

f(g(x))\equiv f(h(x))+f'(h(x))(g(x)-h(x))\pmod{x^n}

由于 f(g(x))\bmod x^n=0,通过移项就能把 g(x) 求出来,得到了下面的这个最终式:

g(x)\equiv h(x)-\frac{f(h(x))}{f'(h(x))}\pmod{x^n}

像求逆一样迭代一下就好啦!

7.3 推式子

最初的式子:g(x)\equiv \mathrm{e}^{ f(x)}\pmod{x^n}

就按常规操作把原来的式子两边同时取对数:\ln g(x)\equiv f(x)\pmod{x^n}

令: P(g(x))=\ln g(x)-f(x)

那么 P(g(x))\equiv0\pmod{x^n}

这就可以套刚才讲的牛顿迭代啦!

假设已经求出来了 h(x) 满足:P(h(x))\equiv0\pmod{x^{\frac{n}{2}}}

g(x)\equiv h(x)-\frac{P(h(x))}{P'(h(x))}\pmod{x^n}

$P'(h(x))=\frac{1}{h(x)}

分别把 h(x) 带到式子里面就得到了:g(x)\equiv h(x)-\frac{\ln h(x)-f(x)}{\frac{1}{h(x)}}\pmod{x^n}

化简一下:

g(x)\equiv h(x)\cdot(1-\ln h(x)+f(x))\pmod{x^n}

然后就是套多项式 \ln 的板子了

时间复杂度:\mathcal{O}(n\log n)

7.4 code here

8.多项式开根

(题目) && (题目-加强版)

问题:给定 f(x),求 g(x) 满足 g^2(x)\equiv f(x)\pmod{x^n}

注: 加强版没有保证 f[0]=1,此时需要二次剩余。

8.1 二次剩余

(题目)

问题:给定 x,p,保证 p奇质数,解方程:x^2\equiv n\pmod{p}

先说一说什么是二次剩余。

概念:

对于任意的一个自然数 n ,若对于一个质数 p,满足关于 x 的方程 x^2\equiv n\pmod{p} 有整数解,则称 n 是模 p 的二次剩余。

判定方法:

先引入一个勒让德记号

1&(\,p\nmid n\ \text{且}\,n\,\text{是模}\,p\,\text{的二次剩余})\\ -1&(\,p\nmid n\ \text{且}\,n\,\text{是模}\,p\,\text{的二次非剩余})\\ 0&(\,p\mid n) \end{cases}

根据勒让德记号又有个欧拉判别准则

\left(\frac{n}{p}\right)\equiv n^{(p-1)/2}\pmod{p}

通过费马小定理珂以证明,这里不详细说了qwq。

现在就可以判断是否有解啦!

求解:

如果有解肿么办?

下面讲 Cipolla算法

先随机找到一个自然数 a\in [0,p),令 \omega^2=a^2-n,使得 \left(\frac{\omega^2}{p}\right)=-1

其实期望随机2次就出来了。

然后其中一个解就是: (a+\omega)^{(p+1)/2} (\,另一个解是 p-(a+\omega)^{(p+1)/2}\,)

这个玩意等于 n ? 不可思议,下面证明一下:

别急 首先有个结论: (a+b)^p\equiv a^p+b^p\pmod{p} 证明很简单,拿二项式定理展开一下就好了,不啰嗦啦。

x\equiv (a+\omega)^{(p+1)/2}\pmod{p}

则: x^2\equiv(a+\omega)^{p+1}\pmod{p}

根据上面的结论得到:x^2\equiv(a^p+\omega^p)(a+\omega)\pmod{p}

由于:\left(\frac{\omega^2}{p}\right)\equiv\omega^{p-1}\equiv -1\pmod{p},还有费马小定理:a^{p-1}\equiv1\pmod{p}

所以:x^2\equiv(a^p-\omega)(a+\omega)\equiv(a-\omega)(a+\omega)\equiv a^2-w^2\pmod{p}

根据前面的定义:\omega^2=a^2-n 就得到了:

x^2\equiv n \pmod{p}

证毕!!

现在唯一的问题就是如何求:(a+\omega)^{(p+1)/2}

由于 \omega^2 的值是已知的,我们珂以把 \omega 看成虚部,两个形如 a+b\omega 相乘会得到:

(a_1+b_1\omega)\cdot(a_2+b_2\omega)=(a_1b_1+a_2b_2)\omega^2+(a_1b_2+a_2b_1)\omega

重载一波运算符即可,code here

8.2 推式子

(需要用到7.2讲的牛顿迭代,一定要先把那个式子记好喽!)

g(x)^2\equiv f(x)\pmod{x^n}

沿用多项式exp的思路,我们设:P(g(x))=g^2(x)-f(x)

则:P(g(x))\equiv0\pmod{x^n}

然后套牛顿迭代:假设已经求出来了 h(x) 满足:P(h(x))\equiv0\pmod{x^{\frac{n}{2}}}

g(x)\equiv h(x)-\frac{P(h(x))}{P'(h(x))}\pmod{x^n}

和exp一样,在这里把 f(x) 看成是常数,那么套多项式求导的公式得到:P'(h(x))=2h(x)

h(x) 带到式子里面得:g(x)\equiv h(x)-\frac{h(x)^2-f(x)}{2h(x)}\pmod{x^n}

化简之后就得到了最终的结果:

g(x)\equiv \frac{h(x)^2+f(x)}{2h(x)}\pmod{x^n}

套板子好啦。但还有一个边界就是要求 g[0],它满足:

g[0]^2\equiv f[0]\pmod{p}

对于非加强版满足 f[0]=1,显然 g[0]=1

对于加强版直接套刚才讲的二次剩余即可。

8.3 code here && here(加强版)

常数巨大!不开O2最慢点要跑1s+0.05s左右qwq。

其实还有个exp+ln的做法很简便,不过代码难写+常数更大,不建议。

思考题: 如果是开 k 次根该怎么搞?提示:还是要用牛顿迭代。

9.多项式快速幂

(题目) & (加强版)

问题:给定 f(x),k\ (k\leq 10^{10^5}),求 g(x) 满足:g(x)\equiv f^k(x)\pmod{x^n}

式子很好推,但要看你的码力+卡常神功了23333。。。

9.1 推式子:

原式:g(x)\equiv f^k(x)\pmod{x^n}

两边同时取对数 ( ln 一下):\ln g(x)\equiv k\ln f(x)\pmod{x^n}

两边同时取指数 ( exp 一下):g(x)\equiv \mathrm{e}^{k\ln f(x)}\pmod{x^n}

然后就是套板子,注意虽然 k 很大,但是它是系数,珂以边读边取余。

珂能大家会有个疑问:

Q: 根据欧拉定理/费马小定理,k 不应该对 998244352 (即 p-1 ) 取模么?为什么要对 998244353 取模?

A: 这是在模 998244353 的意义下对 k\ln f(x) 做exp,而 k 是在系数部分的,所以显然 k 要对 998244353 取模。

9.2 code here

为毛我这道题没卡常就跑的非常快^_^

10.多项式多点求值

(题目)

问题:给定 f(x)m 个整数 a_i ,分别求出 f(a_i)

秦九韶的暴力 \mathcal{O}(nm) 理论上过不了 (然鹅循环展开卡卡常然后夜深人静的时候提交竟然可以过,这时候时间复杂度是 \mathcal{O}(wyh) ??

本题要用分治的思想。(终于轮到分治上场了,之前都是倍增)

10.1 推式子

我们设两个多项式 g_0(x)=\prod\limits_{i=1}^{\lfloor\frac{m}{2}\rfloor}{(x-a_i)}g_1(x)=\prod\limits_{i=\lfloor\frac{m}{2}\rfloor+1}^{m}{(x-a_i)}

方便起见,我们假设 m 是偶数。

显然,g_0(x),g_1(x) 都是 \lfloor\frac{m}{2}\rfloor 次的多项式。

并且,g_0(a_1)=g_0(a_2)=\cdots=g_0(a_{\frac{m}{2}})=g_1(a_{\frac{m}{2}+1})=\cdots=g_1(a_m)=0

我们令 r_0(x)=f(x)\mod g_0(x)r_1(x)=f(x)\mod g_1(x)

f(x)=q_0(x)*g_0(x)+r_0(x)=q_1(x)*g_1(x)+r_1(x)

当把 a_{1 \dots \frac{m}{2}} 带入左边的式子中,由于 g_0(a_i)=0

所以 q_0(x)*g_0(x)=0,则 f(x)=r_0(x)

此时 r_0(x)\lfloor\frac{m}{2}\rfloor 次的,我们珂以把它当作新的 f(x),类似于上面的过程继续分治下去即可。

右边的同理。

最初珂以先分治NTT求出所有的 g(x),这时候用指针存 g(x) 就非常方便了 (当然也可以用vector)。

求所有的 r(x) 直接套多项式除法的板子即可。

时间复杂度:

T(n)&=2\cdot T(\frac{n}{2})+\mathcal{O}(n\log n) \\ &=\mathcal{O}(n\log^2 n) \end{aligned}

这回真的是两只 \log 了!

10.2 code here

这回的代码跑得不慢^_^(但开O2反而更慢了是smg?)

11.多项式快速插值

(题目)

问题:给定 n\ (n\leq 10^5) 个点 (x_i,y_i),求经过这 n 的点的 n-1 次多项式。

原题?BUT,这回拉格朗日插值的 \mathcal{O}(n^2) 复杂度显然行不通。

这时候,我们大力分治!不过先要知道洛必达法则

11.1 洛必达法则

若函数 f(x)g(x) 满足以下条件:

则:

\lim\limits_{x\to x_0}{\frac{f(x)}{g(x)}}=\lim\limits_{x\to x_0}{\frac{f'(x)}{g'(x)}}

11.2 推式子:

还记得开头拉格朗日插值么,实际上快插就是对拉格朗日插值的优化,忘了回去复习一下!

为了不使式子看起来这么恶心,当 j 没写上下界约定的时候默认 j 的取值范围和它前面的 i 的取值范围相同。( 此时 \prod 下面只保留 j\neq i\, )

原式:

f(x)=\sum\limits_{i=1}^n{y_i \prod\limits_{j\neq i}{\frac{x-x_j}{x_i-x_j}}}

先把这一坨拆成两部分:

f(x)=\sum\limits_{i=1}^n\left({\frac{y_i}{\prod\limits_{j\neq i}{(x_i-x_j)}}\cdot\prod\limits_{j\neq i}{(x-x_j)}}\right)

我们可以把 \frac{y_i}{\prod\limits_{j\neq i}{(x_i-x_j)}} 当作常数;

v_i=\frac{y_i}{\prod\limits_{j\neq i}{(x_i-x_j)}}

g(x)=\prod\limits_{j=1}^n{(x-x_j)} ,那么 \prod\limits_{j\neq i}{(x_i-x_j)}=\frac{g(x)}{x-x_j}(x_i)

**洛必达!!** 所以 $\frac{g(x)}{x-x_j}=\frac{g'(x)}{(x-x_j)'}=g'(x)

那么 v_i=\frac{y_i}{g'(x_i)}

这个分治NTT求出 g'(x),然后直接套多点求值的板子即可求出 v_i

同样可以把开始的那个大式子给分治NTT了。

f_{l,r}x_{l\cdots r} 这一段插出来的多项式。

g_{l,r}\prod\limits_{i=l}^{r}(x-x_i) (实际上 g_{l,r} 已经在求 g'(x) 的分治NTT的过程中求出)

mid=\frac{l+r}{2}

则:

&=\sum\limits_{i=l}^{mid}\left({v_i\cdot\prod\limits_{j=l,j\neq i}^{r}{(x-x_j)}}\right)+\sum\limits_{i=mid+1}^{r}\left({v_i\cdot\prod\limits_{j=l,j\neq i}^{r}{(x-x_j)}}\right) \\ &=\left[\prod\limits_{j=mid+1}^{r}{(x-x_j)}\right]\cdot\sum\limits_{i=l}^{mid}\left({v_i\cdot\prod\limits_{j=l,j\neq i}^{mid}{(x-x_j)}}\right)+ \left[\prod\limits_{j=1}^{mid}{(x-x_j)}\right]\cdot\sum\limits_{i=mid+1}^{r}\left({v_i\cdot\prod\limits_{j=mid+1,j\neq i}^{r}{(x-x_j)}}\right) \\ &=g_{mid+1,r}\cdot f_{l,mid}+g_{l,mid}\cdot f_{mid+1,r} \end{aligned}

直接递归求 f_{1,n} 即可。最后 f(x)=f_{1,n}

时间复杂度:\mathcal{O}(n\log^2 n)

11.3 code here

这是一道卡常好题!还需预处理原根次幂才能过。(具体方法见文末)

这回不开O2真过不了了qwq。卡常卡到我心态爆炸。。。

12.多项式三角函数

(题目)

问题:给定 f(x),求 g(x),满足:g(x)\equiv\sin f(x)\pmod{x^n}g(x)\equiv\cos f(x)\pmod{x^n}

emm..又是玄学东西..别说三角函数不知道

首先要把多项式exp写得很熟练。

12.1 推式子

相信大家都熟知一个欧拉公式

\mathrm{e}^{i\theta}=\cos\theta+i\sin\theta

神马?不知道?现在知道了就好。

不过大家至少一定都明白:(初学三角函数就学的内容吧)

\cos\theta=\cos(-\theta),\sin\theta=-\sin(-\theta)

那么把 -\theta 带入到欧拉公式得到:

\mathrm{e}^{i(-\theta)}=\cos\theta-i\sin\theta

这个式子和最初的欧拉公式相加得:

\cos\theta=\frac{\mathrm{e}^{i\theta}+\mathrm{e}^{-i\theta}}{2}

两式相减得:

\sin\theta=\frac{\mathrm{e}^{i\theta}-\mathrm{e}^{-i\theta}}{2i}

\theta 换成 f(x) 即可变成要求的式子:

g(x)\equiv\cos f(x)\equiv\frac{\mathrm{e}^{i\cdot f(x)}+\mathrm{e}^{-i\cdot f(x)}}{2}\pmod{x^n} g(x)\equiv\sin f(x)\equiv\frac{\mathrm{e}^{i\cdot f(x)}-\mathrm{e}^{-i\cdot f(x)}}{2i}\pmod{x^n}

那个虚数 i 是smg?

别急,先从 i 的定义入手:i^2=-1

所以 i^2\equiv-1\equiv998244352\pmod{998244353}

想到了什么?对!二次剩余

恰好是有解的,其中较小解为:i=86583718

然后把所有 i 换成 86583718,再套个exp的板子就好啦!

12.2 code here

无任何卡常依然飞快(雾

13.多项式反三角函数

(题目)

问题:给定 f(x),求 g(x) 满足:g(x)\equiv\arcsin f(x)\pmod{x^n}g(x)\equiv\arctan f(x)\pmod{x^n}

上一个板子用欧拉定理,但本题就用不了了qwq。不过珂以求导锕!求导大法吼!

13.1 推式子

相信大家在数学课上都学过:(雾

\frac{\mathrm{d}\arcsin x}{\mathrm{d}x}=\frac{1}{\sqrt{1-x^2}} \frac{\mathrm{d}\arctan x}{\mathrm{d}x}=\frac{1}{1+x^2}

先看反正弦函数

f(x) 带入进去得:

\frac{\mathrm{d}\arcsin f(x)}{\mathrm{d}f(x)}=\frac{1}{\sqrt{1-f^2(x)}}

再根据链式法则 (忘了回6.1复习一下):

\cdot \frac{\mathrm{d}f(x)}{\mathrm{d}x}=\frac{f'(x)}{\sqrt{1-f^2(x)}}

然后把原式 g(x)\equiv\arcsin f(x)\pmod{x^n} 两边同时求导:

g'(x)\equiv\frac{\mathrm{d}\arcsin f(x)}{\mathrm{d}x}\equiv\frac{f'(x)}{\sqrt{1-f^2(x)}}\pmod{x^n}

再把两边积分回去就得到要求的:

g(x)\equiv\arcsin f(x)\equiv\int{\frac{f'(x)}{\sqrt{1-f^2(x)}}\mathrm{d}x}\pmod{x^n}

然后疯狂套板子:微分,开根,求逆,积分!

处理反正切函数同理,当成思考题,大家照着上边的过程推一推。

这里就放个结论:

g(x)\equiv\arctan f(x)\equiv\int{\frac{f'(x)}{1+f^2(x)}\mathrm{d}x}\pmod{x^n}

13.2 code here

没开O2没卡常就飞快到卡进最优解第二页^_^

14.MTT (FFT进阶)

(题目1) && (题目2)

是我们的老熟人--多项式乘法和多项式求逆

BUT,模数不再是NTT模数,这下子普通的NTT不再能用qwq。

但是珂以选三个NTT模数,然后拿CRT(中国剩余定理)合并。

不过那样要9次NTT,常数太大而且代码不好写qwq。(wtcl,勿喷)

14.1 FFT的精度优化

由于 a_i\leq 10^9 ,所以常规的FFT的精度会爆炸!

考虑如何优化精度。

可以把系数拆开吖! 把每个系数拆成 aw+b 的形式。

通常 w2^{15}=32768

类似的,我们就把一个多项式 f(x) 拆成了两个。

f_0(x)f(x) 的‘低位’,f_1(x)f(x) 的‘高位’。

那么,对于每个多项式: f(x)=f_0(x)+f_1(x)\cdot w

此时把两个多项式相乘就得到了:

&=f_1(x)g_1(x)\cdot w^2+(f_0(x)g_1(x)+f_1(x)g_0(x))w+f_0(x)g_0(x) \end{aligned}

这样分别对 f_0(x),f_1(x),g_0(x),g_1(x) 做DFT,

然后再对 f_1(x)g_1(x),(f_0(x)g_1(x)+f_1(x)g_0(x)),f_0(x)g_0(x) 做IDFT即可

最后合并起来就得到了要求的卷积!

共计7次FFT。

其实还有4次FFT的优秀做法,不过我太菜了qwq。

14.2 code code1 && code2

这两道道题还卡精度,要开long double (这样似乎比NTT跑得更慢啦?qwq),不过MTT的代码极为好写是真的!

114514.TIPS

qwq,众所周知,多项式的常数极大,并且不是很好写,下面介绍一下我知道tips。

结语

终于完啦^_^

希望大家多向我指出细节错误和不足(比如哪里写得不是很详细),我会很感激的!

谢谢大家! (到这里源码有1000+行了qwq)