多項式小結(求逆、求導、積分、ln、牛頓迭代、exp)


求逆

\(A(x)B(x)\equiv 1(mod\;x^n)\),下文為了方便表述把n/2

已知\(A(x)C(x)\equiv 1(mod\;x^n)\),倍增求\(A(x)B(x)\equiv 1(mod\;x^{2n})\)下文為了方便把(x)省掉

\(A(B-C)\equiv 0(mod\;x^n)\)

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

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

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

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

\(B\equiv 2C-AC^2(mod\;x^{2n})\)

求導

\(A'(x)\),本質是求一點的斜率,定義\(A_i=A(x)[x^i]\)

\(A'(x)=\frac{\sum A_i((x+\Delta)^i-x^i}{\Delta}\)

\(\Delta\rightarrow0\)\(A'(x)=\sum_{i>=1} iA_ix^{i-1}\)

復合函數求導

\(G(x)=F(A(x)),G'(x)=F'(A(x))A'(x)\)

這個的意義是以x為自變量的導,而\(F'(A(x))\)的意義是以A(x)為自變量的導

所以ln中求的是G'(x),牛頓迭代把A當作自變量求的是F'(A(x))

積分

是不定積分,即求導的逆運算

\(\int A(x)=\sum_{i>=1} \frac{1}{i}x^iB_{i-1}\)

積分后再求導常數項會消失

ln

\(G(x)=ln(F(x))\)具體的意義是什么並不是很知道

復合函數求導:\(G(x)=F(A(x))\)\(G'(x)=F'(A(x))A'(x)\),並且有\(ln'(x)=\frac{1}{x}\)

\(G'(x)=\frac{F'(x)}{F(x)}\)

\(G(x)=\int \frac{F'(x)}{F(x)}\)

牛頓迭代

對於普通的多項式\(A(x)\)求零點x,設上一次求得的是x0(不一定是零點)

\((x0,f(x0))\)處做切線,斜率為\(f'(x0)\),則根據簡單三角函數有

\(\frac{f(x0)}{x0-x}=f'(x0)\)

\(x=x0-\frac{f(x0)}{f'(x0)}\)

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

exp

先不考慮精度,求\(B(x)=e^{A(x)}(mod\;x^n)\)

\(ln(B(x))=A(x)(mod\;x^n)\)

\(ln(B(x))-A(x)=0(mod\;x^n)\)

把B當作變量,A當作常數並不知道為什么可以這樣

\(F(B(x))=ln(B(x))-A(x)\),則\(F'(B(x))=\frac{1}{B(x)}\)(函數相加的導數=分別的導數和,A看作常數了所以是0)

把n/2,代牛頓迭代的式子,設\(C(x)\)是模\(x^n\)下的解,求模\(x^{2n}\)下的解\(B(x)\)

\(B(x)=C(x)-C(x)(ln(C(x))-A(x))(mod\;x^{2n})\)

\(B(x)=C(x)(1-ln(C(x))+A(x))(mod\;x^{2n})\)

正確性證明

泰勒展開:

x即\(B(x)\),x0即\(C(x)\),相減之后前n項是0,平方之后前2n項是0,於是被模掉了(

所以只剩前兩項了,根據牛頓迭代的定義

發現取的就是前兩項,現在只剩前兩項了,所以是精確解(

其實模完之后變成了一條直線,所以可以求解

n^2lnexp

模數不是ntt模數時可以用,也有學的必要

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

exp

\(g(x)=e^{f(x)}\),設\(g_n\)\(g(x)[x^n]\)\(f_n\)\(f(x)[x^n]\)

\(g(x)=e^{f(x)}\)

\(g'(x)=e^{f(x)}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

原因是f(x)沒有常數項(否則求不出來),所以\(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 fo(a,b,c) for (a=b; a<=c; a++)
#define fd(a,b,c) for (a=b; a>=c; a--)
#define C(n,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 qpower(ll 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 dft(ll *a,int tp,int N,int len)
{
	int i,j,k,l,S=N,s1=2,s2=1;
	ll u,v,w,W;
	
	fo(i,0,N-1) A2[i]=a[a2[len][i]];
	memcpy(a,A2,N*8);
	
	fo(i,1,len)
	{
		W=(tp==1)?qpower(G,(mod-1)/s1):qpower(G,(mod-1)-(mod-1)/s1);S>>=1;
		fo(j,0,S-1)
		{
			w=1;
			fo(k,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 mul(ll *a,ll *b,ll *c,int N,int len)
{
	ll N2=qpower(N,Mod);
	int i,j,k,l;
	
	memcpy(Mul::a,a,N*8),memcpy(Mul::b,b,N*8);
	dft(Mul::a,1,N,len);dft(Mul::b,1,N,len);
	fo(i,0,N-1) c[i]=Mul::a[i]*Mul::b[i]%mod;
	dft(c,-1,N,len);
	fo(i,0,N-1) c[i]=c[i]*N2%mod;
}
namespace Ny{ll a[524288],b[524288];}
void ny(ll *a,ll *b,int N,int len)
{
	int i,j,k,l;
	
	memset(b,0,N*8);
	if (N==1) {b[0]=qpower(a[0],Mod); return;}
	ny(a,b,N/2,len-1);
	
	memset(Ny::a,0,N*16);
	mul(b,b,Ny::a,N,len);
	memset(Ny::b,0,N*16);
	fo(i,0,N-1) Ny::b[i]=a[i];
	mul(Ny::a,Ny::b,Ny::a,N*2,len+1);
	fo(i,0,N-1) b[i]=(b[i]*2-Ny::a[i])%mod;
}
void dao(ll *a,int N)
{
	int i;
	fo(i,0,N-2) a[i]=a[i+1]*(i+1)%mod;a[N-1]=0;
}
void ji(ll *a,int N)
{
	int i;
	fd(i,N-1,1) a[i]=a[i-1]*w[i]%mod;a[0]=0;
}
namespace LN{ll a[524288];}
void Ln(ll *a,int N,int len)
{
	int i;
	memset(LN::a,0,N*16);memcpy(LN::a,a,N*8);
	ny(LN::a,a,N,len);dao(LN::a,N);
	mul(a,LN::a,a,N*2,len+1);
	ji(a,N);fo(i,N,N+N-1) a[i]=0;
}
namespace EXP{ll a[524288];}
void Exp(ll *a,ll *b,int N,int len)
{
	int i,j,k,l;
	
	memset(b,0,N*8);
	if (N==1) {b[0]=1;return;}
	Exp(a,b,N/2,len-1);
	
	memset(EXP::a,0,N*8);memcpy(EXP::a,b,N*4);
	Ln(EXP::a,N,len);
	fo(i,0,N-1) EXP::a[i]=(-EXP::a[i]+a[i])%mod;++EXP::a[0];fo(i,N,N+N-1) EXP::a[i]=0;
	mul(b,EXP::a,b,N*2,len+1);fo(i,N,N+N-1) b[i]=0;
}
namespace Mi{ll a[524288];}
void mi(ll *a,ll k,int N,int len)
{
	ll s=qpower(a[0],k);
	int i;
	
	memcpy(Mi::a,a,N*8);
	Ln(Mi::a,N,len);
	fo(i,0,N-1) Mi::a[i]=Mi::a[i]*(k%mod)%mod;
	Exp(Mi::a,a,N,len);
	fo(i,0,N-1) a[i]=a[i]*s%mod;
}

void init()
{
	int I,s=1;
	w[1]=1;
	fo(i,2,200000) w[i]=mod-w[mod%i]*(mod/i)%mod;
	
	fo(I,0,19)
	{
		fo(i,0,s-1)
		{
			j=i;k=0;
			fo(l,1,I) k=k*2+(j&1),j>>=1;
			a2[I][i]=k;
		}
		s*=2;
	}
}
void work()
{
	s=1;fo(i,1,m-n+1) s=s*((T-i+1)%mod)%mod*w[i]%mod,a[i-1]=s;
	s=1;fo(i,0,m-n) b[i]=s,s=s*(((S-n*T)-(i+1)+1)%mod)%mod*w[i+1]%mod;
	mi(a,n,N,len);
	mul(a,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=ceil(log2(m-n+1));N=qpower(2,len);
	init();
	
	work();
	printf("%lld\n",(a[m-n]+mod)%mod);
	
	fclose(stdin);
	fclose(stdout);
	return 0;
}


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM