zoukankan      html  css  js  c++  java
  • 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])));
    }
    
  • 相关阅读:
    7.21 高博教育 数组 内存
    【基础扎实】Python操作Excel三模块
    PAT 甲级 1012 The Best Rank
    PAT 甲级 1011  World Cup Betting
    PAT 甲级 1010 Radix
    链式线性表——实验及提升训练
    循环程序设计能力自测
    链表应用能力自测
    PAT 甲级 1009 Product of Polynomials
    1008 Elevator (20分)
  • 原文地址:https://www.cnblogs.com/weiyanpeng/p/13125163.html
Copyright © 2011-2022 走看看