多项式 ln exp

多项式 ln exp

单独写篇博客...

用的vector写板子,比数组慢一倍左右,安全性较高

多项式ln

[ln (F(x))=G(x) ]

[G'(x) = frac {dln (F(x))} {dF(x)} frac{d(F(X))}{dx}=frac{F'(x)}{F(x)} ]

[G(x)=intfrac{F'(x)}{F(x)} dx ]

#include<bits/stdc++.h>
using namespace std;
int read(){
	int x=0,pos=1;char ch=getchar();for(;!isdigit(ch);ch=getchar()) if(ch=='-') pos=0;
	for(;isdigit(ch);ch=getchar()) x=(x<<1)+(x<<3)+ch-'0';return pos?x:-x;
}
#define FOR(i,a,b) for(int i=a;i<=b;++i)
const int N = 200201;
#define ROF(i,a,b) for(int i=b;i>=a;--i)
int n,m;
int a[N*2];
typedef vector<int> poly;poly b[N];
const int mod = 998244353;
const int LSZ=20,SZ=1<<LSZ;
int ksm(int a,int b){
	int res=1;
	while(b){
		if(b&1) res=1ll*res*a%mod;
		a=1ll*a*a%mod,b>>=1;
	}return res;
}
int tr[SZ],omg[SZ][2];
void init(int n){
	int lim=0;while((1<<lim)<n) lim++;
	FOR(i,0,n-1) tr[i]=tr[i>>1]>>1|((i&1)<<(lim-1));
	omg[0][0]=omg[0][1]=1;
	omg[1][0]=ksm(3,(mod-1)/n);
	omg[1][1]=ksm(omg[1][0],mod-2);
	FOR(i,2,n) omg[i][0]=1ll*omg[i-1][0]*omg[1][0]%mod,omg[i][1]=1ll*omg[i-1][1]*omg[1][1]%mod;
}
void dft(poly &f,int n,int opt){
	//assert((int)f.size()<=n);
	f.resize(n,0);
	for(int i=0;i<n;++i) if(i<tr[i]) swap(f[i],f[tr[i]]);
	for(int l=2;l<=n;l<<=1){ int m=l/2;
		for(int g=0;g!=n;g+=l){
			for(int i=0;i<m;i++){
				int t=1ll*omg[n/l*i][opt]*(f[g+i+m])%mod;
				f[g+i+m]=(f[g+i]-t+mod)%mod;f[g+i]=(f[g+i]+t)%mod;
			}
		}
	}
	if(opt) for(int i=0,iv=ksm(n,mod-2);i<n;i++) f[i]=1ll*f[i]*iv%mod;
}
void pop(poly &a){
	if(!a.size()) return; 
	while(a[a.size()-1]==0) a.pop_back();
}
void dao(poly& a){
	for(int i=1;i<a.size();i++) a[i-1]=1ll*a[i]*i%mod;a.pop_back();
}
void jifen(poly &a){
	for(int i=a.size()-1;i>=0;--i) a[i]=1ll*a[i-1]*ksm(i,mod-2)%mod;a[0]=0;
}
poly operator+(const poly &a,const poly &b){
	poly res(a);
	res.resize(max(a.size(),b.size()),0);
	for(int i=0;i<b.size();++i) res[i]=(res[i]+b[i])%mod;
	return res;pop(res);
}
poly operator*(const poly &a,const poly &b){
	poly A=a,B=b;//pop(A),pop(B);
	int len=1;while(len<A.size()+B.size()) len*=2;init(len);
	poly C(len,0);
	dft(A,len,0),dft(B,len,0);
	FOR(i,0,len-1) C[i]=1ll*A[i]*B[i]%mod;
	dft(C,len,1);pop(C);return C;
} 
void inv(const poly &a,poly &b,int n){
	if(n==1) return b[0]=ksm(a[0],mod-2),void();
	inv(a,b,(n+1)/2);
	int len=1;while(len<n*2) len*=2;init(len);
	if(b.size()<len) b.resize(len,0);FOR(i,(n+1)/2,len-1) b[i]=0;
	poly c=a;if(c.size()<len) c.resize(len,0);FOR(i,n,len-1) c[i]=0;
	dft(c,len,0);dft(b,len,0);
	FOR(i,0,len-1) b[i]=(2ll-1ll*c[i]*b[i]%mod+mod)%mod*b[i]%mod;
	dft(b,len,1);
}
void getln(const poly &a,poly &ans,int n){
	poly b;b.resize(n,0);
	inv(a,b,n);pop(b);
	poly c=a;dao(c);
	ans=b*c;
	jifen(ans);pop(ans);
}
int main(){
	//freopen("P4725_8.in","r",stdin);
	//freopen("8.txt","w",stdout);
	n=read();poly a;poly ans;
	a.resize(n,0);ans.resize(n,0);
	for(int i=0;i<n;i++){
		a[i]=read();
	}
	getln(a,ans,n);
	for(int i=0;i<n;i++){
		printf("%lld ",ans[i]);
	}
	return 0;
}

多项式exp

[G(x)=e^{F(x)} ]

[ln(G(x))=F(x) ]

[H(G(x))=ln(G(x))-F(x)=0 ]

由牛顿迭代,

[G(x)=G_0(x)-frac{ln(G_0(x))-F(x)}{G_0(x)^{-1}} ]

[=G_0(x)(1-ln(G_0(x))+F(x)) ]

边界,(G(x)[x^0]=e^{F(x)[X^0]})

 #include<bits/stdc++.h>
using namespace std;
int read(){
	int x=0,pos=1;char ch=getchar();for(;!isdigit(ch);ch=getchar()) if(ch=='-') pos=0;
	for(;isdigit(ch);ch=getchar()) x=(x<<1)+(x<<3)+ch-'0';return pos?x:-x;
}
#define FOR(i,a,b) for(int i=a;i<=b;++i)
const int N = 200201;
#define ROF(i,a,b) for(int i=b;i>=a;--i)
int n,m;
int a[N*2];
typedef vector<int> poly;poly b[N];
const int mod = 998244353;
const int LSZ=20,SZ=1<<LSZ;
int ksm(int a,int b){
	int res=1;
	while(b){
		if(b&1) res=1ll*res*a%mod;
		a=1ll*a*a%mod,b>>=1;
	}return res;
}
int tr[SZ],omg[SZ][2];
void init(int n){
	int lim=0;while((1<<lim)<n) lim++;
	FOR(i,0,n-1) tr[i]=tr[i>>1]>>1|((i&1)<<(lim-1));
	omg[0][0]=omg[0][1]=1;
	omg[1][0]=ksm(3,(mod-1)/n);
	omg[1][1]=ksm(omg[1][0],mod-2);
	FOR(i,2,n) omg[i][0]=1ll*omg[i-1][0]*omg[1][0]%mod,omg[i][1]=1ll*omg[i-1][1]*omg[1][1]%mod;
}
void dft(poly &f,int n,int opt){
	//assert((int)f.size()<=n);
	f.resize(n,0);
	for(int i=0;i<n;++i) if(i<tr[i]) swap(f[i],f[tr[i]]);
	for(int l=2;l<=n;l<<=1){ int m=l/2;
		for(int g=0;g!=n;g+=l){
			for(int i=0;i<m;i++){
				int t=1ll*omg[n/l*i][opt]*(f[g+i+m])%mod;
				f[g+i+m]=(f[g+i]-t+mod)%mod;f[g+i]=(f[g+i]+t)%mod;
			}
		}
	}
	if(opt) for(int i=0,iv=ksm(n,mod-2);i<n;i++) f[i]=1ll*f[i]*iv%mod;
}
void pop(poly &a){
	if(!a.size()) return; 
	while(a[a.size()-1]==0) a.pop_back();
}
void dao(poly& a){
	for(int i=1;i<a.size();i++) a[i-1]=1ll*a[i]*i%mod;a.pop_back();
}
void jifen(poly &a){
	if(!a.size()) return; for(int i=a.size()-1;i>=0;--i) a[i]=1ll*a[i-1]*ksm(i,mod-2)%mod;a[0]=0;
}
poly operator+(const poly &a,const poly &b){
	poly res(a);
	res.resize(max(a.size(),b.size()),0);
	for(int i=0;i<b.size();++i) res[i]=(res[i]+b[i])%mod;
	return res;pop(res);
}
poly operator*(const poly &a,const poly &b){
	poly A=a,B=b;//pop(A),pop(B);
	int len=1;while(len<A.size()+B.size()) len*=2;init(len);
	poly C(len,0);
	dft(A,len,0),dft(B,len,0);
	FOR(i,0,len-1) C[i]=1ll*A[i]*B[i]%mod;
	dft(C,len,1);pop(C);return C;
} 
void inv(const poly &a,poly &b,int n){
	if(n==1) return b[0]=ksm(a[0],mod-2),void();
	inv(a,b,(n+1)/2);
	int len=1;while(len<n*2) len*=2;init(len);
	if(b.size()<len) b.resize(len,0);FOR(i,(n+1)/2,len-1) b[i]=0;
	poly c=a;if(c.size()<len) c.resize(len,0);FOR(i,n,len-1) c[i]=0;
	dft(c,len,0);dft(b,len,0);
	FOR(i,0,len-1) b[i]=(2ll-1ll*c[i]*b[i]%mod+mod)%mod*b[i]%mod;
	dft(b,len,1);
}
void getln(const poly &a,poly &ans,int n){
	poly b;b.resize(n,0);
	inv(a,b,n);pop(b);
	poly c=a;dao(c);
	ans=b*c;
	jifen(ans);
}
void getexp(const poly &a,poly &ans,int n){
	if(n==1) return ans[0]=1,void();
	getexp(a,ans,(n+1)/2);
	ans.resize(n,0);
	poly tmp;getln(ans,tmp,n);tmp.resize(n,0);
	tmp[0]=(1ll+a[0]-tmp[0]+mod)%mod;
	FOR(i,1,n-1) tmp[i]=(0ll+a[i]-tmp[i]+mod)%mod;
	ans=ans*tmp;
	ans.resize(n,0);
} 
int main(){
	n=read();poly a;poly ans;
	a.resize(n,0);ans.resize(n,0);
	for(int i=0;i<n;i++){
		a[i]=read();
	}
	getexp(a,ans,n);
	for(int i=0;i<n;i++){
		printf("%d ",ans[i]);
	}
	return 0;
}

练习

以后有其他的这类题也会放在这里。

集训队作业2013 城市规划

考虑ln和exp的本质。

[G(F(x))=e^{F(x)} =sum_{ige 0} frac{F^i(x)}{i!} ]

也就是说,G(F(x))是i个F(x)的EGF乘起来的和,也就是把G的方案分成了i个F的方案。

本题中,G表示无向图个数的EGF,F表示无向连通图的EGF,显然满足上式。G是很好求的,F直接ln即可。

那么要把i个F组合成G,相当于反过来,直接exp。这两个技巧可用于许多计数题。

 #include<bits/stdc++.h>
using namespace std;
int read(){
	int x=0,pos=1;char ch=getchar();for(;!isdigit(ch);ch=getchar()) if(ch=='-') pos=0;
	for(;isdigit(ch);ch=getchar()) x=(x<<1)+(x<<3)+ch-'0';return pos?x:-x;
}
#define FOR(i,a,b) for(int i=a;i<=b;++i)
const int N = 200201;
#define ROF(i,a,b) for(int i=a;i>=b;--i)
int n,m;
typedef vector<int> poly;
const int mod = 1004535809;
const int LSZ=19,SZ=1<<LSZ;
int ksm(int a,int b){
	int res=1;
	while(b){
		if(b&1) res=1ll*res*a%mod;
		a=1ll*a*a%mod,b>>=1;
	}return res;
}
int tr[SZ],omg[SZ][2];
void init(int n){
	int lim=0;while((1<<lim)<n) lim++;
	FOR(i,0,n-1) tr[i]=tr[i>>1]>>1|((i&1)<<(lim-1));
	omg[0][0]=omg[0][1]=1;
	omg[1][0]=ksm(3,(mod-1)/n);
	omg[1][1]=ksm(omg[1][0],mod-2);
	FOR(i,2,n) omg[i][0]=1ll*omg[i-1][0]*omg[1][0]%mod,omg[i][1]=1ll*omg[i-1][1]*omg[1][1]%mod;
}
void dft(poly &f,int n,int opt){
	//assert((int)f.size()<=n);
	f.resize(n,0);
	for(int i=0;i<n;++i) if(i<tr[i]) swap(f[i],f[tr[i]]);
	for(int l=2;l<=n;l<<=1){ int m=l/2;
		for(int g=0;g!=n;g+=l){
			for(int i=0;i<m;i++){
				int t=1ll*omg[n/l*i][opt]*(f[g+i+m])%mod;
				f[g+i+m]=(f[g+i]-t+mod)%mod;f[g+i]=(f[g+i]+t)%mod;
			}
		}
	}
	if(opt) for(int i=0,iv=ksm(n,mod-2);i<n;i++) f[i]=1ll*f[i]*iv%mod;
}
void dao(poly& a){
	for(int i=1;i<a.size();i++) a[i-1]=1ll*a[i]*i%mod;a.pop_back();
}
void jifen(poly &a){
	if(!a.size()) return; for(int i=a.size()-1;i>=0;--i) a[i]=1ll*a[i-1]*ksm(i,mod-2)%mod;a[0]=0;
}
poly operator+(const poly &a,const poly &b){
	poly res(a);
	res.resize(max(a.size(),b.size()),0);
	for(int i=0;i<b.size();++i) res[i]=(res[i]+b[i])%mod;
	return res;
}
poly operator*(const poly &a,const poly &b){
	poly A=a,B=b;//pop(A),pop(B);
	int len=1;while(len<A.size()+B.size()) len*=2;init(len);
	poly C(len,0);
	dft(A,len,0),dft(B,len,0);
	FOR(i,0,len-1) C[i]=1ll*A[i]*B[i]%mod;
	dft(C,len,1);return C;
} 
void inv(const poly &a,poly &b,int n){
	if(n==1) return b[0]=ksm(a[0],mod-2),void();
	inv(a,b,(n+1)/2);
	int len=1;while(len<n*2) len*=2;init(len);
	if(b.size()<len) b.resize(len,0);FOR(i,(n+1)/2,len-1) b[i]=0;
	poly c=a;if(c.size()<len) c.resize(len,0);FOR(i,n,len-1) c[i]=0;
	dft(c,len,0);dft(b,len,0);
	FOR(i,0,len-1) b[i]=(2ll-1ll*c[i]*b[i]%mod+mod)%mod*b[i]%mod;
	dft(b,len,1);
}
void getln(const poly &a,poly &ans,int n){
	poly b;b.resize(n,0);
	inv(a,b,n);
	poly c=a;dao(c);
	ans=b*c;
	jifen(ans);
}
void getexp(const poly &a,poly &ans,int n){
	if(n==1) return ans[0]=1,void();
	getexp(a,ans,(n+1)/2);
	ans.resize(n,0);
	poly tmp;getln(ans,tmp,n);tmp.resize(n,0);
	tmp[0]=(1ll+a[0]-tmp[0]+mod)%mod;
	FOR(i,1,n-1) tmp[i]=(0ll+a[i]-tmp[i]+mod)%mod;
	ans=ans*tmp;
	ans.resize(n,0);
} 
int f[N],g[N];
void init1(){
	f[0]=f[1]=1;
	FOR(i,2,n) f[i]=1ll*f[i-1]*i%mod;
	g[n]=ksm(f[n],mod-2);ROF(i,n-1,0) g[i]=1ll*g[i+1]*(i+1)%mod; 
}
int main(){
	n=read();poly a;a.resize(n+1,0);init1();
	int iv2=ksm(2,mod-2);
	a[0]=a[1]=1;
	for(int i=2;i<=n;i++){
		a[i]=1ll*ksm(2,1ll*i*(i-1)/2%(mod-1))*g[i]%mod;
	}
	poly ans;getln(a,ans,n+1);
	printf("%lld",1ll*ans[n]*f[n]%mod);
	return 0;
}
原文地址:https://www.cnblogs.com/lcyfrog/p/12801206.html