zoukankan      html  css  js  c++  java
  • P4491 [HAOI2018]染色

    P4491 [HAOI2018]染色

    出题人用心险恶啊。。。特意空格强调了“恰好”二字

    如果恰 好出现了 (S) 次的颜色有 (K)

    这题里有两个“恰好”,一个是恰好出现了 (S) 次,一个是恰好 (K) 种颜色。

    这种一般都是转成至多至少反演或者容斥。

    容斥掉的根本不是出题人强调的那个恰好。。。

    一开始我就进坑了,不过还好没几分钟就跳了出来并想到了下面的做法

    稍微想一想就知道我们应该算的是:(f(k)= m 至少 k 种颜色恰好出现 S 次) 的方案数。

    然而其实会算重的,并不严格是方案数,只不过重复是可以计算的。事实上二项式反演好像都是这么干的,非常巧妙。

    这个 (f) 非常好算。

    [f(k)=inom{m}{k}dfrac{inom{n}{Sk}(Sk)!}{(S!)^{k}}(m-k)^{n-Sk} ]

    首先选出 (k) 种颜色,钦定出现 (S) 次,其余乱选,这样必然至少 (k) 种颜色出现 (S) 次。

    这样会钦定 (Sk) 个位置,那么直接选好即为 (inom{n}{Sk}) 。同时这些元素可以互换位置,那么就是全排列除掉每一种颜色内部换顺序的情况。剩下 (n-Sk) 个位置,每一个位置都有 (m-k) 种颜色可选。

    把分子那个二项式系数拆成阶乘,稍微化简一下得

    [f(k)=inom{m}{k}dfrac{n!}{(S!)^{k}(n-Sk)!}(m-k)^{n-Sk} ]

    设恰好选 (k) 种颜色的方案数设为 (g(k))(up=min(m,lfloordfrac{n}{S} floor)) 是选的颜色种数上界。

    这个 (g) 是真的方案数了

    再来说上面插的那句话

    然而其实会算重的,并不严格是方案数

    (这个其实是翻command_block的blog时发现的问题,刚开始学二项式反演被那个错误的叫法迷惑了好久,到现在才解决,虽然记个式子就能做题了。向cmd表示感谢。)

    那么究竟重复了多少呢?

    首先肯定,重复的部分是乱选与钦定出现了相同的方案。

    那么对于 (i(ige k)) 种颜色相同,会在钦定 (k) 种相同时重复 (dbinom{i}{k}) 次。

    比如 (i=4,k=2)({1,2,3,4}) 这个集合会在钦定 ({1,2},{1,3},{1,4},{2,3},{2,4},{3,4}) 的时候分别被计算一次。

    所以:

    [f(k)=sum_{i=k}^{up}inom{i}{k}g(i) ]

    二项式反演得

    [g(k)=sum_{i=k}^{up}(-1)^{i-k}inom{i}{k}f(i) ]

    这个式子得想个办法卷起来,没法化了,而且看着很可卷的样子。

    卷王能一眼看出这个怎么卷,我这种菜鸡只能慢慢试。。。

    如果您是擅长卷可以直接跳过下面这部分。

    组合数肯定得拆

    [g(k)=dfrac{1}{k!}sum_{i=k}^{up}f(i)i!dfrac{(-1)^{i-k}}{(i-k)!} ]

    (A(i)=f(i)i!,B(i)=dfrac{(-1)^{i}}{i!})

    [g(k)=dfrac{1}{k!}sum_{i=k}^{up}A(i)B(i-k) ]

    必然要翻转一个。翻转哪个好呢?试一试就知道应该翻转 (A) ,因为这样下标从 (0) 开始才能卷。

    那么

    [g(k)=dfrac{1}{k!}sum_{i=k}^{up}A'(up-i)B(i-k)\ =dfrac{1}{k!}sum_{i=0}^{up-k}A'(up-i-k)B(i) ]

    总算可以卷了/ll,卷完 ([x^k]g(x))(up-k) 处提取。

    出题人挺良心的,样例看着很强的样子=_=

    #include<bits/stdc++.h>
    using namespace std;
    #define fi first
    #define se second
    #define mkp(x,y) make_pair(x,y)
    #define pb(x) push_back(x)
    #define sz(v) (int)v.size()
    typedef long long LL;
    typedef double db;
    template<class T>bool ckmax(T&x,T y){return x<y?x=y,1:0;}
    template<class T>bool ckmin(T&x,T y){return x>y?x=y,1:0;}
    #define rep(i,x,y) for(int i=x,i##end=y;i<=i##end;++i)
    #define per(i,x,y) for(int i=x,i##end=y;i>=i##end;--i)
    inline int read(){
    	int x=0,f=1;char ch=getchar();
    	while(!isdigit(ch)){if(ch=='-')f=0;ch=getchar();}
    	while(isdigit(ch))x=x*10+ch-'0',ch=getchar();
    	return f?x:-x;
    }
    #define mod 1004535809
    const int V=10000005;
    const int N=100005;
    const int M=N<<2;
    
    namespace math{
    int fac[V],ifc[V];
    inline int qpow(int n,int k){int res=1;for(;k;k>>=1,n=1ll*n*n%mod)if(k&1)res=1ll*n*res%mod;return res;}
    inline void fmod(int&x){x-=mod,x+=x>>31&mod;}
    int binom(int n,int m){return n<m?0:1ll*fac[n]*ifc[m]%mod*ifc[n-m]%mod;}
    void initmath(const int&n=V-1){
    	fac[0]=1;for(int i=1;i<=n;++i)fac[i]=1ll*i*fac[i-1]%mod;
    	ifc[n]=qpow(fac[n],mod-2);for(int i=n-1;i>=0;--i)ifc[i]=1ll*ifc[i+1]*(i+1)%mod;
    }
    
    }
    using math::qpow;
    using math::fmod;
    
    namespace poly{
    int rev[M],lg,lim;
    void poly_init(const int&n){
    	for(lg=0,lim=1;lim<n;++lg,lim<<=1);
    	for(int i=0;i<lim;++i)rev[i]=(rev[i>>1]>>1)|((i&1)<<(lg-1));
    }
    void NTT(int*a,int op){
    	for(int i=0;i<lim;++i)
    		if(i>rev[i])swap(a[i],a[rev[i]]);
    	const int g=op?3:qpow(3,mod-2);
    	for(int i=1;i<lim;i<<=1){
    		const int wn=qpow(g,(mod-1)/(i<<1));
    		for(int j=0;j<lim;j+=i<<1){
    			int w0=1;
    			for(int k=0;k<i;++k,w0=1ll*w0*wn%mod){
    				const int X=a[j+k],Y=1ll*a[i+j+k]*w0%mod;
    				fmod(a[j+k]=X+Y),fmod(a[i+j+k]=mod+X-Y);
    			}
    		}
    	}
    	if(op)return;const int ilim=qpow(lim,mod-2);
    	for(int i=0;i<lim;++i)a[i]=1ll*a[i]*ilim%mod;
    }
    #define clr(a,n) memset(a,0,sizeof(int)*(n))
    #define cpy(a,b,n) memcpy(a,b,sizeof(int)*(n))
    void poly_mul(int*f,int*g,int*ans,int n,int m){
    	static int A[M],B[M];poly_init(n+m);
    	cpy(A,f,n),clr(A+n,lim-n),NTT(A,1);
    	cpy(B,g,m),clr(B+m,lim-m),NTT(B,1);
    	for(int i=0;i<lim;++i)ans[i]=1ll*A[i]*B[i]%mod;
    	NTT(ans,0);
    }
    
    }
    
    int n,m,S,up,f[M],g[M],A[M],B[M],C[M],W[N],ans;
    signed main(){
    	math::initmath();
    	n=read(),m=read(),S=read(),up=min(n/S,m);
    	rep(i,0,m)W[i]=read();
    	for(int i=0;i<=up;++i)f[i]=1ll*math::binom(m,i)*math::fac[n]%mod*qpow(math::ifc[S],i)%mod*math::ifc[n-S*i]%mod*qpow(m-i,n-S*i)%mod;
    	for(int i=0;i<=up;++i)A[i]=1ll*math::fac[up-i]*f[up-i]%mod;
    	for(int i=0;i<=up;++i)B[i]=(i&1)?mod-math::ifc[i]:math::ifc[i];
    	poly::poly_mul(A,B,C,up+1,up+1);
    	for(int i=0;i<=up;++i)g[i]=1ll*math::ifc[i]*C[up-i]%mod;
    	for(int i=0;i<=up;++i)fmod(ans+=1ll*g[i]*W[i]%mod);
    	printf("%d
    ",ans);
    	return 0;
    }
    
  • 相关阅读:
    (转)消息队列 Kafka 的基本知识及 .NET Core 客户端
    Neo4j学习笔记
    科技论文推荐系统
    下载pubmed数据
    杂项
    Scrapy 知乎验证码
    Scrapy 爬取网站文章
    爬虫基础知识
    Django linux uWsgi Nginx 部署
    DocumentSimilarity
  • 原文地址:https://www.cnblogs.com/zzctommy/p/14225121.html
Copyright © 2011-2022 走看看