luogu P4726 【模板】多项式指数函数 多项式 exp 牛顿迭代 泰勒展开

LINK：多项式 exp

做多项式的题 简直在嗑药。

前置只是 泰勒展开

这个东西用于 对于一个函数f(x) 我们不好得到 其在x处的取值。

所以另外设一个函数g(x) 来在x点处无限逼近f(x).

具体的 (f(x) ≈ g(x)=g(0)+frac{f^1(0)}{1!}x+frac{f^2(0)}{2!}x^2+...+frac{f^n(0)}{n!}x^n)

牛顿迭代：

常用来求一个函数的零点：假设我们已经求得一个近似值x0 那么我们只需要过(x0,f(x0))这个点做函数图像的切线 取切线与x轴的交点作为新的x0.

迭代几次就可以比较精确。

如：现在要求一个函数f(x) 近似值为x0 y=f'(x0)(x-x0)+f(x0);

y=0时 (x=x0-frac{f(x0)}{f'(x0)})

当然也可以放到多项式上 现在要求一个G(x) 我们想要求出F(G(x))=0的零点G(x).

(G(x)=G0(x)-frac{F(G0(x))}{F'(G0(x))})

本质上每迭代一次都可以迅速逼近真实值。

如果(F(G0(x))≡0(mod x^{2n}))那么(F(G(x))equiv 0(mod{x^{n}}))

关于这道题的推导过程：

(B(x)equiv e^{A(x)}(mod x^n))

(InB(x)-A(x)equiv 0(mod x^n))

现在设一个函数(F(G(x))=InG(x)-A(x)equiv 0)

其实在求F这个函数的零点.

两边直接求导可得 ((F(G0(x)))'=frac{G'0(x)}{G0(x)})

带入牛顿迭代的式子里。

(G(x)=frac{G0(x)(1-InG0(x)+A(x))}{G'0(x)})

每次迭代需要 求逆 做多项式In 再来一遍多项式乘法即可。

废什么话 码!

求导 和 积分需要仔细熟悉一下.

const int MAXN=600010,G=3;
int n;
int A[MAXN],B[MAXN],E[MAXN],F[MAXN],C[MAXN],D[MAXN],g[MAXN],rev[MAXN],inv[MAXN],O[MAXN];
inline int mul(int a,int b){return (ll)a*b%mod;}
inline int add(int a,int b){return a+b>=mod?a+b-mod:a+b;}
inline int mus(int a,int b){return a-b<0?a-b+mod:a-b;}
inline int ksm(int b,int p)
{
	int cnt=1;
	while(p)
	{
		if(p&1)cnt=mul(cnt,b);
		p=p>>1;b=mul(b,b);
	}
	return cnt;
}
inline void NTT(int *a,int op,int ww)
{
	int lim=1;while(lim<ww)lim=lim<<1;
	rep(0,lim-1,i)
	{
		rev[i]=rev[i>>1]>>1|((i&1)?lim>>1:0);
		if(i<rev[i])swap(a[i],a[rev[i]]);
	}
	for(int len=2;len<=lim;len=len<<1)
	{
		int mid=len>>1;
		int wn=ksm(G,op==1?(mod-1)/len:mod-1-(mod-1)/len);
		O[0]=1;rep(1,mid-1,i)O[i]=(ll)O[i-1]*wn%mod;
		for(int j=0;j<lim;j+=len)
		{
			for(int i=0;i<mid;++i)
			{
				int x=a[i+j],y=(ll)a[i+j+mid]*O[i]%mod;
				a[i+j]=(x+y)%mod;a[i+j+mid]=(x-y+mod)%mod;
			}
		}
	}
	if(op==-1)for(int i=0,inv=ksm(lim,mod-2);i<lim;++i)a[i]=(ll)a[i]*inv%mod;
}
inline void Direv(int *a,int *b,int len)//求导
{
	rep(0,len-2,i)b[i]=mul(a[i+1],i+1);b[len-1]=0;
}
inline void Inv(int *a,int *b,int len)
{
	if(len==1)return b[0]=ksm(a[0],mod-2),void();
	Inv(a,b,len>>1);rep(0,len-1,i)C[i]=a[i],D[i]=b[i];
	NTT(C,1,len<<1);NTT(D,1,len<<1);
	rep(0,(len<<1)-1,i)D[i]=mul(D[i],mul(D[i],C[i]));
	NTT(D,-1,len<<1);
	rep(0,len-1,i)b[i]=((ll)2*b[i]-D[i]+mod)%mod;
	rep(0,(len<<1)-1,i)C[i]=D[i]=0;
}
inline void Inter(int *a,int *b,int len)
{
	rep(1,len-1,i)b[i]=mul(a[i-1],inv[i]);b[0]=0;
}
inline void Ln(int *a,int *b,int len)
{
	Inv(a,E,len);Direv(a,F,len);
	NTT(E,1,len<<1);NTT(F,1,len<<1);
	rep(0,(len<<1)-1,i)E[i]=mul(E[i],F[i]);
	NTT(E,-1,len<<1);Inter(E,b,len);
	rep(0,(len<<1)-1,i)E[i]=F[i]=0;
}
inline void Exp(int *a,int *b,int len)
{
	if(len==1)return b[0]=1,void();
	Exp(a,b,len>>1);Ln(b,g,len);
	g[0]=(a[0]+1-g[0]+mod)%mod;
	rep(1,len-1,i)g[i]=mus(a[i],g[i]);
	NTT(g,1,len<<1);NTT(b,1,len<<1);
	rep(0,(len<<1)-1,i)b[i]=mul(b[i],g[i]);
	NTT(b,-1,len<<1);rep(len,(len<<1)-1,i)g[i]=b[i]=0;
}
int main()
{
	//freopen("1.in","r",stdin);
	get(n);rep(0,n-1,i)get(A[i]);
	int len=1;while(len<n)len=len<<1;
	inv[1]=1;rep(2,(len<<1),i)inv[i]=(ll)inv[mod%i]*(mod-mod/i)%mod;
	Exp(A,B,len);rep(0,n-1,i)put_(B[i]);return 0;
}