[学习笔记]拉格朗日反演

拉格朗日反演可以互推一对多项式复合逆，拉反通常起到将函数自、因变量交换，简化柿子的作用。

定义

• 复合逆：$F(G(x))\equiv G(F(x))\equiv x\pmod {x^n}$
  则称$F$和$G$互为复合逆，记作$G(x)=F^{(-1)}(x)$

• 引理：$[x^{-1}]F^k(x)F'(x)=[k=-1]$

证明

当$k\ne-1$，有$F'(x)F(x)^k=(\frac{1}{k+1}F(x){k+1})'$,是整式显然无$[x^{-1}]$项；

当$k=-1$,有$[x^{-1}]\frac{F'(x)}{F(x)}=[x^0]\frac{F'(x)}{F(x)/x}$ 上下的常数项均为$[x^1]F(x)$所以值为$1$。

• 拉格朗日反演：若$G(F(x))=x$，且满足常数项为$0$,一次项可逆，则$[x^n]G^k(x)=\dfrac{k}{n}[x^{-k}]F^{-n}(x)=\dfrac{k}{n}[x^{n-k}](\dfrac{x}{F(x)})^n$

通常使用最右边形式，因为要常数项可逆。

证明

$G^k(F(x))=x^k$

两侧求导得: $(G^k)'(F(x))F'(x)=(\sum\limits_{i}i([x^i]G^k)F^{i-1})F'(x)=kx^{k-1}$

同除$F^n(x)$得：$\sum\limits_{i}i([x^i]G^k)F^{i-n-1}F(x)=\dfrac{kx^{k-1}}{F^n(x)}=k[x^{-k}]F^{-n}(x)$
根据引理1: $n[x^n]G^k(x)=k[x^{-k}]F^{-n}(x)$

• 推论(常用)：$[x^n]G(x)=\dfrac{1}{n}[x^{-1}]F^{-n}(x)=\dfrac{1}{n}[x^{n-1}](\dfrac{x}{F(x)})^n$

• 扩展拉格朗日反演：给定$H(x)\in F((x)),G(F(x))=x$，则$[x^n]H(G(x))=\dfrac{1}{n}[x^{-1}]\dfrac{H'(x)}{F^n(x)}=\dfrac{1}{n}[x^{n-1}]H'(x)(\dfrac{x}{F(x)})^n$
  大概意义是复合逆的复合函数。

证明

设$A(x)=H(G(x))$

$H(x)=H(G(F(x)))=A(F(x))=\sum\limits_{i\ge0}([x^i]A(x))F^i(x)$

求导可得：$H'=\sum\limits_{i\ge1}([x^i]A(x))i\cdot F^{i-1}F'$

同除$F^n$得：$H'F^{-n}=\sum\limits_{i\ge1}([x^i]A(x))i\cdot F^{i-1-n}F'$

根据引理1：$[x^{-1}]H'F^{-n}=n[x^n]A(x)$

当 $n=0$ 时，分母为 $0$ 就需要用到下面的另类方法

• 另类拉格朗日反演：
  $[x^n]G^k(x)=[x^{-k-1}]\dfrac{F'(x)}{F^{n+1}(x)}=[x^{n-k}]F'(x)(\dfrac{x}{F(x)})^{n+1}$

证明

$G^k(F(x))=x^k$

展开可得：$\sum\limits_{i}([x^i]G^k)F^i=x^k$

两侧同除$x^n$：$\sum\limits_{i}([x^i]G^k)F'F^{i-n-1}=x^kF'F^{-n-1}$

根据引理1:$[x^n]G^k=[x^{-1}]x^kF'F^{-n-1}=[x^{-k-1}]F'F^{-n-1}$

• 扩展另类拉格朗日反演：$[x^n]H(G(x))=[x^n]H(x)F'(x)(\dfrac{x}{F(x)})^{n+1}$

应用

【BZOJ3684】大朋友和多叉树

• 题意：
  给你$s$和$m$,在给你$m$个数的集合$D$。问你叶子数为$s$且每个度数大于$1$的点的儿子个数在$D$中的方案数。
• 思路：
  $F(x)$为叶子数为$x$的方案数。
  $F(x)=x+\sum\limits_{i\in D}F^i(x)$
  发现$F(x)$项比$x$项复杂。考虑拉反（代入$x=F^{(-1)}(x)$）。
  $x=F^{(-1)}(x)+\sum\limits_{i\in D}x^i$
  $F^{(-1)}(x)=x-\sum\limits_{i\in D}x^i$
  由拉反得：$Ans=[x^s]F(x)=[x^{s-1}](\dfrac{x}{F^{(-1)}(x)})^s$
  因此贺一下多项式幂函数+多项式求逆的板子就可以了（板子太烂，心累ing~）。
• code:

贺的

#include<bits/stdc++.h>
using namespace std;
typedef double db;
typedef long long ll;
const int N=1e6+5;
const ll mod=950009857;
ll inv[N],inv7,gen[2][N];

ll ksm(ll a,ll b) {ll mul=1;for(;b;b>>=1,a=a*a%mod)if(b&1)mul=mul*a%mod;return mul;}
int rev[N],L,up;
void NTT(ll *a,int op) {
	for(int i=0;i<up;i++) {
		if(rev[i]>i)swap(a[i],a[rev[i]]);
	}
	for(int mid=1;mid<up;mid<<=1) {
		int len=mid<<1;ll w1=gen[op][len];
 		for(int l=0;l<up;l+=len) {
 			ll W=1;
			for(int i=0;i<mid;i++,W=W*w1%mod) {
				int p=l+i,q=p+mid;
				ll x=a[p],y=W*a[q];
				a[p]=(x+y)%mod;a[q]=(x-y)%mod;
			}
		}
	}
}

void gt_up(int len) {
	up=1,L=0;
	while(up<len) {up<<=1,L++;}
	for(int i=1;i<up;i++) {rev[i]=(rev[i>>1]>>1)|((i&1)<<(L-1));}
}

ll a[N];
void poly_inv(int dep,ll *f,ll *b) {	// b = f mod(x^dep) 
	if(dep==1) {b[0]=ksm(f[0],mod-2);return;}
	poly_inv((dep+1)>>1,f,b);
	gt_up(dep<<1);
	for(int i=0;i<dep;i++)a[i]=f[i];for(int i=dep;i<up;i++)b[i]=a[i]=0;
	NTT(b,0);NTT(a,0);
	for(int i=0;i<up;i++) {a[i]=(2-a[i]*b[i]%mod)*b[i]%mod;}
	NTT(a,1);
	ll i_up=ksm(up,mod-2);
	for(int i=0;i<dep;i++) {b[i]=a[i]*i_up%mod;} for(int i=dep;i<up;i++)b[i]=0;
}

void poly_dao(int n,ll *f,ll *g) {
	for(int i=1;i<n;i++) g[i-1]=f[i]*i%mod;g[n-1]=0;
}
void poly_jf(int n,ll *f,ll *g) {
	for(int i=1;i<n;i++) g[i]=f[i-1]*inv[i]%mod;g[0]=0;
}

ll d[N],g[N];
void poly_ln(int n,ll *f) {
	poly_dao(n,f,d);
	for(int i=0;i<n;i++)g[i]=0;poly_inv(n,f,g);		//****b must be 1,0....0 at first 
	gt_up(n<<1);
	for(int i=n;i<up;i++)g[i]=d[i]=0;
	NTT(d,0),NTT(g,0);
	for(int i=0;i<up;i++)g[i]=d[i]*g[i]%mod;
	NTT(g,1);
	for(int i=0,i_up=ksm(up,mod-2);i<up;i++)g[i]=(i<n)?(g[i]*i_up%mod):0;
	poly_jf(n,g,f);
}
ll c[N];
void poly_exp(int dep,ll *f,ll *b) {
	if(dep==1) {b[0]=1;return;}
	poly_exp((dep+1)>>1,f,b);
	for(int i=0;i<dep;i++)c[i]=b[i];
	poly_ln(dep,c);
	gt_up(dep<<1);
	for(int i=0;i<dep;i++)c[i]=(f[i]-c[i])%mod;for(int i=dep;i<up;i++)b[i]=c[i]=0;
	c[0]++;

	NTT(b,0);NTT(c,0);
	for(int i=0;i<up;i++)b[i]=c[i]*b[i]%mod;
	NTT(b,1);
	for(int i=0,i_up=ksm(up,mod-2);i<up;i++) {
		b[i]=(i<dep)?(b[i]*i_up%mod):0;
	}
}

void poly_pow(int n,ll *f,ll *g,ll k) {
	poly_ln(n,f);
	for(int i=0;i<n;i++) {f[i]=f[i]*k%mod;}
	poly_exp(n,f,g);
}

void gt_gen(int len) {
	inv7=ksm(7,mod-2);
	gt_up(len);
	up=min(up,1<<21);
	gen[0][up]=ksm(7,(mod-1)/up);gen[1][up]=ksm(inv7,(mod-1)/up);
	for(int i=up;i;i>>=1) {
		gen[0][i>>1]=gen[0][i]*gen[0][i]%mod;
		gen[1][i>>1]=gen[1][i]*gen[1][i]%mod;
	}
	inv[1]=1;for(int i=2;i<up;i++)inv[i]=(mod-mod/i)*inv[mod%i]%mod;
}

ll f[N],oo[N],ans[N];
int main() {
	int s,m;
	scanf("%d%d",&s,&m);
	for(int i=1;i<=m;i++) {int k;scanf("%d",&k);f[k-1]=-1;}
	f[0]=1;
	gt_gen(s<<1);
	poly_inv(s,f,oo);
	poly_pow(s,oo,ans,s);
	printf("%lld\n",(ans[s-1]*ksm(s,mod-2)%mod+mod)%mod);
	return 0;
}

树的数量

传送门
转载