多项式小结(求逆、求导、积分、ln、牛顿迭代、exp)

目录求逆求导复合函数求导积分ln牛顿迭代exp正确性证明n^2lnexpexpln快速幂时间复杂度调试方法&注意事项例题题解code

求逆

Ax)Bx)equiv 1mod;x^n)),下文为了方便表述把n/2

已知Ax)Cx)equiv 1mod;x^n)),倍增求Ax)Bx)equiv 1mod;x^{2n}))下文为了方便把x)省掉

AB-C)equiv 0mod;x^n))

A^2B-C)^2equiv 0mod;x^{2n}))

AB-C)^2equiv 0mod;x^{2n}))

AB^2-2ABC+AC^2equiv 0mod;x^{2n}))

B-2C+AC^2equiv 0mod;x^{2n}))

Bequiv 2C-AC^2mod;x^{2n}))

求导

A’x)),本质是求一点的斜率,定义A_i=Ax)[x^i])

A’x)=frac{sum A_ix+Delta)^i-x^i}{Delta})

Deltaightarrow0)A’x)=sum_{i>=1} iA_ix^{i-1})

复合函数求导

Gx)=FAx)),G’x)=F’Ax))A’x))

这个的意义是以x为自变量的导,而F’Ax)))的意义是以Ax)为自变量的导

所以ln中求的是G’x),牛顿迭代把A当作自变量求的是F’Ax))

积分

是不定积分,即求导的逆运算

int Ax)=sum_{i>=1} frac{1}{i}x^iB_{i-1})

积分后再求导常数项会消失

ln

Gx)=lnFx)))具体的意义是什么并不是很知道

复合函数求导:Gx)=FAx)))G’x)=F’Ax))A’x)),并且有ln’x)=frac{1}{x})

G’x)=frac{F’x)}{Fx)})

Gx)=int frac{F’x)}{Fx)})

牛顿迭代

对于普通的多项式Ax))求零点x,设上一次求得的是x0(不一定是零点)

x0,fx0)))处做切线,斜率为f’x0)),则根据简单三角函数有

frac{fx0)}{x0-x}=f’x0))

x=x0-frac{fx0)}{f’x0)})

牛顿迭代求的是近似解,但

exp

先不考虑精度,求Bx)=e^{Ax)}mod;x^n))

lnBx))=Ax)mod;x^n))

lnBx))-Ax)=0mod;x^n))

把B当作变量,A当作常数并不知道为什么可以这样

FBx))=lnBx))-Ax)),则F’Bx))=frac{1}{Bx)})(函数相加的导数=分别的导数和,A看作常数了所以是0)

把n/2,代牛顿迭代的式子,设Cx))是模x^n)下的解,求模x^{2n})下的解Bx))

Bx)=Cx)-Cx)lnCx))-Ax))mod;x^{2n}))

Bx)=Cx)1-lnCx))+Ax))mod;x^{2n}))

正确性证明

泰勒展开:

x即Bx)),x0即Cx)),相减之后前n项是0,平方之后前2n项是0,于是被模掉了(

所以只剩前两项了,根据牛顿迭代的定义

发现取的就是前两项,现在只剩前两项了,所以是精确解(

其实模完之后变成了一条直线,所以可以求解

n^2lnexp

模数不是ntt模数时可以用,也有学的必要

https://www.cnblogs.com/gmh77/p/13162153.html

exp

gx)=e^{fx)}),设g_n)gx)[x^n])f_n)fx)[x^n])

gx)=e^{fx)})

g’x)=e^{fx)}f’x))

因为f’x)=sum{i*x^{i-1}f_i})

所以xf’x)=sum_{i>=1}{i*x^if_i})

ng_n=sum_{i=0}^{n}{g_{n-i}if_i})

硬点g0=1,然后即可n^2求得exp

原因是fx)没有常数项(否则求不出来),所以e^x)展开后只有frac{x^0}{0!}=1)

ln

ng_n=sum_{i=0}^{n-1}{g_{n-i}if_i}+nf_n)

nf_n=ng_n-sum_{i=0}^{n-1}{g_{n-i}if_i})

快速幂

先ln,乘上系数后exp

lnexp来搞是线性卷积,dft再乘k是循环卷积

时间复杂度

上面的那一坨都是nlogn的,但是常数略大

调试方法&注意事项

一定要把不用的位置清空,相乘长度开到长度和

两个长度为N的多项式相乘时一定要把[N,2N-1]清空

调试小技♂巧:

快速幂->exp->ln->求导积分求逆->NTT,从后往前调试

要测小数据和中数据,多测几遍+改n和len看有没有变

ln再exp和exp再ln结果不变,可以利用这点来查错

可以在不用的位上加一些数看答案是否会改变

exp的组合意义:e^x=sum frac{x^i}{i!}),所以可以手玩判断exp是否写对

如果exp对了而ln+exp/exp+ln错了那就是数组没有清空

例题

6712.【2020.06.09省选模拟】题3

题解

推式子省略,变成快速幂

code

#include <bits/stdc++.h>
#define foa,b,c) for a=b; a<=c; a++)
#define fda,b,c) for a=b; a>=c; a--)
#define Cn,m) jc[n]*Jc[m]%998244353*Jc[n)-m)]%998244353)
#define mod 998244353
#define Mod 998244351
#define G 3
#define ll long long
#define file
using namespace std;

ll A2[524288],a[524288],b[524288],c[524288],w[524288],S,T,n,m,s;
int a2[20][524288],i,j,k,l,N,len;

ll qpowerll a,int b) {ll ans=1; while b) {if b&1) ans=ans*a%mod;a=a*a%mod;b>>=1;} return ans;}

//static ll a[maxn];
void dftll *a,int tp,int N,int len)
{
	int i,j,k,l,S=N,s1=2,s2=1;
	ll u,v,w,W;
	
	foi,0,N-1) A2[i]=a[a2[len][i]];
	memcpya,A2,N*8);
	
	foi,1,len)
	{
		W=tp==1)?qpowerG,mod-1)/s1):qpowerG,mod-1)-mod-1)/s1);S>>=1;
		foj,0,S-1)
		{
			w=1;
			fok,0,s2-1)
			{
				u=a[j*s1+k],v=a[j*s1+k+s2]*w;
				a[j*s1+k]=u+v)%mod;
				a[j*s1+k+s2]=u-v)%mod;
				w=w*W%mod;
			}
		}
		s1<<=1,s2<<=1;
	}
}
namespace Mul{ll a[524288],b[524288];}
void mulll *a,ll *b,ll *c,int N,int len)
{
	ll N2=qpowerN,Mod);
	int i,j,k,l;
	
	memcpyMul::a,a,N*8),memcpyMul::b,b,N*8);
	dftMul::a,1,N,len);dftMul::b,1,N,len);
	foi,0,N-1) c[i]=Mul::a[i]*Mul::b[i]%mod;
	dftc,-1,N,len);
	foi,0,N-1) c[i]=c[i]*N2%mod;
}
namespace Ny{ll a[524288],b[524288];}
void nyll *a,ll *b,int N,int len)
{
	int i,j,k,l;
	
	memsetb,0,N*8);
	if N==1) {b[0]=qpowera[0],Mod); return;}
	nya,b,N/2,len-1);
	
	memsetNy::a,0,N*16);
	mulb,b,Ny::a,N,len);
	memsetNy::b,0,N*16);
	foi,0,N-1) Ny::b[i]=a[i];
	mulNy::a,Ny::b,Ny::a,N*2,len+1);
	foi,0,N-1) b[i]=b[i]*2-Ny::a[i])%mod;
}
void daoll *a,int N)
{
	int i;
	foi,0,N-2) a[i]=a[i+1]*i+1)%mod;a[N-1]=0;
}
void jill *a,int N)
{
	int i;
	fdi,N-1,1) a[i]=a[i-1]*w[i]%mod;a[0]=0;
}
namespace LN{ll a[524288];}
void Lnll *a,int N,int len)
{
	int i;
	memsetLN::a,0,N*16);memcpyLN::a,a,N*8);
	nyLN::a,a,N,len);daoLN::a,N);
	mula,LN::a,a,N*2,len+1);
	jia,N);foi,N,N+N-1) a[i]=0;
}
namespace EXP{ll a[524288];}
void Expll *a,ll *b,int N,int len)
{
	int i,j,k,l;
	
	memsetb,0,N*8);
	if N==1) {b[0]=1;return;}
	Expa,b,N/2,len-1);
	
	memsetEXP::a,0,N*8);memcpyEXP::a,b,N*4);
	LnEXP::a,N,len);
	foi,0,N-1) EXP::a[i]=-EXP::a[i]+a[i])%mod;++EXP::a[0];foi,N,N+N-1) EXP::a[i]=0;
	mulb,EXP::a,b,N*2,len+1);foi,N,N+N-1) b[i]=0;
}
namespace Mi{ll a[524288];}
void mill *a,ll k,int N,int len)
{
	ll s=qpowera[0],k);
	int i;
	
	memcpyMi::a,a,N*8);
	LnMi::a,N,len);
	foi,0,N-1) Mi::a[i]=Mi::a[i]*k%mod)%mod;
	ExpMi::a,a,N,len);
	foi,0,N-1) a[i]=a[i]*s%mod;
}

void init)
{
	int I,s=1;
	w[1]=1;
	foi,2,200000) w[i]=mod-w[mod%i]*mod/i)%mod;
	
	foI,0,19)
	{
		foi,0,s-1)
		{
			j=i;k=0;
			fol,1,I) k=k*2+j&1),j>>=1;
			a2[I][i]=k;
		}
		s*=2;
	}
}
void work)
{
	s=1;foi,1,m-n+1) s=s*T-i+1)%mod)%mod*w[i]%mod,a[i-1]=s;
	s=1;foi,0,m-n) b[i]=s,s=s*S-n*T)-i+1)+1)%mod)%mod*w[i+1]%mod;
	mia,n,N,len);
	mula,b,a,N*2,len+1);
}

int main)
{
	freopen"sum.in","r",stdin);
	#ifdef file
	freopen"sum.out","w",stdout);
	#endif
	
	scanf"%lld%lld%lld%lld",&S,&T,&n,&m);len=ceillog2m-n+1));N=qpower2,len);
	init);
	
	work);
	printf"%lld
",a[m-n]+mod)%mod);
	
	fclosestdin);
	fclosestdout);
	return 0;
}

Published by

风君子

独自遨游何稽首 揭天掀地慰生平

发表回复

您的邮箱地址不会被公开。 必填项已用 * 标注