多项式Ⅱ
前言
因为之前的 Ⅰ 写的太屑了，所以又写了篇更屑的Ⅱ，其实是不同时间的个人垃圾理解了..
牛顿迭代
- 要求解\(f(x)=0\)\[\frac{f(x_0)}{x_0-x}=f'(x_0)\\ x=x_0-\frac{f(x_0)}{f'(x_0)} \]不断往下迭代可以找到近似解。
多项式牛顿迭代
和普通的牛顿迭代其实一样
\[F=F_0-\frac{G(F_0)}{G'(F_0)}
\]
如何理解近似解？
从\(\bmod x^i\)开始向上\(\bmod x^{i+1}...\)
- 多项式求逆\[A*B\equiv 1\pmod {x^n} \]构造\[G(B)=A*B-1 \]牛顿迭代找根
\[ \begin{aligned}
B&=B_0-\frac{AB_0-1}{A}\\
&=B_0-B_0(AB_0-1)\\
&=2B_0-AB_0^2
\end{aligned}
\]
向上倍增即可。
-
多项式开根
\[B^2\equiv A\pmod {x^n} \]构造
\[G(B)=A-B^2 \]牛顿迭代找根
\[\begin{aligned} B&=B_0-\frac{A-B_0^2}{-2B_0}\\ &=\frac{A}{2B_0}+\frac{B_0}{2} \end{aligned} \]向上倍增即可
-
多项式exp
\[B\equiv e^A \pmod {x^n} \]并没有弄清楚这个咋定义的，因此还不会直接用。
\[\ln B=A \]构造
\[G(B)=\ln B-A \]搞
\[\begin{aligned} B&=B_0-\frac{\ln B_0-A}{\frac{1}{B_0}}\\ &=B_0(1-\ln B_0+A) \end{aligned} \]
多项式取ln
-
多项式积分和求导
\[\begin{aligned} F(x)&=\sum_{i=0}^na_ix^i\\ \frac{\mathrm dF}{\mathrm dx}&=\sum_{i=0}^{n-1}a_{i+1}(i+1)x^i\\ \int F(x)&=\sum_{i=1}^n\frac{a_{i-1}}{i}x^i \end{aligned} \] -
多项式取ln
\[\begin{aligned} B&=\ln A\pmod {x^n}\\ B'&=\ln'(A)A'\pmod {x^n}\\ &=\frac{A'}{A}\pmod {x^n} \end{aligned} \] -
多项式快速幂
\[F^k=e^{k\ln n}\pmod {x^n} \]常数无敌的\(n\log n\)
拉格朗日反演和快速插值什么的，我还没学...
放一个exp的代码
提几个细节
- 不知道为什么递归版效率优于非递归版，而递归版又好写一点，因此以后就写递归版了
- 每个需要倍增的，每次要外面开两个自己的数组，不许改传进来的值，注意这两个数组一定需要清空
- 要注意一下每次NTT的长度什么的
- 不要封装卷积乘法，因为慢，可以把多个多项式乘起来全部先转换成点表达，然后一起乘，再转回系数，这样需要注意长度。
#include <cstdio>
#include <algorithm>
const int N=(1<<18)+10;
const int mod=998244353,Gi=332748118;
#define mul(a,b) (1ll*(a)*(b)%mod)
#define add(a,b) ((a+b)%mod)
int invA[N],invB[N],lnA[N],lnB[N],exA[N],exB[N],ans[N],f[N],turn[N],n;
int qp(int d,int k){int f=1;while(k){if(k&1)f=mul(f,d);d=mul(d,d),k>>=1;}return f;}
void NTT(int *a,int typ,int len)
{
int L=-1;for(int i=1;i<len;i<<=1) ++L;
for(int i=1;i<len;i++)
{
turn[i]=turn[i>>1]>>1|(i&1)<<L;
if(i<turn[i]) std::swap(a[i],a[turn[i]]);
}
for(int le=1;le<len;le<<=1)
{
int wn=qp(typ?3:Gi,(mod-1)/(le<<1));
for(int p=0;p<len;p+=le<<1)
{
int w=1;
for(int i=p;i<p+le;i++,w=mul(w,wn))
{
int tx=a[i],ty=mul(w,a[i+le]);
a[i]=add(tx,ty);
a[i+le]=add(tx,mod-ty);
}
}
}
if(!typ)
{
int inv=qp(len,mod-2);
for(int i=0;i<len;i++) a[i]=mul(a[i],inv);
}
}
void polyinv(int *a,int *b,int len)
{
if(len==1) {b[0]=qp(a[0],mod-2);return;}
polyinv(a,b,len>>1);
for(int i=0;i<len<<1;i++) invA[i]=invB[i]=0;//
for(int i=0;i<len;i++) invA[i]=a[i],invB[i]=b[i];
NTT(invA,1,len<<1),NTT(invB,1,len<<1);
for(int i=0;i<len<<1;i++) invA[i]=mul(invB[i],add(2,mod-mul(invA[i],invB[i])));
NTT(invA,0,len<<1);
for(int i=0;i<len;i++) b[i]=invA[i];
}
void direv(int *a,int n){for(int i=0;i<n;i++)a[i]=mul(a[i+1],i+1);a[n]=0;}
void inter(int *a,int n){for(int i=n;i;i--)a[i]=mul(a[i-1],qp(i,mod-2));a[0]=0;}
void polyln(int *a,int *b,int len)
{
for(int i=0;i<len;i++) lnA[i]=a[i],lnB[i]=b[i];
polyinv(lnA,lnB,len);
direv(lnA,len-1);
NTT(lnA,1,len<<1),NTT(lnB,1,len<<1);
for(int i=0;i<len<<1;i++) lnA[i]=mul(lnA[i],lnB[i]);
NTT(lnA,0,len<<1);
inter(lnA,len-1);
for(int i=0;i<len;i++) b[i]=lnA[i];
}
void polyexp(int *a,int *b,int len)
{
if(len==1){b[0]=1;return;}//默认a[0]=0
polyexp(a,b,len>>1);
for(int i=0;i<len<<1;i++) exA[i]=exB[i]=0;//
polyln(b,exA,len);
for(int i=0;i<len;i++) exA[i]=add(a[i]+(i==0),mod-exA[i]);
for(int i=0;i<len;i++) exB[i]=b[i];
NTT(exA,1,len<<1),NTT(exB,1,len<<1);
for(int i=0;i<len<<1;i++) exA[i]=mul(exA[i],exB[i]);
NTT(exA,0,len<<1);
for(int i=0;i<len;i++) b[i]=exA[i];
}
int main()
{
scanf("%d",&n);
for(int i=0;i<n;i++) scanf("%d",f+i);
int len=1;while(len<n) len<<=1;
polyexp(f,ans,len);
for(int i=0;i<n;i++) printf("%d ",ans[i]);
return 0;
}
2018.12.29
