• BZOJ3456城市规划


    BZOJ3456

    http://www.lydsy.com/JudgeOnline/problem.php?id=3456

    #include<cstdio>
    #include<cstdlib>
    #include<algorithm>
    #include<cstring>
    using namespace std;
    const int mod=1004535809,maxn=1<<18|1;
    int p[maxn],gn[233];
    int n;
    int fac[maxn],inv[maxn],tp[maxn];
    inline int fp(int a,int b){
    	int res=1;
    	while(b){
    		if(b&1)res=1ll*res*a%mod;
    		a=1ll*a*a%mod;b>>=1;
    	}
    	return res;
    }
    inline void init(){
    	for(register int i=0;i<=23;++i)gn[i]=fp(3,(mod-1)/(1<<(i+1)));
    	fac[0]=1;for(register int i=1;i<=n+1;++i)fac[i]=1ll*fac[i-1]*i%mod;
    	inv[1]=1;for(register int i=2;i<=n+1;++i)inv[i]=mod-1ll*mod/i*inv[mod%i]%mod;
    	inv[0]=1;for(register int i=1;i<=n+1;++i)inv[i]=1ll*inv[i]*inv[i-1]%mod;
    	for(register int i=0;i<=n+1;++i)tp[i]=fp(2,(1ll*i*(i-1)/2)%(mod-1));
    }
    struct poly{
    	int a[maxn],n;
    	poly(){memset(a,0,sizeof(a));}
    	inline int &operator[](const int x){return a[x];}
    	inline void ntt(int d,int f){
    		for(register int i=0;i<d;++i)if(i<p[i])swap(a[i],a[p[i]]);
    		for(register int i=1,t=0,w,v;i<d;i<<=1,++t){
    			for(register int j=0;j<d;j+=(i<<1)){
    				w=1;
    				for(register int k=j;k<i+j;++k,w=1ll*w*gn[t]%mod)
    					v=1ll*w*a[i+k]%mod,a[i+k]=(a[k]-v+mod)%mod,a[k]=(a[k]+v)%mod;
    			}
    		}
    		if(f==1)return;
    		reverse(a+1,a+d);
    		register int ny=fp(d,mod-2);
    		for(register int i=0;i<d;++i)a[i]=1ll*a[i]*ny%mod;
    	}
    	friend inline poly operator*(poly &A,poly &B){
    		register poly res;res.n=A.n+B.n;int d,lg2;
    		for(d=1,lg2=0;d<=res.n;d<<=1)++lg2;
    		for(register int i=0;i<d;++i)p[i]=(p[i>>1]>>1)^((i&1)<<(lg2-1));
    		A.ntt(d,1);B.ntt(d,1);
    		for(register int i=0;i<d;++i)res[i]=1ll*A[i]*B[i]%mod;
    		res.ntt(d,-1);A.ntt(d,-1);B.ntt(d,-1);
    		return res;
    	}
    }A,B,tmp,invA,Ans;
    inline void poly_inv(int deg){
    	if(deg==1){
    		invA[0]=fp(A[0],mod-2);
    		return ;
    	}
    	poly_inv((deg+1)>>1);
    	register int d=1,lg2=0;for(d=1;d<=(deg<<1);d<<=1)++lg2;
    	for(register int i=0;i<d;++i)p[i]=(p[i>>1]>>1)^((i&1)<<(lg2-1));
    	for(register int i=0;i<d;++i)tmp[i]=i<deg?A[i]:0;
    	for(register int i=(deg+1)>>1;i<d;++i)invA[i]=0;
    	invA.ntt(d,1);tmp.ntt(d,1);
    	for(register int i=0;i<d;++i)invA[i]=(mod+2ll*invA[i]%mod-1ll*invA[i]*invA[i]%mod*tmp[i]%mod)%mod;
    	invA.ntt(d,-1);invA.n=deg-1;
    }
    int main(){
    	scanf("%d",&n);init();A.n=B.n=n;
    	for(register int i=0;i<=n;++i)A[i]=1ll*tp[i]*inv[i]%mod;
    	for(register int i=1;i<=n;++i)B[i]=1ll*tp[i]*inv[i-1]%mod;
    	poly_inv(n+1);
    	printf("%d
    ",1ll*(invA*B)[n]*fac[n-1]%mod);	
    	return 0;
    }
    

      

  • 相关阅读:
    group_concat的长度限制
    mb_strlen默认字符集问题
    &符号导致的一个bug
    python面向对象编程系列
    python基础之面向过程编程系列
    RPA流程自动化
    什么是DevOps?
    ansible详解
    saas和paas的区别
    CPT/cpt接口
  • 原文地址:https://www.cnblogs.com/Stump/p/8454369.html
Copyright © 2020-2023  润新知