【LOJ6363】「地底蔷薇」【点双】【指数型生成函数】【扩展拉格朗日反演】【多项式幂函数】

传送门

题意：给定$n$和集合$S$,求出$n$个点的「所有极大点双连通分量的大小都在$S$ 内」的不同简单无向连通图的个数 模 $998244353$。

$n,{\sum }_{i\in S}i\le 10^5$

道理我都懂，可为啥我百度搜地灵殿ex终符一半都是cookie☆同人OI题

思路挺清奇的

首先注意题目中的点双定义是不存在割点，也就是说只有两个点的连通图也是点双

同时一个点可能属于多个点双

显然要设出答案的EGF

因为无向图关系很复杂，所以我们要强行整点东西上去

设$[x^n]F(x)$表示 $n$个点 带标号 有根 满足题目条件 的方案数 的EGF

有根是指：每张图钦定一个根，如果两张图边完全相同但根不同也视为不同的图

其实就是乘一个$n$

考虑怎么表示出$F(x)$

对于一张有根图，我们把根删掉，这样会出现多个连通块

对于每个连通块，会发现都会有若干点和根处于同一个点双内

但是不会有不同连通块的点处于同一点双，因为根结点删去后它们分开了，与定义矛盾

考虑一个连通块的情况

先考虑根结点所在的点双大小，发现同一组点组成的点双也有多种连边方案

设$a_n$表示$n+1$个点构成的带标号点双个数

再设出满足题目要求的点双的EGF

$$A(x)=\sum _{i+1\in S}\frac{a_i}{i!}{x}^i$$

因为没有跨连通块，所以这个连通块的点双一定在$S$内

然后将点双中的所有边删掉，这样点双中的点两两间都不会连通 因为如果连通，由于内部边都断完了，只能走外面的点，而经过的点都会在这个点双内

也就是说，这个大连通块被分成了若干小连通块，且连通块个数等于点双的大小（不含根）,小连通块的根是原来的点双中的点

这样可以写出一个大连通块的EGF

$$\sum _{i=1}^{\infty }\frac{a_i}{i!}{F}^i(x)$$

也就是$A(F(x))$

注意：这个过程并没有选出点双中有哪些点，只是把连通块中的点分组并将每组的根按点双的方式连接起来

然后$F(x)$可以表示成任意个数连通块组合再加上根

也就是

$$F(x)=x{e}^{A(F(x))}$$

假设我们已经求出了$A(x)$，我们尝试算出$F(x)$

$$\frac{F(x)}{e^{A(F(x))}}=x$$

设

$$f(F(x))=\frac{F(x)}{e^{A(F(x))}}=x$$

即$f$和$F$互为复合逆

而

$$f(x)=\frac{x}{e^{A(x)}}$$

显然$[x^0]A(x)=0$,所以可以算出$f(x)$，进而用拉格朗日反演算出$[x^n]F(x)$

---

现在的问题是如何求$A(x)$

直接求不好求，考虑构造一个相似的问题用其他方法算一次，再把$A(x)$反推回来

那考虑减少限制 发现减掉在$S$中的限制，也就是任意连通图，我们也可以用这个思路计算

并且这是一道原题

城市规划

直接取对数就可以轻松算出答案，岂不美哉？

好，我们仿（fu）照（zhi）这道题，算出任意无向连通图的EGF $G(x)$,并乘一个$n$上个根

现在$G(x)$ 已知了

设

$$B(x)=\sum _{i=1}^{\infty }\frac{a_i}{i!}{x}^i$$

相同的思路，得到

$$G(x)=x{e}^{B(G(x))}$$

$$\frac{G(x)}{e^{B(G(x))}}=x$$

$$g(G(x))=\frac{G(x)}{e^{B(G(x))}}=x$$

$$g(x)=\frac{x}{e^{B(x)}}$$

变一下可以得到

$$B(x)=\ln (\frac{x}{g(x)})$$

用拉格朗日反演算出$g(x)$……

拉个鬼啊，一次拉格朗日反演就是$O(nlogn)$的，你要算出$g(x)$就要跑$n$次，再怎么也有$O(n^2)$了吧

但我们注意到$A(x)$有值的地方很少，而且下标和不超过$1e5$，所以如果我们直接拉出$A(x)$就没有复杂度问题啦

好，整理一下已知条件

$$g(G(x))=\frac{G(x)}{e^{B(G(x))}}=x$$

$$B(x)=\ln (\frac{x}{g(x)})$$

闲着没事把上面代到下面

（注：这里为了统一代码改一下）

$$B(G(g(x)))=\ln (\frac{x}{g(x)})$$

好多括号啊，拆掉里面那个

$$B(G(x))=\ln (\frac{G(x)}{x})$$

---

停，我说一句

有个叫扩展拉格朗日反演的东西,长这样:

$$[x^n]H(G(x))=\frac{1}{n}[x^{-1}]\frac{H^{\prime }(x)}{F^n(x)}$$

证明的话在$G(F(x))=x$两边同时复合一个$H(x)$,并令$H(G(x))=T(x)$，后面就一模一样了

---

好继续

我们强行构造上面的结论，令

$$H(g(x))=B(x)$$

这样

$$[x^n]B(x)=[x^n]H(g(x))=\frac{1}{n}[x^{-1}]\frac{H^{\prime }(x)}{G^n(x)}$$

就可以把$B(x)$算出来

现在只需要求出$H(x)$

由于

$$H(g(x))=B(x)$$

得到

$$H(x)=H(g(G(x)))=B(G(x))=\ln (\frac{G(x)}{x})$$

并不需要手动算，除出来$O(n)$挪一下就可以了

把读入的位置拉出来直接丢到$A(x)$的对应位置，然后推出$F(x)$即可

#include <iostream>
#include <cstdio>
#include <cstring>
#include <cctype>
#include <ctime>
#define MAXN 400005
using namespace std;
const int MOD=998244353;
inline int read()
{
	int ans=0;
	char c=getchar();
	while (!isdigit(c)) c=getchar();
	while (isdigit(c)) ans=(ans<<3)+(ans<<1)+(c^48),c=getchar();
	return ans;
}
typedef long long ll;
inline int add(const int& x,const int& y){return x+y>=MOD? x+y-MOD:x+y;}
inline int dec(const int& x,const int& y){return x<y? x-y+MOD:x-y;}
inline int qpow(int a,int p)
{
	int ans=1;
	while (p)
	{
		if (p&1) ans=(ll)ans*a%MOD;
		a=(ll)a*a%MOD;p>>=1;
	}
	return ans;
}
#define inv(x) qpow(x,MOD-2)
int fac[MAXN],finv[MAXN];
int rt[2][24];
int l,r[MAXN];
inline void init(){for (int i=0;i<(1<<l);i++) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));}
inline void NTT(int* a,int type)
{
	int lim=1<<l;
	for (int i=0;i<lim;i++) if (i<r[i]) swap(a[i],a[r[i]]);
	for (int L=0;L<l;L++)
	{
		int mid=1<<L,len=mid<<1,Wn=rt[type][L+1];
		for (int s=0;s<lim;s+=len)
			for (int k=0,w=1;k<mid;k++,w=(ll)w*Wn%MOD)
			{
				int x=a[s+k],y=(ll)w*a[s+mid+k]%MOD;
				a[s+k]=add(x,y);a[s+mid+k]=dec(x,y);
			}
	}
	if (type)
	{
		int t=inv(lim);
		for (int i=0;i<lim;i++) a[i]=(ll)a[i]*t%MOD;
	}
}
void getinv(int* A,int* B,int n)
{
	static int f[MAXN],t[MAXN];
	if (n==1) return (void)(*B=inv(*A));
	getinv(A,t,(n+1)>>1);
	for (l=0;(1<<l)<(n<<1);l++);
	init();
	for (int i=0;i<n;i++) f[i]=A[i];
	for (int i=n;i<(1<<l);i++) f[i]=t[i]=0;
	NTT(f,0);NTT(t,0);
	for (int i=0;i<(1<<l);i++) B[i]=(ll)t[i]*dec(2,(ll)f[i]*t[i]%MOD)%MOD;
	NTT(B,1);
	for (int i=n;i<(1<<l);i++) B[i]=0;
}
inline void deriv(int* A,int* B,int n){for (int i=0;i<n-1;i++) B[i]=(ll)A[i+1]*(i+1)%MOD;B[n-1]=0;}
inline void integ(int* A,int* B,int n){for (int i=1;i<n;i++) B[i]=(ll)A[i-1]*fac[i-1]%MOD*finv[i]%MOD;B[0]=0;}
void getln(int* A,int* B,int n)
{
	static int f[MAXN],g[MAXN];
	deriv(A,f,n);getinv(A,g,n);
	for (int i=n;i<(1<<l);i++) f[i]=g[i]=0;
	NTT(f,0);NTT(g,0);
	for (int i=0;i<(1<<l);i++) f[i]=(ll)f[i]*g[i]%MOD;
	NTT(f,1);
	integ(f,B,n);
	for (int i=n;i<(1<<l);i++) B[i]=0;
}
void getexp(int* A,int* B,int n)
{
	static int f[MAXN],g[MAXN];
	if (n==1) return (void)(*B=1);
	getexp(A,g,(n+1)>>1);getln(g,f,n);
	for (int i=0;i<n;i++) f[i]=dec(A[i],f[i]);
	++f[0];
	for (int i=n;i<(1<<l);i++) f[i]=g[i]=0;
	NTT(f,0);NTT(g,0);
	for (int i=0;i<(1<<l);i++) B[i]=(ll)f[i]*g[i]%MOD;
	NTT(B,1);
	for (int i=n;i<(1<<l);i++) B[i]=0;
}
void getpow(int* A,int* B,int n,int k)
{
	static int t[MAXN];
	getln(A,t,n);
	for (int i=0;i<n;i++) t[i]=(ll)t[i]*k%MOD;
	getexp(t,B,n);
}
int G[MAXN],H[MAXN],dH[MAXN],A[MAXN],t[MAXN],f[MAXN];
int LangInv(int* F,int n,bool ex=false)
{
	static int f[MAXN],g[MAXN];
	for (int i=0;i<n;i++) f[i]=F[i+1];
	f[n]=0;
	getinv(f,g,n);
	getpow(g,f,n,n);
	int ans=0;
	if (ex)	for (int i=0;i<n;i++) ans=add(ans,(ll)dH[i]*f[n-i-1]%MOD);
	else ans=f[n-1];
	return (ll)ans*fac[n-1]%MOD*finv[n]%MOD;
}
time_t START;
int main()
{
	fac[0]=1;
	for (int i=1;i<MAXN;i++) fac[i]=(ll)fac[i-1]*i%MOD;
	finv[MAXN-1]=inv(fac[MAXN-1]);
	for (int i=MAXN-2;i>=0;i--) finv[i]=(ll)finv[i+1]*(i+1)%MOD;
	rt[0][23]=qpow(3,119);rt[1][23]=inv(rt[0][23]);
	for (int i=22;i>=0;i--)
	{
		rt[0][i]=(ll)rt[0][i+1]*rt[0][i+1]%MOD;
		rt[1][i]=(ll)rt[1][i+1]*rt[1][i+1]%MOD;
	}
	int n=read();
	for (int i=0;i<=n;i++) t[i]=(ll)qpow(2,((ll)i*(i-1)/2)%(MOD-1))*finv[i]%MOD;
	getln(t,G,n+1);
	for (int i=0;i<=n;i++) G[i]=(ll)G[i]*i%MOD;
	for (int i=0;i<n;i++) t[i]=G[i+1];
	t[n]=0;
	getln(t,H,n+1);deriv(H,dH,n+1);
	for (int s=read();s;s--)
	{
		int p=read()-1;
		A[p]=LangInv(G,p,true);
	}
	getexp(A,t,n+1);getinv(t,A,n+1);
	for (int i=1;i<=n;i++) f[i]=A[i-1];
	f[0]=0;
	printf("%d\n",(ll)LangInv(f,n)*fac[n-1]%MOD);
	return 0;
}

P.S. 如果你用LOJ的C++（NOI） T成了狗，试试C++……