多项式入门
前言
作者会讲的尽量详细。
当然了，因为作者不太行，有些东西无法给出证明。同时，可能会有很多不严谨。
但是写的应该会对新手非常友好。（毕竟作者也是新手）
代码在文末，有需要的可以参考。
作者的板子常数非常魔鬼。建议弄一个常数小点的板子……不要对着作者的写……
前置知识
\(\text{FFT} \& \text{NTT}\)。没了。
记号：
\([x^t]f(x)\) 表示 \(f(x)\) \(t\) 次项的系数。
\(\varepsilon\) 表示一个趋近于零的极小量。（也就是无限小）\(\varepsilon>0\)，且小于任何正数。有限个 \(\varepsilon\) 的和仍然满足上述性质。
\(e=\lim_{x \to \infty}\limits (1+\dfrac{1}{x})^x,\ln(x)=log_e(x),\exp(x)=e^x\)
\(f'(x)\) 表示 \(f(x)\) 的导数。不会求导没关系，会讲的。
\(f^{(i)}\) 表示 \(f(x)\) 的 \(i\) 阶导数。\(f^{(0)}(x)=f(x)\) ，\(f^{(2)}(x)\) 就是 \(f'(x)\) 的导数，以此类推。
下文默认系数对 \(998244353\) 取模。
关于生成函数：
生成函数与多项式有着千丝万缕的联系。因此在这里写一点生成函数的基础知识。
生成函数每一项的系数存放着一些信息。（如方案数）。生成函数每一项会带上一个“形式”，类似于下标。例如：
\(\text{OGF}\) ，第 \(t\) 项会带上 \(x^t\)
\(\text{EGF}\) ，第 \(t\) 项会带上 \(\dfrac{x^t}{t!}\)
其中， \(\text{OGF}\) 常常表示 “数量”，即 \(x^t\) 表示 \(t\) 个物品；
例如描述用一个大小为 \(2\) 时，内部方案有 \(3\) 种，大小为 \(1\) 时，内部方案有 \(1\) 种的组 \(\alpha\) 和一个大小为 \(1\) 时，内部方案有 \(4\) 种，大小为 \(0\) 时，内部方案有 \(1\) 种的组 \(\beta\) 中挑选一些物品得到组 \(\gamma\) 的方案数（\(\alpha\) 的物品与 \(\beta\) 的物品之间有区别）：
我们可以给它们 \(3\) 个\(\text{OGF}\) ：\(f_\alpha(x),f_\beta(x),f_\gamma(x)\)
也就是说，选 \(3\) 个物品的方案数为 \(12\)，选 \(2\) 个物品的方案数为 \(7\) ，选 \(1\) 个物品的方案数为 \(1\)。
你也许可以感受到，我们可以用乘法来处理生成函数之间的“合并”。
\(\text{EGF}\) 则常常与“排列”有关。例如，对 \(n\) 个物品进行分组，在确定了一个方案中分成了几组，及每个组的大小后，还要确定每个物品分到哪个组里面。
而这个可以用 \(\dfrac{n!}{\prod_{i}\limits S_i!}\) 描述。（\(S_i\) 表示第 \(i\) 组的大小）也就是所有物品的全排列，除掉每组内部的全排列。
那么干脆把这个 \(\dfrac{1}{\prod_i S_i!}\) 放到每个组的生成函数里面去好了，也就是给大小为 \(t\) 的带上一个 \(\dfrac{1}{t!}\) 。
于是 \(x^t \to \dfrac{x^t}{t!}\)。统计方案数的时候记得把 \(n!\) 给乘上啊。（可以理解为用 \(\text{EGF}\) 运算，得到的也是 \(\text{EGF}\)，所以在计算系数的时候要把 \(n!\) 乘回去）
关于\(\text{NTT}\)
一般刚敲完板子的话，对 \(\text{NTT}\) 的理解都会比较奇怪。这里一并讲了。
\(\text{NTT}\) 允许我们把多项式转到一些特定的点值上去。
点值表示的多项式可以与常数或者其它点值表示的多项式直接进行加/减/乘。
\(\text{NTT}\) 允许我们把用一些特定的点值表示的多项式转回其系数表示。
多项式求逆
给定一个最高次为 \(n-1\) 的多项式 \(f(x)\) ，求一个多项式 \(g(x)\) ，其满足
你可能会对 \(\pmod{x^n}\) 产生困惑。这个意思是舍弃次数大于等于 \(n\) 的项数。（它们都是 \(x^n\) 的倍数）
解法：
考虑倍增。假如我们已经求出了 \(g_0(x)\) ，其满足 \(f(x)\cdot g_0(x)\equiv 1 \pmod{x^{\lceil \frac{n}{2} \rceil}}\)
左右同时平方的话，会得到什么？
右边显然是 \(0\)。左边的话，既然第 \({{\lceil \frac{n}{2} \rceil}}\) 项以下都是 \(0\) ，那么其自乘之后第 \(n\) 项以下也全是 \(0\)。
即：
展开并移项，可以得到
我们知道 \(f(x)\cdot g(x)\equiv 1 \pmod{x^n}\) 。将两边同时乘 \(g(x)\) ，得到：
这个运算可以用 \(\text{NTT}\) 实现。因此，可以递归求解。边界是在 \(n=1\) 的时候把 \([x^0]g(x)\) 写成 \([x^0]f(x)\) 的逆元。
（所以 \([x^0]f(x)=0\) 时 \(f(x)\) 没有逆元）
细节：
此处 \(f(x)\) 只需要前 \(n\) 项。\(g_0(x)\) 长度为 \(\lceil \dfrac{n}{2} \rceil\) ，平方后长度为 \(n\)
所以 \(\text{NTT}\) 的时候要记得总长度为大于等于 \(2n\) 的第一个二次幂。
同时，最后的结果只需要前 \(n\) 项就好了，因此记得把后面的清成 \(0\) 。
然后虽然每个人写的求逆都不一样，但是有很大概率你写的求逆需要保证 \(g(x)\) 初始为空。
复杂度
例题
算是生成函数入门？P4841 [集训队作业2013]城市规划
分治 \(\text{NTT}\)
有时候你会得到一个跟自己有关的式子。以模板为例：
你愿意的话可以多项式求逆带走它：
但是并不是所有这样的式子你都能用求逆带走的，所以这里还是讲一讲分治 \(\text{NTT}\) 这个重要的思想。
解法
法如其名，考虑分治。假如当前所处理的区间为 \([l,r]\)，中间点为 \(mid\)
那么先处理掉 \([l,mid]\) 这一部分。接下来，考虑 \([l,mid]\) 对 \([mid+1,r]\) 的影响。然后去处理 \([mid+1,r]\) 这一部分。
放在这道题上，就是把 \([l,mid]\) 这一部分在递归后已经完全处理好的 \(f(x)\)，和 \(g(x)\) 卷起来，将其结果加到 \([mid+1,r]\) 上面去就好了。可以注意到 \(g(x)\) 有用的只有前 \(r-l+1\) 项。
递归边界为 \(l=r\) ，此时直接返回。
举个例子，\(n=8\) 时：递归到 \([0,0]\) 的时候返回（此时 \([x^0]f(x)\) 已经可以确定），退回到 \([0,1]\) ，将 \([0,0]\) 范围内的 \(f(x)\) 与 \(g(x)\) 相卷，用其结果更新 \([1,1]\)；然后递归到 \([1,1]\)，直接返回（此时 \([x^1]f(x)\) 已经可以确定）；退回到 \([0,3]\) ，将 \([0,1]\) 范围内的 \(f(x)\) 与 \(g(x)\) 相卷，用其结果更新 \([2,3]\) 范围内的 \(f(x)\)；递归 \([2,3]\)，递归 \([2,2]\) ，返回（此时 \([x^2]f(x)\) 已经可以确定），退回到 \([2,3]\) ，用 \([2,2]\) 范围内的 \(f(x)\) 更新 \([3,3]\) 范围内的 \(f(x)\)，递归 \([3,3]\)，返回（此时 \([x^3]f(x)\) 已经可以确定），返回到 \([0,7]\) ，用 \([0,3]\) 范围内的 \(f(x)\) 与 \(g(x)\) 相卷更新 \([4,7]\) 范围内的 \(f(x)\)，递归 \([4,7]\) ，继续类似的过程。
也许你会发现，我们每递归到一个满足 \(l=r\) 的叶子，它的值就已经可以被完全确定。同时，每次我们用左边的 \(f(x)\) 去更新右边的的 \(f(x)\) 时，左边的 \(f(x)\) 也已经被完全确定。
同时，这样的更新不会带来重复。如果把递归的过程看成一棵以 \([0,n-1]\)为根，以 \([i,i]\) 为叶子的类似线段树的二叉树，对于两个点 \(i,j,i<j\) ，\(i\) 只会在二者的 \(\text{LCA}\) 处贡献到 \(j\) 上。
（实际上这就是 \(\text{CDQ}\) 分治的一个应用）
在处理更为复杂的式子时，细节可能会非常诡异。
复杂度
（这个常数挺小的）
例题
二项式定理
也就是先枚举 \(\alpha\) 被乘了几次（\(i\)），而这一项的系数就是从 \(t\) 个 \((\alpha+\beta)\) 中选出 \(i\) 个 \(\alpha\) 的方案数，即 \(\binom{t}{i}\)
求导&积分
在开始后面的内容之前，你可能需要一些关于这方面的知识。学过的话可以跳过。
先来看看导数的定义：
例如，对路程求导就可以得到瞬时速度。
再例如，尝试对 \(f(x)=x^2+1\) 求导：
由于 \(\varepsilon\) 是一个趋近于 \(0\) 的极小数，\(\varepsilon\) 可以忽略。即，\(f'(x)=2x\)
可以发现，常数项对求导没有影响。
接下来，我们看看一个多项式要怎么求导：一个多项式可以拆成多个单项式，将这些单项式的导数相加就是多项式的导数
那么，一个单项式的求导：
把 \((x+\varepsilon)^t\) 用二项式定理展开。其中， \(\varepsilon\) 的次数如果大于 \(1\)，那么这一项的 \(\varepsilon\) 就除不干净，结果里会带有一个 \(\varepsilon\)。由于 \(\varepsilon\) 是一个趋近于 \(0\) 的极小数，这一项可以忽略。
而如果 \(\varepsilon\) 的次数为 \(0\) ，它就会和 \(x^t\) 相消。因此，可以得到
验证一下，把 \(f(x)=x^2\) 带进去，可以得到 \(f'(x)=2x\)
于是多项式的求导可以简单实现。
积分，姑且可以理解成求导的逆运算。多项式的积分，通过上面的式子倒推回去，也可以简单实现。
注意求导再积回来会导致常数项信息丢失。需要自己看情况补回去。
复合函数求导法则
证明：
由于 \(g'(x)\) 有限，而有限个 \(\varepsilon\) 的和仍然是无限小，分子分母的两个 \(\varepsilon\cdot g'(x)\) 等价，所以上式等于
泰勒展开
后面的内容需要这个。这是我无法给出证明的东西。
对于一个函数 \(f\) ，一个在其定义域内，且在其可导的位置的 \(x_0\)，有
\(e\) 的一些性质
1，\(f(x)=\ln(x),f'(x)=\frac{1}{x}\)
证明：
注意到 \(\frac{\varepsilon}{x}\) 趋向 \(0\) ，且左边式子的次方趋向于\(\infty\)，再回顾 \(e\) 的定义，可以得到 \(\frac{\varepsilon}{x}\cdot \frac{1}{\varepsilon}\cdot \frac{1}{f'(x)}=1\)
即，\(f'(x)=\frac{1}{x}\)
2，\(\ln{(\varepsilon+1)}=\varepsilon\)
证明：
3，\(f(x)=e^x,f'(x)=e^x\)（真神奇）
证明：
最后一步是性质2。
多项式 \(\ln\)
给出一个最高次为 $ n-1$ 的多项式 \(f(x)\) ，求 \(g(x)\) 满足 \(\ln(f(x))\equiv g(x) \pmod{x^n}\)
保证 \([x^0]f(x)=1\)
关于 \(e\)
你可能会一脸蒙蔽，\(e\) 不是个数吗？它为什么可以有多项式次方？
我们可以把 \(e^x\) 在 \(x_0=0\) 的位置泰勒展开，得到
我们发现就算 \(x\) 是个多项式它也是有意义的，于是 \(e\) 的多项式次方就是合法的。
但是有一个限制，对于 \(\exp(f(x))\) 合法时 \([x^0]f(x)\) 必须为 \(0\) 。不然 \([x^0]\exp(f(x))\) 无法计算。原因显然。
于是 \([x^0]\exp(f(x))\) 的结果必然为 \(1\)。（只在 \(i=0\) 的时候会算到常数项）\([x^0]\ln(f(x))\) 则必然为 \(0\) 。
（当然你要是比较特立独行，你可以去找指导谈谈，谈谈 \([x^0]\exp(f(x))\) 在 \([x^0]f(x)\) 不为 \(0\) 的时候是什么东西）
解法
对两边一起求导，得到
（\(\ln(f(x))\) 是个复合函数，要使用复合函数的求导法则）
然后我们发现求逆和求导我们都会，上面的式子我们会算，于是可以得到 \(g'(x)\)
然后给他积回去就好了。复杂度为 \(O(n\log n)\) 。
牛顿迭代
给定一个多项式 \(f(x)\) ，求 \(g(x)\) 使 \(f(g(x)) \equiv 0 \pmod {x^n}\) 。
解法
假如我们已经成功求出了 \(f(g_0(x)) \equiv 0 \pmod {x^{\lceil \frac{n}{2} \rceil}}\) ，
那么将其泰勒展开：
看着非常的阴间，但是有一个好消息是， \(g(x)\) 跟 \(g_0(x)\) 在前 \(\lceil \dfrac{n}{2} \rceil\) 项都是相同的。
（\(\lceil \dfrac{n}{2} \rceil\) 及以上的项数不会对 前 \(\lceil \dfrac{n}{2} \rceil\) 项有影响）
那么我们发现上面那个式子在 \(i \ge 2\) 的时候 \(\pmod {x^n}\) 就都是 \(0\) 了。
好的，所以原式等价于
这样就可以愉快的递归了。当然，递归到底层的时候，要看情况处理掉 \([x^0]f(g_0(x)) = 0\)
细节
注意，\(f(g(x))\) 是关于 \(g(x)\) 的函数，而不是关于 \(x\) 的复合函数。
例如，如果有 \(f(g(x))=g(x)+x^2\)，那么 \(f'(g(x))=1\)
多项式 \(\exp\)
给出一个最高次为 \(n-1\) 的多项式 \(\alpha(x)\) ，求 \(g(x)\) 使 \(\exp(\alpha(x))=g(x)\)
保证 \([x^0]\alpha(x)=0\)
解法
两边一起取 \(\ln\) 之后移项，得到 \(\ln(g(x))-\alpha(x)=0\) 。
那么令 \(f(g(x))=\ln(g(x))-\alpha(x)\) ，只要求出 \(g(x)\) 满足 \(f(g(x))\equiv 0\pmod{x^n}\) 就好了。
\(\alpha(x)\) 与 \(g(x)\) 无关。
把它喂到牛迭的式子里面，得到
这个我们就会了。所以递归就好。边界如上文所述，\([x^0]\exp(f(x))\) 必然为 \(1\)。
复杂度是 \(T(n)=T(\dfrac{n}{2})+O(n\log n)=O(n\log n)\)
警告
注意，做 \(\exp\) 需要递归，每次递归需要 \(\ln\) 和 \(\text{NTT}\)。
然后，做 \(\ln\) 需要求逆。
接着，做求逆需要递归，每次递归需要\(\text{NTT}\)。
所以这个 \(O(n\log n)\) 常数非常大。
另外，跟求逆一样，你的 \(\exp\) 也很有可能需要保证初始为空。
应用
你现在也许会好奇，这有什么用啊？
首先，这允许我们化乘为加：
推式子的时候常用。例如，P4389 付公主的背包
也可以拿去 \(n\log n\) 多项式快速幂。
你以为它就只能干这点事了？不不不， \(\exp\) 比这猛多了。
考虑这样一个问题：P5748 集合划分计数
简单的描述一下，就是把 \(n\) 个不同的球分到任意多个集合里面，集合之间无序，求方案数。
由于每个集合内部球的排列没有顺序，也就是说，对于每个已经确定大小，且已经确定里面有哪些球的集合，其方案数为 \(1\)。
于是我们可以写出每个集合的 \(\text{EGF}\) ：
为什么是 \(\text{EGF}\) ？因为我们要确定每个球分到哪个集合里面。为什么第 \(0\) 项是 \(0\) ？因为是非空集合。
那么接下来，我们枚举一共有几个集合，接下来用乘法合并这些集合的生成函数，并除去集合数量的阶乘，得到集合划分数的\(\text{EGF}\) \(g\)。（除以集合数量的阶乘是因为集合之间无序）
这好像没法做？但是……有点眼熟。
再次观察 \(\exp\) 的展开式。
那么！\(g(x)=\exp(f(x))\) 。于是你可以去把上面的题切了。
总结一下，\([x^t]\exp(f(x))\) 表示将 \(t\) 个不同物品分成任意多组，其中大小为 \(k\) 的一组有 \([x^k]f(x)\) 种方案的方案数。
（\(f(x)\) 是 \(\text{EGF}\)）
\(\text{Euler}\) 变换
现在我们有了 \(\exp\) 这个强有力的工具来数有标号的计数问题，那么怎么数无标号的呢？
例题：无标号有根树计数。解法在这
解法写的是无标号无根树计数，包括了无标号有根树计数，看跟无标号有根树计数相关的部分就好。
多项式开根
这个感觉不太重要，但还是讲了吧
给定 \(n-1\) 次多项式 \(f(x)\)，求 \(g(x)\) ，满足 \(f(x)\equiv g^2(x)\pmod{x^n}\)
解法
还是倍增，考虑已经求得了 \(g_0(x)\)，其满足 \(f(x)\equiv g_0^2(x)\pmod{x^{\lceil \frac{n}{2} \rceil}}\)
然后是在推导求逆时类似的操作，有兴趣的话可以自己推推看
这个式子可以用求逆跟 \(\text{NTT}\) 算出来，然后就可以递归带走了。
边界
如果 \([x^0]f(x)=1\) 就是 \(1\)，\([x^0]f(x)=0\) 就是 \(0\)。
否则的话，需要保证 \([x^0]f(x)\) 是 \(\text{998244353}\) 下的二次剩余（就是 \(([x^0]f(x))^2\equiv 1 \pmod {\text{998244353}}\)要有解），接着用一个奇怪的算法求解。
考虑到不太重要，实际应用也不大，我不打算讲解那个奇怪的算法。
但是，跟求逆一样，你写的开根也很有可能需要保证初始为空。
也许现在你可以去这个题单里面玩了？去应用一下。 感谢x义x神仙
感受数数的毒瘤吧……
至此，你大约可以宣称自己多项式入门了？
最后是说好的代码。
#include <vector>
#include <stdexcept>
#include <cstring>
#include <iostream>
#include <algorithm>
const int maxn = 3e5+5 , mod = 998244353;
template<typename T>
inline T max(const T &a,const T &b){
return a>b?a:b;
}
template<typename T>
inline void swap(T &a,T &b){
T temp=a;a=b;b=temp;
}
int fastpow(int x,int y){
if(y==0) return 1;int tmp=fastpow(x,y>>1);
return y&1?1ll*tmp*tmp%mod*x%mod:1ll*tmp*tmp%mod ;
}
inline int add(int x,int y){return x+y>=mod?x+y-mod:x+y;}
inline int sub(int x,int y){return x-y<0?x-y+mod:x-y;}
//Alpha -> Beta
inline void copy(int *Alpha,const int *Beta,int len){
for(int i=0;i<len;i++) Alpha[i]=Beta[i];
}
namespace Poly{
const int gate = 3,invg = fastpow(gate,mod-2);
int lim,lg,sing[maxn];int Gate[maxn],Invg[maxn];int inv[maxn];
void init(int len){
lim = 1,lg = 0;while(lim<len) lim<<=1,lg++;
for(int i=1;i<lim;i++) sing[i]=(sing[i>>1]>>1)|((i&1)<<(lg-1));
}
int Hacking_to_the(){
for(int i=1;i<maxn;i<<=1) Gate[i] = fastpow(gate,(mod-1)/(i<<1)),Invg[i] = fastpow(invg,(mod-1)/(i<<1));
inv[1] = 1 ;for(int i=2;i<maxn;i++) inv[i] = (mod-1ll*(mod/i)*inv[mod%i]%mod) ;
return 0 ;
}
int GATE = Hacking_to_the();
void NTT(int *Steins,int type){
for(int i=0;i<lim;i++) if(i<sing[i]) swap(Steins[i],Steins[sing[i]]);
for(int len=1;len<lim;len<<=1){
int unit = type==1?Gate[len]:Invg[len];
for(int i=0;i<lim;i+=(len<<1)){
int w = 1;
for(int k=0;k<len;k++,w=1ll*w*unit%mod){
int x=Steins[i+k],y=1ll*w*Steins[i+k+len]%mod;
Steins[i+k]=add(x,y),Steins[i+k+len]=sub(x,y);
}
}
}
if(type!=1){int Okabe=fastpow(lim,mod-2);for(int i=0;i<lim;i++) Steins[i]=1ll*Steins[i]*Okabe%mod;}
}
int ALPHA[maxn],BETA[maxn];
void convolution(const int *Alpha,int la,const int *Beta,int lb,int *Zeta){
init(la+lb-1);copy(ALPHA,Alpha,la),copy(BETA,Beta,lb);for(int i=la;i<lim;i++) ALPHA[i]=0;for(int i=lb;i<lim;i++) BETA[i]=0;
NTT(ALPHA,1),NTT(BETA,1);
for(int i=0;i<lim;i++) Zeta[i]=1ll*ALPHA[i]*BETA[i]%mod;return NTT(Zeta,-1) ;
}
int g[maxn];
void get_inv(const int *Alpha,int len,int *Beta){
//使用前清空Beta数组
if(len==1){Beta[0]=fastpow(Alpha[0],mod-2);return ;}
get_inv(Alpha,(len+1)>>1,Beta);init(len+len-1);
copy(g,Alpha,len);for(int i=len;i<lim;i++) g[i] = 0;
NTT(Beta,1),NTT(g,1);
for(int i=0;i<lim;i++) Beta[i]=1ll*Beta[i]*(mod+2-1ll*Beta[i]*g[i]%mod)%mod;
NTT(Beta,-1);
for(int i=len;i<lim;i++) Beta[i]=0;
}
const int inv2=fastpow(2,mod-2);
int t[maxn];
void get_sqrt(const int *Alpha,int len,int *Beta){
//使用前清空Beta数组
if(len==1){Beta[0]=1;return ;}get_sqrt(Alpha,(len+1)>>1,Beta);
get_inv(Beta,len,t);
init(len+len-1);copy(g,Alpha,len);for(int i=len;i<lim;i++) g[i]=0;
convolution(g,len,t,len,g);
for(int i=0;i<lim;i++) Beta[i]=1ll*inv2*(Beta[i]+g[i])%mod;
for(int i=len;i<lim;i++) Beta[i]=0;for(int i=0;i<lim;i++) t[i]=0;
}
void derivative(const int *Alpha,int len,int *Beta){
for(int i=1;i<len;i++) Beta[i-1]=1ll*i*Alpha[i]%mod;
Beta[len-1] = 0;
}
void integral(const int *Alpha,int len,int *Beta){
for(int i=len-1;i>=1;i--) Beta[i]=1ll*Alpha[i-1]*inv[i];
Beta[0] = 0;
}
int f1[maxn],f2[maxn];
void get_ln(const int *Alpha,int len,int *Beta){
memset(f2,0,sizeof(f2));
derivative(Alpha,len,f1),get_inv(Alpha,len,f2);
convolution(f1,len,f2,len,Beta);
return integral(Beta,len,Beta) ;
}
int e1[maxn],e2[maxn];
void get_exp(const int *Alpha,int len,int *Beta){
//使用前清空Beta数组
if(len==1) return void(Beta[0]=1);
get_exp(Alpha,(len+1)>>1,Beta);get_ln(Beta,len,e1);
init(len+len-1);
for(int i=0;i<len;i++) e2[i] = sub(Alpha[i],e1[i]);for(int i=len;i<lim;i++) e2[i]=0;
NTT(e2,1),NTT(Beta,1);
for(int i=0;i<lim;i++) Beta[i]=1ll*Beta[i]*(e2[i]+1)%mod;
NTT(Beta,-1);
for(int i=len;i<lim;i++) Beta[i]=0;for(int i=0;i<lim;i++) e1[i]=0;
}
int g1[maxn],g2[maxn],g3[maxn];
void division(const int *Alpha,int la,const int *Beta,int lb,int *Zeta,int *Gamma){
//Zeta's length requires 2la-lb
if(la<lb){
for(int i=0;i<=la-lb;i++) Zeta[i] = 0;
for(int i=la;i<lb;i++) Gamma[i] = 0;return copy(Gamma,Alpha,la);
}
for(int i=0;i<la;i++) g1[i]=Alpha[la-i-1];for(int i=0;i<lb;i++) g2[i]=Beta[lb-i-1];
for(int i=lb;i<=la-lb;i++) g2[i]=0;
lim=1;while(lim<=(la-lb)<<1) lim<<=1;for(int i=0;i<lim;i++) g3[i]=0;
get_inv(g2,la-lb+1,g3);
convolution(g1,la,g3,la-lb+1,Zeta);for(int i=0;i<la-lb-i;i++) swap(Zeta[i],Zeta[la-lb-i]);
for(int i=la-lb+1;i<lim;i++) Zeta[i]=0;
copy(g1,Beta,lb);copy(g2,Zeta,la-lb+1);
convolution(g1,lb,g2,la-lb+1,g1);for(int i=0;i<lb;i++) Gamma[i]=sub(Alpha[i],g1[i]);
}
}
struct poly{
std::vector<int> f;
poly():f(){}
poly(const std::vector<int> &F):f(F){}
poly(const std::initializer_list<int>& lit)
: poly(std::vector<int>(lit)){}
poly(const int &length): f(length){}
int operator[](const int &id) const {return f[id];}
int& operator[](const int &id){return f[id];}
int* data(){return f.data();}
const int* data() const {return f.data();}
std::vector<int>::size_type size() const {return f.size();}
poly copy() const {return poly(f);}
poly resize(const int &length){f.resize(length);return *this;}
poly reverse() {std::reverse(f.begin(),f.end());return *this;}
poly auto_resize(){while(f.size() > 1 && f.back() == 0 ) f.pop_back() ; return *this;}
void print() const {
for(int i = 0 ; i < f.size() ; i ++ ) printf("%d%c",f[i]," \n"[i==f.size()-1]) ;
}
//左闭右开
poly SUB(const int &l,const int &r = -1) const {
if(r == -1) return poly(std::vector<int>(f.begin()+l,f.end()));
return poly(std::vector<int>(f.begin()+l,f.begin()+r));
}
poly operator<<(const int &k) const {
if(f.size()+k<0) throw std::runtime_error("f.size()+k<0 but calculate f << k") ;
poly r(f.size()+k) ;
for(int i = 0 ; i < f.size() ; i ++ ) r[i+k] = f[i] ;
return r ;
}
poly operator>>(const int &k) const {
if(f.size() < k) throw std::runtime_error("f.size() < k but calculate f >> k") ;
poly r(f.size()-k) ;
for(int i = k ; i < f.size() ; i ++ ) r[i-k] = f[i] ;
return r ;
}
poly operator+(const poly &g) const {
poly r(max(f.size(),g.size())) ; int ms = std::min(f.size(),g.size()) ;
for(int i = 0;i < ms; i ++ ) r[i] = add(f[i],g[i]) ;
for(int i = ms;i < f.size(); i ++ ) r[i] = f[i] ;
for(int i = ms;i < g.size(); i ++ ) r[i] = g[i] ;
return r ;
}
poly operator-(const poly &g) const {
poly r(max(f.size(),g.size())) ; int ms = std::min(f.size(),g.size()) ;
for(int i = 0;i < ms; i ++ ) r[i] = sub(f[i],g[i]) ;
for(int i = ms;i < f.size(); i ++ ) r[i] = f[i] ;
for(int i = ms;i < g.size(); i ++ ) r[i] = mod-g[i] ;
return r ;
}
poly operator*(const poly &g) const {
poly r((f.size()+g.size()) << 1) ;
Poly::convolution(f.data(),f.size(),g.data(),g.size(),r.data()) ;
return r.auto_resize(),r ;
}
poly inv() const {
poly r( f.size() << 2 ) ;
Poly::get_inv(f.data(),f.size(),r.data()) ;
return r.auto_resize(),r ;
}
};
