zoukankan      html  css  js  c++  java
  • [LibreOJ 3120]【CTS2019】珍珠 【生成函数】【计数】

    Description

    在这里插入图片描述
    在这里插入图片描述

    Solution

    有一个直观的思路是考虑每种颜色个数的奇偶性,奇数个数的颜色不能超过(n-2m)
    因此若(n-2mgeq D)则答案一定是(D^n)

    否则由于每种颜色其实没有区别,我们考虑一种颜色为奇数和为偶数的指数型生成函数
    奇数是(e^x-e^{-x}over 2),偶数是(e^x+e^{-x}over 2)

    我们枚举有多少个奇数的颜色
    容易得到

    [Ans=n!sumlimits_{i=0}^{n-2m}{Dchoose i}left({e^x-e^{-x}over 2} ight)^ileft({e^x+e^{-x}over 2} ight)^{D-i}[x^n] ]

    我考场上写的是这个式子(和题解的本质相同,不过化起来比较麻烦)

    提出一个(2^{-D}),把后面的东西二项式展开

    [=2^{-D}n!sumlimits_{i=0}^{n-2m}{Dchoose i}sumlimits_{p=0}^{i}e^{(2p-i)x}{ichoose p}(-1)^{i-p}sumlimits_{q=0}^{D-i}e^{(2q+i-D)x}{D-ichoose q}[x^n] ]

    此处我们可以枚举(T=p+q),并移到最外层

    [=2^{-D}n!sumlimits_{T=0}^{D}e^{2T-D}sumlimits_{i=0}^{n-2m}{Dchoose i}sumlimits_{p+q=T}(-1)^{i-p}{ichoose p}{D-ichoose q}[x^n] ]

    容易知道(n!e^{px}[x^n]=p^n)
    把q换成T-p

    [=2^{-D}sumlimits_{T=0}^{D}(2T-D)^nsumlimits_{i=0}^{n-2m}sumlimits_{p=0}^{T}(-1)^{i-p}{Dchoose i}{ichoose p}{D-ichoose T-p} ]

    考场上的时候我就卡在这里推不动了
    实际上那三个组合数可以化开成({D!over i!(D-i)!}{i!over (i-p)!p!}{(D-i)!over (T-p)!(D-i-T+p)!})

    分配阶乘,约分,补上一个((D-T)!T!over (D-T)!T!)
    可以得到

    [=2^{-D}sumlimits_{T=0}^{D}(2T-D)^n{Dchoose T}sumlimits_{i=0}^{n-2m}sumlimits_{p=0}^{T}(-1)^{i-p}{Tchoose p}{D-Tchoose i-p} ]

    可以发现后面是两个二项式卷积的形式
    实际上就是(left((1+y)^T(1-y)^{D-T} ight)[y^i])
    后面的就和题解是一样的了。

    题解的推法是这样的

    [Ans=n!sumlimits_{i=0}^{n-2m}left(y{e^x-e^{-x}over 2}+{e^x+e^{-x}over 2} ight)^D[x^n][y^i] ]

    (e^x)(e^{-x})分开

    [Ans=2^{-D}n!sumlimits_{i=0}^{n-2m}left(e^x(1+y)+e^{-x}(1-y) ight)^D[x^n][y^i] ]

    二项式展开

    [Ans=2^{-D}n!sumlimits_{i=0}^{n-2m}sumlimits_{T=0}^{D}e^{(2T-D)x}(1+y)^T(1-y)^{D-T}{Dchoose T}[x^n][y^i] ]

    (i)放到后面去,就是

    [=2^{-D}sumlimits_{T=0}^{D}(2T-D)^n{Dchoose T}sumlimits_{i=0}^{n-2m}left((1+y)^T(1-y)^{D-T} ight)[y^i] ]

    (看来是我的推法太蠢了)

    记后面的东西(F(T,D)=sumlimits_{i=0}^{n-2m}left((1+y)^T(1-y)^{D-T} ight)[y^i])
    接下来就是高端操作时间
    把一个((1+y))拆成(-(1-y)+2),式子就可以化开成两边,具体略去
    立刻可以得到(F(T,D)=-F(T-1,D)+2F(T-1,D-1))

    然后(F(0,D)=sumlimits_{i=0}^{n-2m}{Dchoose i}(-1)^i)
    通过杨辉三角的性质,化到上一行去发现都消掉了,结果就是({D-1choose n-2m}(-1)^{n-2m})

    然后就可以通过NTT加速递推过程求出(F(T,D))了。

    实际上网上似乎有一种很简单的推法,利用容斥,变成至少i个为奇数,然后式子就好化很多,结果就是两次直接的卷积,就不需要后面的高端操作了。
    时间复杂度(O(nlog n))

    Code

    #include <bits/stdc++.h>
    #define fo(i,a,b) for(int i=a;i<=b;++i)
    #define fod(i,a,b) for(int i=a;i>=b;--i)
    #define N 100005
    #define M 262144
    #define T 18
    #define mo 998244353
    #define LL long long
    using namespace std;
    LL l,n,m,js[M+1],ns[M+1],ny[M+1];
    LL ksm(LL k,LL n)
    {
    	k=(k+mo)%mo;
    	LL s=1;
    	for(;n;n>>=1,k=k*k%mo) if(n&1) s=s*k%mo;
    	return s;
    }
    LL C(int n,int m)
    {
    	if(n<m||n<0||m<0) return 0;
    	return js[n]*ns[m]%mo*ns[n-m]%mo;
    }
    int a[M+1],b[M+1];
    LL F(int w)
    {
    	if(w==0) return 1;
    	return ((n-2*m)&1)?(mo-C(w-1,n-2*m))%mo:C(w-1,n-2*m);
    }
    
    int bit[M+1];
    int wi[M+1];
    namespace polynomial
    {
    	void prp()
    	{
    		LL v=ksm(3,(mo-1)/M);
    		wi[0]=1;
    		fo(i,1,M) 
    		{
    			wi[i]=(LL)wi[i-1]*v%mo;
    			bit[i]=(bit[i>>1]>>1)|((i&1)<<(T-1));
    		}
    	}
    	int inc(int a,int b)
    	{
    		return (a+=b)>=mo?a-mo:a;
    	}
    	int dec(int a,int b)
    	{
    		return (a-=b)<0?a+mo:a;
    	}
    	void DFT(int *a)
    	{
    		fo(i,0,M-1) if(bit[i]<i) swap(a[i],a[bit[i]]);
    		for(int h=1,l=(M>>1),v;h<M;h<<=1,l>>=1) 
    		{
    	   		for(int j=0;j<M;j+=h<<1) 
    	   		{
    	        	int *x=a+j,*y=x+h,*w=wi;
    				for(int i=0;i<h;++i,++x,++y,w+=l)
    				{
    	        		v=((LL)*y * *w)%mo;
    	        		*y=dec(*x,v);
    	        		*x=inc(*x,v);
    	        	}
    	      	}
    	    }
    	}
    	void IDFT(int *a)
    	{
    		DFT(a);
    		fo(i,0,M-1) a[i]=a[i]*ny[M]%mo;
    		reverse(a+1,a+M);
    	}
    }
    using polynomial::prp;
    using polynomial::DFT;
    using polynomial::IDFT;
    int main()
    {
    	cin>>l>>n>>m;
    	if(n-2*m>=l) 
    	{
    		printf("%lld
    ",ksm(l,n));
    		return 0;
    	} 
    	js[0]=js[1]=ny[1]=ns[0]=ns[1]=1;
    	fo(i,2,M) 
    	{
    		js[i]=js[i-1]*(LL)i%mo;
    		ny[i]=(-ny[mo%i]*(LL)(mo/i)%mo+mo)%mo;
    		ns[i]=ns[i-1]*ny[i]%mo; 
    	}
    	LL v=1;
    	fo(i,0,l) 
    	{
    		a[i]=v*ns[i]%mo*F(l-i)%mo,b[i]=(LL)((i&1)?mo-1:1)*ns[i]%mo;
    		v=v*(LL)2%mo; 
    	}
    	prp();
    	DFT(a),DFT(b);
    	fo(i,0,M-1) a[i]=(LL)a[i]*(LL)b[i]%mo;
    	IDFT(a);
    	LL ans=0;
    	fo(i,0,l) ans=(ans+a[i]*js[i]%mo*ksm(2*i-l,n)%mo*C(l,i))%mo;
    	printf("%lld
    ",ans*ksm(ksm(2,l),mo-2)%mo);
    }
    
  • 相关阅读:
    「从零单排canal 03」 canal源码分析大纲
    「从零单排canal 02」canal集群版 + admin控制台 最新搭建姿势(基于1.1.4版本)
    「从零单排canal 01」 canal 10分钟入门(基于1.1.4版本)
    实时数据订阅与分发系统概述
    使用phoenix踩的坑与设计思考
    「从零单排HBase 12」HBase二级索引Phoenix使用与最佳实践
    「从零单排HBase 11」HBase二级索引解决方案
    「从零单排HBase 10」HBase集群多租户实践
    「从零单排HBase 09」HBase的那些数据结构和算法
    Netty源码分析之自定义编解码器
  • 原文地址:https://www.cnblogs.com/BAJimH/p/10902115.html
Copyright © 2011-2022 走看看