基础多项式

形式幂级数(sumlimits_{i=0}^{+infty}a_ix^i)和多项式具有很强的相似性，这里不做区分。
至于收敛不收敛的问题，我们不讨论。

多项式加法

略。

多项式减法

略。

多项式求相反数

略。

多项式积分

(int f(x)mathrm dx=sumlimits_{n=1}^{+infty}frac{[x^{n-1}]f(x)}nx^n)

多项式求导

(f^{(0)}=f(x),f^{(n)}(x)=sumlimits_{i=0}^{+infty}(i+1)[x^{i+1}]{f^{(n-1)}(x)}x^i)

多项式差分

(Delta^0f(x)=f(x),Delta^nf(x)=Delta^{n-1}f(x+1)-Delta^{n-1}f(x))
(Delta^nf(x)=sumlimits_{k=0}^n(-1)^{n-k}{nchoose k}f(x+k))

多项式乘法

给定(f(x)=sum_{i=0}^{n-1}a_ix^i,g(x)=sum_{i=0}^{m-1}b_ix^i)，求(h(x)=f(x)g(x))。

做法有可以进行浮点数运算的FFT和可以进行取模运算的NTT。

FFT

基本流程

DFT（系数表示法转点值表示法）( ightarrow)直接相乘( ightarrow)IDFT（点值表示法转系数表示法）。

DFT

我们现在假如你对(A(x)=sum_{i=0}^{n-1}a_ix^i)进行DFT。其中(n=2^t)。
按下标的奇偶性分类并把右边提一个(x)出来。
(A_0(x)=sum_{i=0}^{frac n2-1}a_{2i}x^i,A_1(x)=sum_{i=0}^{frac n2-1}a_{2i+1}x^i)
那么(A(x)=A_0(x^2)+xA_1(x^2))。
对于(k<frac n2)，将(w_n^k)以及(w_n^{k+frac n2})分别代入(A(x))。
(A(w_n^k)=A_0(w_n^{2k})+w_n^kA_1(w_n^{2k})=A_0(w_{frac n2}^k)+w_n^kA_1(w_{frac n2}^k))
(A(w_n^{k+frac n2})=A_0(w_{frac n2}^k)-w_n^kA_1(w_{frac n2}^k))
发现只有后面一项的符号不同，所以我们可以分治处理。

迭代求法

手玩以后发现每个位置分治后的位置就是其原本位置二进制反转后的结果。
所以我们可以先把原序列先变换成分治后的序列，然后一层层向上合并。
预处理(i)的二进制翻转(rev[i])，如果(i<rev[i])就交换序列元素。

IDFT

把DFT中(A_1(x))前面的系数由(w_n^k)变为(overline{w_n^k})即可。
最后还要除一个(n)。

小优化

假如要求(H=FG)，令(P=F+iG)，求出(P^2=(F^2-G^2)+2FGi)。
这样我们可以只需要给(P)做一次DFT。

NTT

DFT

用(g^frac{P-1}n)来代替(w_n)。
必须满足(P=rcdot2^k+1)。这样的(P)称为费马素数。
一般是(P=998244353,g=3)。

IDFT

用(a^{-1})代替(overline{a})。
（因为(w_n^k=1)，所以单位根的共轭实质上就是求逆）

多项式乘法逆

给定(F(x))，求(G(x)equiv F(x)^{-1}pmod {x^n})即(F(x)G(x)equiv1pmod {x^n})。

有解的充要条件：([x^0]F(x))存在(pmod P)意义下的逆元。
考虑倍增。
首先我们有(([x^0]F(x))^{-1}F(x)equiv1pmod{x})。
现在我们考虑已知(f(x)F(x)equiv1pmod{x^{lceilfrac n2 ceil}})如何求(g(x)F(x)equiv1pmod{x^n}))。
显然首先有(g(x)F(x)equiv1equiv f(x)F(x)pmod{x^{lceilfrac n2 ceil}})。
移项并且把(F(x))除掉可以得到(g(x)-f(x)equiv0pmod{x^{lceilfrac n2 ceil}})。
平方可以看做自己卷自己，而(forall i+j<n,min(i,j)<lceilfrac n2 ceil)，因此((g(x)-f(x))^2equiv0pmod{x^n})即(g(x)^2+f(x)^2+2f(x)g(x)equiv0pmod{x^n})。
两边同乘一个(F(x))可以得到(g(x)+f(x)^2F(x)+2f(x)equiv0pmod{x^n})即(g(x)equiv f(x)(2-f(x)F(x))pmod{x^n})。

多项式除法

分为带余除法和除法，除法就是乘一个乘法逆，所以我们现在只讨论带余除法。
给定(F(x),G(x))，求(P(x),Q(x))使得(F(x)=P(x)G(x)+Q(x)(deg(Q)<deg(G)))。保证(deg(G)<deg(F))。

设(deg(F)=n,deg(G)=m)，易知(deg(P)=n-m,deg(Q)<m)。
记(F_R(x)=x^{deg(F)}F(x^{-1}))，即将(F(x))的系数反位。
要求的(P(x),Q(x))满足(F(x)=P(x)G(x)+Q(x))。
那么我们有(F(x^{-1})=P(x^{-1})G(x^{-1})+Q(x^{-1}))。
两边同乘(x^n)得到(F_R(x)=G_R(x)P_R(x)+x^{n-m+1}Q_R(x))。
加个(pmod{x^{n-m+1}})，(F_R(x)equiv G_R(x)P_R(x)pmod{x^{n-m+1}})。
那么(P_R(x)equiv G_R(x)F_R(x)^{-1}pmod{x^{n-m+1}})，恰好(deg(P)=n-m)。
(Q(x)=F(x)-G(x)P(x))

多项式开方

给定(F(x))，求(G(x)=sqrt{F(x)}pmod{x^n})即(G(x)^2equiv F(x)pmod{x^n})。

有解的条件：看下面。
先不考虑直接(Ln+Exp)的作法，考虑倍增。
首先我们有((sqrt{[x^0]F(x)})^2equiv F(x)pmod{x})。
需要Cipolla算法。
现在考虑已知(f(x)^2equiv F(x)pmod{x^{lceilfrac n2 ceil}})如何求(g(x)^2equiv F(x)pmod{x^n})。
类似于多项式求逆，我们可以得到(f(x)equiv g(x)pmod{x^{lceilfrac n2 ceil}})即
(f(x)-g(x)equiv0pmod{x^{lceilfrac n2 ceil}})。
根据多项式求逆那里的证明我们可以得到(f(x)^2+g(x)^2+2f(x)g(x)equiv0pmod{x^n})。
也就是(2f(x)g(x)+f(x)^2+F(x)equiv0pmod{x^n})。
那么(g(x)equiv(f(x)^2+F(x))(2f(x))^{-1}pmod{x^n})。
闲着无聊再拆一步可以得到(g(x)equiv 2^{-1}f(x)+F(x)(2f(x))^{-1}pmod{x^n})。

多项式Newton迭代

已知(G(x))，求(F(x))使得(G(F(x))equiv0pmod{x^n})。

首先我们求出(f(x))使得(G(f(x))equiv0pmod{x})。
然后考虑已知(f(x))使得(G(f(x))equiv0pmod{x^{lceilfrac n2 ceil}})如何求(G(g(x))equiv0pmod{x^n})。
将(G(g(x)))在(f(x))处Taylor展开。
(G(g(x))equivsumlimitsfrac{G^{(i)(f(x))}}{i!}(g(x)-f(x))^ipmod{x^n})
需要注意(G^{(n)}(f(x))=frac{d^nG}{dx^n}|_{x=f(x)})。
我们知道(G(g(x))equiv G(f(x))equiv0pmod{x^{lceilfrac n2 ceil}})。
稍微分析一下可以得到((g(x)-f(x))equiv0pmod{x^{lceilfrac n2 ceil}})。
那么((g(x)-f(x))^2equiv0pmod{x^n})。
也就是说上面那个和式中(i>2)的项全都是(0)。
即(G(g(x))equiv G(f(x))+G'(f(x))(g(x)-f(x))pmod{x^n})。
因为(G(g(x))equiv0pmod{x^n})，所以(g(x)equiv f(x)-G(f(x))(G'(f(x)))^{-1}pmod{x^n})。
看上去好像很美好，实际上有很多限制。

多项式对数函数

给定(F(x))，求(G(x)equivln F(x)pmod{x^n})。

有解的充要条件：([x^0]F(x)=1)。
两边对(x)求导(G'(x)equiv F'(x)F(x)^{-1}pmod{x^n})。
即(G(x)=int F'(x)F(x)^{-1}pmod{x^n}).

多项式指数函数

给定(F(x))，求(G(x)equivexp{F(x)}pmod{x^n})。

有解的充要条件：([x^0]F(x)=0)。
我们知道(ln G(x)equiv F(x)pmod{x^n})。
设(H(x)=ln G(x)-F(x))，要求的就是(H(x))在(pmod{x^n})意义下的零点。
套到Newton迭代的板子里面，得到最后的那个递推式(g(x)equiv f(x)(1-ln f(x)+F(x))pmod{x^n})。

这里说一下两个有解的充要条件。（上面开根和求逆的太显然了就不说了）
(ln)的定义是基于(exp)的，所以我们先看(exp)。
(exp)在多项式意义下的定义是基于Taylor级数的，即(exp f(x)=sumfrac{f(x)^i}{i!})。
若([x^0]f(x) e0)，那么(exp f(x)=+infty)无法在(mathbb{F_p})中表示。
因此(exp f(x))有解的充要条件是([x^0]F(x)=0)。
(ln)的定义是(exp)的反函数，但是由于一些值域上的问题需要一些小小的变化，准确来说(ln (f(x)+1))是((exp f(x))-1)的反函数。
因为要满足(exp f(x))存在所以([x^0]f(x)=0)，那么(ln f(x))存在当且仅当([x^0]f(x)=1)。

多项式幂函数

给定(F(x))和(kin N+)，求(G(x)equiv F(x)^kpmod{x^n})。

一个比较优秀的做法是倍增快速幂，复杂度(O(nlog nlog P))。相比(O(nlog n))的做法常数小很多，最多慢一倍。
(k<0)可以先求(|k|)次幂然后求逆，(k=0)就是(1)，所以只需要考虑(k>0)的情况即(kinmathbb{N_+})。
易知(F(x)^k=exp(kln F(x)))。
注意到此时(k)可以对(P)取模。
但是如果([x^0]F(x) e1)是不能够求(ln)的。
如果([x^0]F(x) e 0)，我们可以考虑乘一个系数(lambda=([x^0]F(x))^{-1})让([x^0]lambda F(x)=1)，然后求((lambda F(x))^k)，最后再乘上(lambda^{-k})。
如果([x^0]F(x)=0)，那么我们可以提一个(x)出来，(F(x)=xF_1(x))，那么(F(x)^k=F_1(x)^kx^k)。一直提到([x^0]F_i(x) e0)为止。

多项式开高次方

给定(F(x))和(kinmathbb{N_+})，求(G(x)equivsqrt[k]{F(x)}pmod{x^n})即(G(x)^kequiv F(x)pmod{x^n})。

求出(t)，满足(forall iin[0,t-1],[x^i]F(x)=0wedge [x^t]F(x) e0)。
如果(k mid t)，那么无解。
否则令(F(x)=x^tf(x))，则(G(x)=x^{frac tk}sqrt[k]{f(x)})。
所以我们现在只要考虑([x^0]F(x) e0)的情况。
利用BSGS求出(pequivsqrt[k]{([x^0]F(x))}pmod P)。
令(F(x)=(p^kf(x)))，其中([x^0]f(x)=1)。
那么(G(x)=psqrt[k]{f(x)})。
也就是说每一个(p)对应一个(G(x))。
注意需要满足([x^0]sqrt[k]{f(x)}=1)。
而显然有(sqrt[k]{f(x)}=exp(frac{ln f(x)}k))。
（实际上最后一步与求出(Kequiv k^{-1}pmod P)然后用多项式快速幂计算(f(x)^K)等价）

多项式三角函数

给定(F(x))，求(G(x)equivsin F(x)pmod{x^n})和(G(x)equivcos F(x)pmod{x^n})。

有解的充要条件：([x^0]F(x)=0)。
考虑Euler公式(exp(iF(x))=cos F(x)+isin F(x))。
同理(exp(-iF(x))=cos F(x)-isin F(x))。
因此(cos F(x)equivfrac{exp(iF(x))+exp(-iF(x))}2pmod{x^n},sin F(x)equivfrac{exp(iF(x))-exp(-iF(x))}{2i}pmod{x^n})。
注意到(i^2equiv-1equiv P-1pmod P)，所以我们可以求出(P-1)的二次剩余来表示(i)。

多项式反三角函数

给定(F(x))，求(G(x)equivarcsin F(x)pmod{x^n},G(x)equivarccos F(x)pmod{x^n},G(x)equivarctan F(x)pmod{x^n})。

有解的充要条件：([x^0]F(x)=0)。
直接链式求导法则。
(arcsin F(x)equivfrac{F'(x)}{sqrt{1-F(x)^2}}pmod{x^n})
(arccos F(x)equiv-frac{F'(x)}{sqrt{1-F(x)^2}}pmod{x^n})
(arctan F(x)equivfrac{F'(x)}{1-F(x)^2}pmod{x^n})

多项式双曲函数

给定(F(x))，求(G(x)equivoperatorname{sh} F(x)pmod{x^n},G(x)equivoperatorname{ch} F(x)pmod{x^n})。

有解的充要条件：([x^0]F(x)=0)。
(operatorname{sh} F(x)=frac{exp F(x)-exp(-F(x))}2pmod{x^n})
(operatorname{ch} F(x)=frac{exp F(x)+exp(-F(x))}2pmod{x^n})

多项式反双曲函数

给定(F(x))，求(G(x)equivoperatorname{argsh}F(x)pmod{x^n},G(x)equivoperatorname{argch}F(x)pmod{x^n},G(x)equivoperatorname{argth}F(x)pmod{x^n})。

先看(operatorname{argsh},operatorname{argth})。
有解的充要条件：([x^0]F(x)=0)。
法一：根据定义，(operatorname{argsh}F(x)=ln(F(x)+sqrt{F(x)^2+1}),operatorname{argth}F(x)=frac12lnfrac{1+F(x)}{1-F(x)})。
法二：使用链式求导法则，(operatorname{argsh}F(x)=intfrac{F'(x)}{sqrt{F(x)^2+1}},operatorname{argth}F(x)=intfrac{F'(x)}{1-F(x)^2})。
再看(operatorname{argch})。
有解的充要条件：
(1.[x^0]F(x)=1)。
(2.)设有(t,forall iin[1,t-1][x^i]F(x)=0wedge[x^t]F(x)=1)。
第二个条件可以结合(operatorname{ch})的Taylor展开、(operatorname{argch})的定义式、(operatorname{argch})的导数理解。
法一：根据定义，(operatorname{argch}F(x)=ln(F(x)+sqrt{F(x)^2-1}))。
法二：使用链式求导法则，(operatorname{argch}=intfrac{F'(x)}{sqrt{F(x)^2-1}})。
实际上法二就是法一往下推一步的得到的式子。

任意模数NTT

给定(A(x),B(x))，求(C(x)=A(x)B(x))，模数不一定是NTT质数。

三模NTT

搞三个NTT模数分别NTT，然后CRT合并。

拆系数FFT

一个很简单的想法是FFT。
但是FFT过程中数字会达到(P^2n)的级别，大概是(10^{23})，无法保证精度。
设(A=aS+c,B(x)=bS+d)，其中(S=sqrt P)。
那么(AB=abS^2+(ad+bc)S+cd)。
现在我们需要(4+3)次DFT，比三模NTT优秀。数字的大小不会超过(Pn)，大概是(10^{14})，精度基本得到了保证。
运用MTT/exMTT可以优化至(2+2)甚至(2+1.5)次DFT。

MTT

myy的论文中讲述了一种将(2)次DFT合并为(1)次DFT的方法。
假如我们要求(DFT(A),DFT(B))，令(P=A+iB,Q=A-iB)。
那么显然有(DFT(A)=frac{DFT(P)+DFT(Q)}{2},DFT(B)=ifrac{DFT(P)-DFT(Q)}{2})。
然后随便手推一下可以发现(DFT(Q)=operatorname{conj}operatorname{rev}DFT(P))。
那么我们就只需要一次DFT了。
运用至拆系数FFT中可以将其优化至(2+2)次DFT。

多项式多点求值

给定(F(x),a_1,cdots,a_m)，(forall iin[1,m])，求出(F(a_i))。

考虑这样一个事实：(F(a)equiv F(x)pmod{x-a})。
这样我们就可以用多项式取模求啦！
用分治优化这个过程：
假如我们已经求出(G(x)equiv F(x)pmod{prodlimits_{i=l}^r(x-a_i)})。
然后计算出(A_0=prodlimits_{i=l}^{mid}(x-a_i),A_1=prodlimits_{i=mid+1}^r(x-a_i))。
再计算(G(x)mod A_0(x),G(x)mod A_1(x))就好了。
复杂度(T(n)=2T(frac n2)+O(nlog n)=O(nlog^2n))。
这个分治过程可以用类似于线段树的结构来表示。

多项式快速插值

给定(n)个点((x_i,y_i))，求一个(n-1)次多项式(F(x))使得(forall iin[1,n],y_i=F(x_i))。

( ext{Part.1})

先看(A_i=prodlimits_{j e i}(x_i-x_j))。
令(g(x)=prodlimits_{i=1}^n(x-x_i))。
那么(A_i=limlimits_{x o x_i}frac{g(x_i)}{x-x_i})。
根据L'Hopital法则(A_i=g'(x_i))。
那么我们可以分治+NTT计算出(g(x))，然后求导得到(g'(x))，再多点求值算出每个(g'(x_i))即(A_i)。

( ext{Part.2})

现在(F(x)=sumlimits_{i=1}^nfrac{y_i}{A_i}prodlimits_{j e i}(x-x_j))。
考虑分治，令(F_{i,j}=sumlimits_{i=l}^rfrac{y_i}{A_i}prodlimits_{j e i}(x-x_j))。
那么有(F_{l,r}=F_{l,mid}prodlimits_{i=mid+1}^r(x-x_i)+F_{mid+1,r}prodlimits_{i=l}^{mid}(x-x_i))。

多项式复合

给定(F(x),G(x))，(deg F=n-1,deg G=m-1)，求(H(x)equiv F(G(x))equivsumlimits_{i=0}^{n-1}([x^i]F(x)G(x))^ipmod{x^n})。

分块暴力

不妨先令(T=lceilsqrt n ceil)。
(H(x)=sumlimits_{i=0}^{lfloorfrac nT floor}sumlimits_{j=0}^{T-1}[iT+j<n]([x^{iT+j}]F(x))G(x)^{iT+j} =sumlimits_{i=0}^{lfloorfrac nT floor}G(x)^{iT}sumlimits_{j=0}^{T-1}[iT+j<n]([x^{iT+j}]F(x))G(x)^j)
先预处理(forall iin[0,lfloorfrac nT floor],G(x)^{iT})和(forall iin[0,T),G(x)^i)，复杂度为(O(n(T+frac nT)log n))处理。
然后对于(sumlimits_{j=0}^{T-1}[iT+j<n]([x^{iT+j}]F(x))G(x)^j)部分暴力计算，复杂度为(O(n^2))。
最后再把刚才计算出来的一个个和(G(x)^{iT})相乘，复杂度为(O(nfrac nTlog n))。
很显然依旧是(T=sqrt n)时复杂度最优，取到(O(n^2+nsqrt nlog n))。
实际上(O(n^2))部分常数小跑得飞快，慢的部分反而是两个(O(nsqrt nlog n))的。

Brent-Kung算法

把(G(x))的前(m)项拆出来记为(G_t(x))，剩下的记为(G_r(x))。
然后把(F(G(x)))在(G_r(x))处Taylor展开，(F(G(x))=sumlimits_{i=0}^{+infty}frac{F^{(i)}(G_t(x))G_r(x)^i}{i!})。
因为我们要求的是(H(x)mod x^n)，注意到(x^t|G_r(x))，所以我们只需要知道这个展开式前(l=lceilfrac nt ceil)次项的值。
即(H(x)equivsumlimits_{i=0}^lfrac{F^{(i)}(G_t(x))G_r(x)^i}{i!}pmod{x^n})。
然后我们考虑利用分治计算(F(G_t(x)))。
假设(deg F=s)，取(F_1(x)equiv F(x)pmod{x^{lceilfrac s2 ceil}},F_2(x)equivfrac{F(x)-F_1(x)}{x^{lceilfrac s2 ceil}}pmod{x^{lceilfrac s2 ceil}})。
那么有(F(G_t(x))equiv F_1(G_t(x))+G_t(x)^{lceilfrac s2 ceil}F_2(G_t(x))pmod{x^n})。
令(A(s)=O(slog s))，(T(s))为计算(G_t^{lceilfrac s2 ceil})和(F(G_t(x)))的时间复杂度。
那么有(T(s)=2T(frac s2)+A(min(n,st)))。
计算得到(T(n)=O(ntlog^2n))。
接下来考虑计算(F^{(k)}(G_t(x)))。
令(h(x)=F(G_t(x)))。
我们知道有(F'(G_t(x))=h'(x)G'_t(x)^{-1})。
那么我们预处理(G_t(x)^{-1})，先用分治计算出(h(x)=F(G_t(x)))并随手算出(h'(x))，然后计算出(F'(G_t(x)))并把(F'(G_t(x)))赋给(h(x))。。。
这里的复杂度为(O(nllog n)=O(frac{n^2log n}t))。
计算每一项(G_r(x)^i)还是(O(frac{n^2log n}t))的，最后求和是(O(frac{n^2}t))的。
总复杂度为(O(frac{n^2log n}t+ntlog^2n))，当(t=sqrt{frac n{log n}})时取到最小值(O((nlog n)^{frac32}))。
然后这东西跑得比(O(n^2))的慢，然后就没有然后了。

多项式复合逆

给定(F(x),deg F=n-1)，求(G(x))使得(deg G=n-1wedge G(F(x))equiv x(mod x^n))。

有解的充要条件：([x^1]F(x))存在(pmod P)意义下的逆元且([x^0]F(x)=0)。
一次多项式需要特殊判断。
根据Lagrange反演，([x^i]G(x)=frac1{i-1}[x^{i-1}](frac x{F(x)})^i)。
类似于多项式复合函数的做法，考虑分块。取(T=sqrt n)。
则(G(x)equivsumlimits_{i=0}^{lceilfrac nT ceil}sumlimits_{j=0}^{T-1}[iT+j<n](frac1{iT+j}[x^{iT+j}](frac x{F(x)})^{iT+j})x^{iT+j}equivpmod{x^n})。
剩下的就跟多项式复合函数的过程一样了。

快速阶乘

计算(n!mod p)。

(O(sqrt nlog^2n))的暴力分块就不讲了。
设(B=sqrt n)，记(f(d,x)=prodlimits_{i=1}^d(x+i))。
我们最后要求的是({f(B,0),cdots,f(B,B^2)})。
我们知道(n+1)个点可以确定一个(n)次多项式，所以如果我们求出了({f(d,0),cdots,f(d,dB)})，那么我们就可以确定(f(d,x))。
考虑倍增。
具体而言我们需要实现两个东西：
(1.)已知({f(d,0),cdots,f(d,dB)})，求({f(2d,0),cdots,f(2d,2dB)})。
(2.)已知({f(d,0),cdots,f(d,dB)})，求({f(d+1,0),cdots,f(d+1,(d+1)B)})。
假如我们能够实现这两个东西的话我们就可以通过(O(log n))轮迭代求出({f(B,0),cdots,f(B,B^2)})。
这里给一个式子，这是下面推式子的基础。
(f(d+t,x)=f(d,x)f(d,x+t))

( ext{Part.1})

我们需要求出({f(d,0),cdots,f(d,2dB)})和({f(d,d),cdots,f(d,2dB+d)})。
构造(h(x)=f(d,Bx))。
那么问题变成已知({h(0),cdots,h(d)})，需要求出({h(d+1),cdots,h(2d+1)})和({h(dB^{-1}),cdots,h(dB^{-1}+2d)})。
我们可以把这个问题分为两步并且归为已知({h(0),cdots,h(k)})，需要求出({h(Delta),cdots,h(Delta+k)})。
考虑Lagrange插值，(f(Delta+n)=sumlimits_{i=0}^kh(i)prodlimits_{j e i}frac{Delta+n-j}{i-j}=prodlimits_{j=0}^n(Delta+n-j)sumlimits_{i=0}^kfrac{h(i)}{(Delta+n-i)i!(k-i)!(-1)^{k-i}})
（注意到(Delta+n-i e0)）
后面这一坨是个卷积。前面的东西是一段连续数字之积，可以双指针。

( ext{Part.2})

(f(d+1,(d+1)B))暴力算，剩下的利用上面的那个式子直接做就好了。

总的复杂度为(O(Blog B)=O(sqrt nlog n))。
一个小trick：根据Wilson定理，(n!equiv(-1)^{n+1}(p-n-1)!pmod p)。所以我们可以把(n)降到(frac p2)级别。

优化

(1. ext{negiizhao})的( ext{Newton})迭代优化Link。
(2. ext{rushcheyo})的转置原理多点求值。
(3. ext{skip2004})的( ext{cdq exp})。