这篇博客写得很好,本文大部分参考此博客。
一些证明的细节在本文中不会提及。
前置知识:多项式求逆,多项式求导、积分 , 泰勒展开。
多项式ln
给定(G(x)),求
[F(x)=ln G(x)
]
两边求导,得
[F'(x)=frac{G'(x)}{G(x)}
]
[F(x)=int frac{G'(x)}{G(x)}
]
求逆,计算即可。
牛顿迭代
以下我们把一个多项式当做一个数来理解。
设(delta(F(x)))是一个关于多项式(F(x))的一个函数,现在我们要求一个(F(x)),满足
[delta(F(x))=0
]
假设我们已经求出了(F(x))的前(frac{n}{2})项,记为(F_0(x)),现在我们要求它的前(n)项
对(delta(F(x)))进行泰勒展开,有
[delta(F(x))=delta(F_0(x))+frac{delta '(F_0(x))}{1!}(F(x)-F_0(x))+frac{delta ''(F_0(x))}{2!}(F(x)-F_0(x))^2+...
]
其中(delta '(F(x)))为(delta (F(x)))关于(F(x))的导数
因为(F_0(x))与(F(x))前(frac{n}{2})项相同,所以((F(x)-F_0(x))^2equiv 0 (mod x^n)),所以上式从第三项开始的值全为0。
也就是
[delta(F(x))=delta(F_0(x))+delta '(F_0(x))(F(x)-F_0(x))=0
]
整理得
[F(x)=F_0(x)-frac{delta(F_0(x))}{delta '(F_0(x))}
]
如果没搞懂的话可以参照下面的例子来理解。
牛顿迭代——多项式求逆
给定(G(x)),求(F(x))满足
[F(x)G(x)equiv1 (mod x^n)
]
就是要使得(F(x)G(x)-1=0)
所以设 (delta(F(x))=F(x)G(x)- 1)
就有
[F(x)=F_0(x)-frac{delta(F_0(x))}{delta '(F_0(x))}
]
[=F_0(x)-frac{F_0(x)G(x)-1}{G(x)}
]
[=F_0(x)-F_0(F_0(x)G(x)-1)
]
[=2F_0(x)-F_0^2(x)G(x) (mod x^n)
]
牛顿迭代——多项式开根
即求(F^2(x)-G(x)=0)
有
[F(x)=F_0(x)-frac{F_0^2(x)-G(x)}{2F_0(x)}
]
[=frac{F_0^2(x)+G(x)}{2F_0(x)}
]
多项式Exp
[F(x)=e^ {G(x)}
]
[ln F(x)=G(x)
]
[ln F(x)-G(x)=0
]
[F(x)=F_0(x)-frac{ln F_0(x)-G(x)}{frac{1}{F_0(x)}}
]
[=F_0(x)-F_0(x)(ln F_0(x)-G(x))
]
[=F_0(x)(1-ln F_0(x)+G(x))
]
最后奉上我丑陋的代码模板:
#include<bits/stdc++.h>
using namespace std;
#define N 600007
#define M 2100007
#define ll long long
const ll mod=998244353;
const int lim=4e5;
ll inv[N];
int rev[M],len;
ll p2(ll x){return x*x%mod;}
ll pw(ll x,ll p)
{
return p?p2(pw(x,p/2))*(p&1?x:1)%mod:1;
}
void getrev()
{
for(int i=0;i<len;i++)
rev[i]=rev[i>>1]>>1|(i&1?len>>1:0);
}
void getlen(int n)
{
for(len=1;len<=n;len<<=1);
getrev();
}
void NTT(ll *a,int op)
{
for(int i=0;i<len;i++)if(i<rev[i])swap(a[i],a[rev[i]]);
for(int i=1;i<len;i<<=1)
{
ll nw=pw(3,(mod-1)/(i<<1));
for(int j=0;j<len;j+=i<<1)
{
ll w=1;
for(int k=j;k<j+i;k++)
{
ll x=a[k],y=a[k+i]*w%mod;
a[k]=(x+y)%mod,a[k+i]=(x-y+mod)%mod;
w=w*nw%mod;
}
}
}
if(op<0)
{
reverse(a+1,a+len);
ll Inv=pw(len,mod-2);
for(int i=0;i<len;i++)a[i]=a[i]*Inv%mod;
}
}
ll inv_c[N];
void getinv(int p,ll *a,ll *b)
{
if(p==1)return a[0]=pw(b[0],mod-2),(void)1;
getinv((p+1)/2,a,b);
getlen(p<<1);
ll *c=inv_c;
for(int i=0;i<p;i++)c[i]=b[i];
for(int i=p;i<len;i++)c[i]=0;
NTT(a,1),NTT(c,1);
for(int i=0;i<len;i++)a[i]=a[i]*(2-a[i]*c[i]%mod+mod)%mod;
NTT(a,-1);
for(int i=p;i<len;i++)a[i]=0;
}
ll mul_c[M],mul_d[M];
void mul(ll *t,ll *a,ll *b)
{
ll *c=mul_c,*d=mul_d;
for(int i=0;i<len;i++)c[i]=a[i],d[i]=b[i];
NTT(c,1),NTT(d,1);
for(int i=0;i<len;i++)c[i]=c[i]*d[i]%mod;
NTT(c,-1);
for(int i=0;i<len;i++)t[i]=c[i];
}
void devir(ll *a,int n)
{
for(int i=1;i<=n;i++)a[i-1]=a[i]*i%mod;
a[n]=0;
}
void inter(ll *a,int n)
{
for(int i=n;i>=0;i--)a[i+1]=a[i]*inv[i+1]%mod;
a[0]=0;
}
ll ln_c[N],ln_d[N];
void getln(int n,ll *a,ll *b)
{
ll *c=ln_c,*d=ln_d;
getlen(2*n);
for(int i=0;i<len;i++)c[i]=b[i],d[i]=0;
getinv(n,d,c),devir(c,n);
getlen(2*n);
mul(a,c,d);
inter(a,n);
for(int i=n;i<len;i++)a[i]=0;
}
ll exp_c[N];
void getexp(int p,ll *a,ll *b)
{
if(p==1)return a[0]=1,(void)1;
getexp((p+1)/2,a,b);
ll *c=exp_c;
getlen(2*p);
for(int i=0;i<len;i++)c[i]=0;
getln(p,c,a);
for(int i=0;i<p;i++)c[i]=(b[i]-c[i]+mod)%mod;
c[0]=(c[0]+1)%mod;
getlen(2*p);
mul(a,a,c);
for(int i=p;i<len;i++)a[i]=0;
}
ll sqrt_c[N],sqrt_d[N];
void getsqrt(int p,ll *a,ll *b)
{
if(p==1)return a[0]=1,(void)1;
getsqrt((p+1)/2,a,b);
ll *c=sqrt_c,*d=sqrt_d;
getlen(2*p);
for(int i=0;i<len;i++)c[i]=d[i]=0;
for(int i=0;i<p;i++)d[i]=b[i];
getinv(p,c,a);
getlen(2*p),mul(c,c,d);
for(int i=0;i<p;i++)a[i]=(a[i]+c[i])*inv[2]%mod;
}
ll a[M],b[M],c[M];
int n;
void work_mul()
{
int n,m;
scanf("%d%d",&n,&m);
for(int i=0;i<=n;i++)scanf("%lld",&a[i]);
for(int i=0;i<=m;i++)scanf("%lld",&b[i]);
getlen(n+m);
mul(c,a,b);
for(int i=0;i<=n+m;i++)printf("%lld ",c[i]);
}
void read()
{
scanf("%d",&n);
for(int i=0;i<n;i++)scanf("%lld",&a[i]);
}
void write(ll *a)
{
for(int i=0;i<n;i++)printf("%lld ",a[i]);
}
void work_inv()
{
read();
getinv(n,b,a);
write(b);
}
void work_sqrt()
{
read();
getsqrt(n,b,a);
write(b);
}
void work_ln()
{
read();
getln(n,b,a);
write(b);
}
void work_exp()
{
read();
getexp(n,b,a);
write(b);
}
void Init()
{
inv[1]=1;
for(int i=2;i<=lim;i++)
inv[i]=mod-mod/i*inv[mod%i]%mod;
}
int main()
{
//freopen("data.in","r",stdin);
Init();
work_exp();
return 0;
}