拉格朗日反演
Alioth_
2020-02-01 11:52:10
# 拉格朗日反演
对于两个至少常数项为0的多项式
$$
g(f(x))=x
$$
则互为复合逆
同时
$$
f(g(x))=x
$$
我们有
$$
[x^n]f(x)=\frac{1}{n} [x^{-1}] (\frac{1}{g(x)})^n
$$
扩域后
$$
[x^n]f(x)=\frac{1}{n} [x^{dn-1}] (\frac{x^d}{g(x)})^n
$$
其中$d$为$g$前缀为$0$项的个数 相当于把$g$左移$d$项使其可逆
这样 我们可以通过$g$求$f$
有其扩展形式
若$F(G(x))=x$ 则
$$
[x^n]H(F(x))=\frac1 n [x^{n-1}]H'(x)(\frac x {G(x)})^n
$$
# 例题
[BZOJ3684 大朋友和多叉树](https://www.lydsy.com/JudgeOnline/problem.php?id=3684)
设答案的生成函数为$f$
则
$$
f(x)=x+\sum\limits_{k\in D}f^k(x)
$$
其中 $+x$为单独一个点也算一棵树
移项 得
$$
f(x)-\sum\limits_{k\in D}f^k(x)=x
$$
令
$$
g(f(x))=f(x)-\sum\limits_{k\in D}f^k(x)=x
$$
则$g(x)$已知 可通过求逆、快速幂反演出$[x^n]f(x)$
```cpp
#include<bits/stdc++.h>
using namespace std;
const int maxn=4e5+7;
const int p=950009857;
const int inv2=475004929;
int g[maxn],n,r[maxn],f[maxn],invg[maxn],ans[maxn];
inline int read(){
int x=0;char ch=getchar();
while(ch<'0'||ch>'9')ch=getchar();
while(ch>='0'&&ch<='9'){x=(x*10ll%p+ch-'0')%p;ch=getchar();}
return x;
}
inline int poww(int a,int b){
int ans=1;
while(b){
if(b&1)ans=1ll*ans*a%p;
b>>=1;
a=1ll*a*a%p;
}
return ans%p;
}
inline int inv(int x){
return poww(x,p-2)%p;
}
int getrev(int deg){
int lim=1,l=0;
while(lim<=deg)++l,lim<<=1;
for(int i=0;i<=lim;++i)r[i]=r[i>>1]>>1|(i&1)<<(l-1);
return lim;
}
inline void NTT(int limit,int *A,int type){
int pr=7;
for(int i=0;i<limit;++i)if(r[i]>i)swap(A[i],A[r[i]]);
for(int mid=1;mid<limit;mid<<=1){
int wn=poww(type==-1?inv(pr):pr,(p-1)/(mid<<1));
for(int j=0,R=mid<<1;j<limit;j+=R){
int w=1;
for(int k=0;k<mid;++k,w=1ll*w*wn%p){
int x=A[j+k],y=1ll*w*A[j+k+mid]%p;
A[j+k]=(x+y)%p;
A[j+k+mid]=(x-y+p)%p;
}
}
}
if(type==-1){
int INV=inv(limit);
for(int i=0;i<limit;++i)A[i]=1ll*A[i]*INV%p;
}
}
void polyinv(int mod,int *A,int *B){
if(mod==1)return B[0]=inv(A[0]),void();
polyinv((mod+1)>>1,A,B);
int lim=getrev(mod<<1);
static int T[maxn];
for(int i=0;i<mod;++i)T[i]=A[i];
for(int i=mod;i<=lim;++i)T[i]=0;
NTT(lim,T,1);
NTT(lim,B,1);
for(int i=0;i<=lim;++i)B[i]=(2ll-1ll*T[i]*B[i]%p+p)%p*B[i]%p,T[i]=0;
NTT(lim,B,-1);
for(int i=mod;i<=lim;++i)B[i]=0;
}
void polysqrt(int mod,int *A,int *B){
if(mod==1)return B[0]=1,void();
polysqrt((mod+1)>>1,A,B);
int lim=getrev(mod<<1);
static int T[maxn],Inv[maxn];
for(int i=0;i<mod;++i)T[i]=A[i],Inv[i]=0;
for(int i=mod;i<=lim;++i)T[i]=0,Inv[i]=0;
polyinv(mod,B,Inv);
NTT(lim,B,1);NTT(lim,T,1);NTT(lim,Inv,1);
for(int i=0;i<=lim;++i)B[i]=(1ll*B[i]*B[i]%p+T[i]%p)%p*(1ll*Inv[i]*inv2%p)%p;
NTT(lim,B,-1);
for(int i=mod;i<=lim;++i)B[i]=0;
}
inline void polyder(int deg,int *A,int *B){
for(int i=1;i<=deg;++i)
B[i-1]=1ll*A[i]*i%p;
B[deg]=0;
}
inline void polyinte(int deg,int *A,int *B){
for(int i=deg+1;i>=1;--i)
B[i]=1ll*A[i-1]*inv(i)%p;
B[0]=0;
}
inline void polyln(int mod,int *A,int *B){//A[0] must be 1
static int D[maxn],Inv[maxn];
polyder(mod-1,A,D);
polyinv(mod,A,Inv);
int lim=getrev(mod<<1);
NTT(lim,Inv,1);NTT(lim,D,1);
for(int i=0;i<=lim;i++)B[i]=1ll*D[i]*Inv[i]%p,D[i]=Inv[i]=0;
NTT(lim,B,-1);
polyinte(mod-1,B,B);
}
inline void polyexp(int mod,int *A,int *B){//A[0] must be 0
if(mod==1)return B[0]=1,void();
polyexp((mod+1)>>1,A,B);
int lim=getrev(mod<<1);
static int Ln[maxn],T[maxn];
for(int i=0;i<mod;++i)T[i]=A[i],Ln[i]=0;
for(int i=mod;i<=lim;++i)T[i]=Ln[i]=0;
polyln(mod,B,Ln);
NTT(lim,B,1);NTT(lim,Ln,1);NTT(lim,T,1);
for(int i=0;i<=lim;++i)B[i]=((B[i]-1ll*B[i]*Ln[i]%p+1ll*T[i]*B[i]%p)%p+p)%p;
NTT(lim,B,-1);
for(int i=mod;i<=lim;++i)B[i]=0;
}
inline void polypow(int mod,int *A,int *B,int k){//A[0] must be 1
static int Ln[maxn];
polyln(mod,A,Ln);
for(int i=0;i<mod;++i)Ln[i]=1ll*Ln[i]*k%p;
for(int i=mod;i<maxn;++i)Ln[i]=0;
polyexp(mod,Ln,B);
}
int main(){
int n,m;
scanf("%d%d",&n,&m);
for(int i=1;i<=m;++i){
int k;scanf("%d",&k);
g[k-1]=p-1;
}
g[0]=1;
polyinv(n,g,invg);
polypow(n,invg,ans,n);
cout<<1ll*inv(n)*ans[n-1]%p;
}
```
[TJOI2015概率论](https://www.luogu.com.cn/problem/P3978)
设答案的生成函数为$f$ 则
$$
f(x)=xf^2(x)+1
$$
$$
xf(x)=x^2f^2(x)+x
$$
令
$$
T(x)=xf(x)
$$
则
$$
g(T(x))=T(x)-T^2(x)=x
$$
套公式得
$$
[x^n]f(x)=\binom{2n-2}{n-1}
$$
再除以卡特兰数就行了