• 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(nlog n+sum x_ilog 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
    ",mul(qpow(n,mod-2),mul(ans, fac[n])));
    }
    
  • 相关阅读:
    【原创】编程题练习:反转字符串中的单词
    【最近的学习安排】
    【转载】判断两个链表是否相交、一个链表是否有环
    获取字符串字节长度
    如何找到GridView里的控件,建立GridViewRow对象
    Flex和.NET协同开发利器FluorineFx Flex与.NET互操作
    两款基于Visual Studio开发Flex的插件
    Mysql:向信号量添加给定计数将导致其超出它的最大计数错误
    Sql Server数据库触发器实例
    国外一些知名ASP.Net开源CMS系统
  • 原文地址:https://www.cnblogs.com/weiyanpeng/p/13125163.html
Copyright © 2020-2023  润新知