多项式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))^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;
}


原文地址:https://www.cnblogs.com/lishuyu2003/p/12123705.html