考虑给一个根。记 \(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])));
}