多项式ln、牛顿迭代学习笔记
这篇博客写得很好,本文大部分参考此博客。
一些证明的细节在本文中不会提及。
前置知识:多项式求逆,多项式求导、积分 , 泰勒展开。
多项式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))^2\equiv 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;
}