Loj #6363. 「地底蔷薇」

考虑给一个根。记 \(B\) 是有根联通图，\(D\) 是点双连通图。

现在考虑有根无向图：

\[B(x) = x*\exp(\sum_i D_{i+1}/i! B^i) \\ \frac{B(x)}{\exp(D'(B(x)))}=x \]

扩展拉格朗日反演：

\[[x^n] H(\frac{x}{\exp(D'(x))}) = \frac{1}{n}[x^{n-1}]H'(x)\frac{x^n}{B(x)^n} \]

取 \(H(x) = \ln(\frac{B(x)}{x})\)，左边即为 \(D'(x)​\)。

现在可以得到一个 \(D_s​\)，表示所有集合内的。

设答案是 \(F\)，有：

\[F(x) = x*\exp(D_s'(F(x))) \\ \frac{F(x)}{\exp(D_s'(F(x)))} = x \]

可以直接拉格朗日反演：

\[[x^n]F(x) = \frac{1}{n}[x^{n-1}](\frac{x}{\frac{x}{\exp(D_s'(x))}})^n \\=\frac{1}{n}[x^{n-1}]\exp(D_s'(x))^n \]

最后除点数就是答案。

复杂度 \(O(n\log n+\sum x_i\log x_i)\)。

代码

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod=998244353;
inline int add(int a,int b){a+=b;return a>=mod?a-mod:a;}
inline int sub(int a,int b){a-=b;return a<0?a+mod:a;}
inline int mul(int a,int b){return 1ll*a*b%mod;}
inline int qpow(int a,int b){int ret=1;for(;b;b>>=1,a=mul(a,a))if(b&1)ret=mul(ret,a);return ret;}
const int inv_2=(mod+1)>>1;
/* math */

typedef vector<int> poly;
namespace polynomial{
const int Ntt_Lim = 4e5+6;
	int rev[Ntt_Lim],_w[Ntt_Lim];
	const int G_mod = 3;
	poly deri(poly a){
		for(int i=0;i+1<(int)a.size();i++)a[i]=mul(a[i+1],i+1);
		a.pop_back();return a;
	}
	poly inte(poly a){
		a.push_back(0);for(int i=(int)a.size()-2;~i;i--)a[i+1]=mul(a[i],qpow(i+1,mod-2));
		a[0]=0;return a;
	}
	poly p_add(poly a,poly b){a.resize(max(a.size(),b.size()));for(size_t i=0;i<b.size();i++)a[i]=add(a[i],b[i]);return a;}
	poly p_sub(poly a,poly b){a.resize(max(a.size(),b.size()));for(size_t i=0;i<b.size();i++)a[i]=sub(a[i],b[i]);return a;}
	inline void _w_init(){
		for(int step=1;step*2<=Ntt_Lim;step<<=1){
			int wn = qpow(G_mod, (mod-1)/(step<<1));
			for(int j=step,w=1;j<step<<1;j++,w=mul(w,wn)){
				_w[j]=w;
			}
		}
	}
	inline void dft(int *f,int len,int type){
		int l=0;while(1<<l<len)++l;
		for(int i=0;i<len;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
		for(int i=0;i<len;i++)if(i<rev[i])swap(f[i],f[rev[i]]);
		for(int step=1;step<len;step<<=1){
//			int wn=_w[step];// int wn = qpow(G_mod, (mod-1)/(step<<1));
			for(int i=0;i<len;i+=step<<1)for(int x,y,j=0;j<step;j++){
				x=f[i+j],y=mul(_w[j+step],f[i+j+step]);
				f[i+j]=add(x,y),f[i+j+step]=sub(x,y);
			}
		}
		if(type==1)return;
		int inv=qpow(len,mod-2);reverse(f+1,f+len);
		for(int i=0;i<len;i++)f[i]=mul(f[i],inv);
	}
	poly ntt(poly a,poly b,int n,int m){
		int l=1;while(l<n+m-1)l<<=1;
		a.resize(l),b.resize(l);dft(&a[0],l,1),dft(&b[0],l,1);
		for(int i=0;i<l;i++)a[i]=mul(a[i],b[i]);
		dft(&a[0],l,-1);a.resize(n+m-1);
		return a;
	}
	poly ntt(poly a,poly b){return ntt(a,b,a.size(),b.size());}
	poly inv(poly &f,int n){
		if(n==1)return poly(1,qpow(f[0],mod-2));
		poly a(&f[0],&f[n]),b=inv(f,(n+1)/2);int l=1;while(l<n<<1)l<<=1;
		a.resize(l),b.resize(l);
		dft(&a[0],l,1),dft(&b[0],l,1);
		for(int i=0;i<l;i++)a[i]=mul(b[i], sub(2,mul(a[i],b[i])));
		dft(&a[0],l,-1);a.resize(n);
		return a;
	}
	poly inv(poly a){return inv(a,a.size());}
	poly sqrt(poly &f,int n){
		if(n==1)return poly(1,1);
		poly a(&f[0],&f[n]),b=sqrt(f,(n+1)/2);
		b.resize(n);a=ntt(a,inv(b));a.resize(n);
		for(int i=0;i<n;i++)a[i]=mul(inv_2, add(a[i],b[i]));
		return a;
	}
	poly sqrt(poly a){return sqrt(a,a.size());}
	poly ln(poly a){
		int l=a.size();a=inte(ntt(deri(a),inv(a)));
		a.resize(l);return a;
	}
	poly exp(poly& f,int n){
		if(n==1)return poly(1,1);//f[0]=1
		poly a(n,0),b=exp(f,(n+1)/2);
		b.resize(n);a=ln(b);
		for(int i=0;i<n;i++)a[i]=sub(f[i],a[i]);a[0]=add(a[0],1);
		a=ntt(a,b);a.resize(n);
		return a;
	}
	poly exp(poly a){return exp(a,a.size());}
	pair<poly,poly> div(poly a,poly b){//assert(a.size()>=b.size())
		if(a.size()<b.size())return make_pair(poly(1,0),a);
		int n=a.size(),m=b.size();
		poly ra=a,rb=b;
		reverse(ra.begin(),ra.end()),reverse(rb.begin(),rb.end());
		ra.resize(n-m+1),rb.resize(n-m+1);
		poly c=ntt(ra,inv(rb));c.resize(n-m+1);reverse(c.begin(),c.end());
		poly d=p_sub(a,ntt(b,c));d.resize(m-1);
		return make_pair(c,d);
	}
}
using namespace polynomial;
int n,t;
const int N = 2e5+5;
int fac[N], ifac[N];
inline void init(int n=2e5){
	fac[0]=ifac[0]=1;for(int i=1;i<=n;i++)fac[i]=mul(fac[i-1],i);
	ifac[n]=qpow(fac[n],mod-2);for(int i=n-1;i;i--)ifac[i]=mul(ifac[i+1],i+1);
}
poly B,H_,_B,Ds_;
inline int Solve_S(int n){
	if(n==1)return 1;
	--n;
	poly a(&H_[0]+0, &H_[0]+n);
	poly b(&_B[0]+0, &_B[0]+n);
	for(int i=0;i<n;i++)b[i]=mul(b[i],n);
	poly ret=ntt(a,exp(b));
	return mul(mul(qpow(n,mod-2),ret[n-1]),fac[n]);
}

int main()
{
	_w_init();
	init();
	cin >> n >> t;
	Ds_.resize(n);
	B.resize(n+2);
	for(int i=0;i<=n+1;i++)
		B[i]=mul(ifac[i],qpow(2,(1ll*i*(i-1)/2)%(mod-1)));
	B=ln(B);
	for(int i=0;i<=n;i++)B[i]=mul(B[i],i);
	_B.resize(n+1);for(int i=0;i<=n;i++)_B[i]=B[i+1];
	H_=ntt(deri(_B),inv(_B));
	H_.resize(n);
	_B=ln(inv(_B));
	while(t--){
		int q;scanf("%d",&q);
		if(q==1)continue;
		int tmp=Solve_S(q);
		// cerr << tmp << endl;
		Ds_[q-1]=mul(ifac[q-1],tmp);
	}
	for(int i=0;i<n;i++)Ds_[i]=mul(Ds_[i],n);
	Ds_=exp(Ds_);
	int ans=mul(qpow(n,mod-2),Ds_[n-1]);
	printf("%d\n",mul(qpow(n,mod-2),mul(ans, fac[n])));
}