【题解】[HAOI2018]染色(NTT+容斥/二项式反演)

【题解】[HAOI2018]染色(NTT+容斥/二项式反演)

可以直接写出式子：

\[f(x)={m \choose x}n!{(\dfrac 1 {(Sx)!})}^x(m-x)^{n-Sx}\dfrac 1 {(n-Sx)!} \]

\(f(x)\) 钦定有\(x\)种颜色出现了恰好\(S\)的方案

然后推一下恰好有\(x\)种颜色出现了恰好\(S\)次的方案\(g(x)\) 。推导在下下面。

最后的答案是\(\sum w_i g(i)\)

推导：

显然颜色种类不会超过\(L=\lfloor \dfrac n s \rfloor\) ，然后限制一下\(L\le n ,L\le m\)

先有关系式

\[f(x)=\sum_{i=x}^L {i\choose x}g(i) \]

(恰好出现多次，会有\({i\choose x}\)被钦定的可能。)

再容斥(或者叫二项式反演?)得到

\[g(x)=\sum_{i=x}^L (-1)^{i-x}{i\choose x}f(i) \]

解释一波这个容斥的含义

由于对于每个有\(u\)种颜色出现\(S\)的方案，有三种状态：

• 当\(u\)少于目标\(x\)时：在\(f(x)\)中被计算了\(0\)次。

• 当\(u\)等于目标\(x\)时：在\(f(x)\)中只被计算了\(1\)次。

• 当\(u\)大于目标\(x\)时：

  考虑这个方案，在\(f(x)\)到\(f(u)\)之中分别计算了

  \[{u\choose x},{u\choose x+1},\dots{u\choose u} \]

  次，很像一个二项式定理展开。

  然后我们考虑一下\({x+i\choose x}f(x+i)\)的情况

  \[{u\choose x}{x\choose x},{u\choose x+1}{x+1\choose x},\dots{u\choose u}{u\choose x} \]

  晓得一个公式，那就是

  \[{a\choose b}{b\choose c}={a\choose c}{a-c\choose b-c} \]

  传送门:【总结】组合模型及其组合意义的阐释

  用这个公式化简一下上面那个式子

  \[\sum_{i=0}^{u-x}{u\choose x}{u-x\choose i} \]

  好家伙，右边就是一个二项式定理啊！

  继续化简

  \[={u\choose x}\sum_{i=0}^{u-x}{u-x\choose i}={u\choose x}(1+1)^{u-x} \]

  现在我们要让等式右边变成

  \[0={u\choose x}\times 0={u\choose x}(1-1)^{u-x} \]

  我们在\(f(i)\)前面乘上一个容斥系数即可。此时，\(u\le x\)的情况不受影响。

  故

  \[\sum_{i=x}^L (-1)^{i-x}{i\choose x}f(i)=g(x) \]

然后我们把\(g(x)\)的式子化简一下

\[g(x)=\sum_{i=x}^L (-1)^{i-x}{i\choose x}f(i)=\sum_{i=x}^L (-1)^{i-x}\dfrac {i!} {(i-x)!x!}f(i) \\ \therefore x!g(x)=\sum_{i=x}^L\dfrac {(-1)^{i-x}}{(i-x)!}\times i!f(i) \]

拿出来

\[x!g(x)=\sum_{i=x}^L\dfrac {(-1)^{i-x}}{(i-x)!}\times i!f(i) \]

我们把f数组带i!给reverse一下

现在问题就是快速求\(x!g(x)\)了，考虑到\(NTT\)实际上就是求\(c_i=\sum_{k+j=i}a_kb_j\)。套进NTT就好了。注意到这不是标准的NTT式子，所以个别地方微调。

人丑代码丑，不保证所有时刻可以AC

//@winlere
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>


using namespace std;  typedef long long ll;
inline int qr(){
      register int ret=0,f=0;
      register char c=getchar();
      while(c<48||c>57)f|=c==45,c=getchar();
      while(c>=48&&c<=57) ret=ret*10+c-48,c=getchar();
      return f?-ret:ret;
}

namespace poly{
      const int maxn=1<<18|1;
      int a[maxn],b[maxn],A[maxn],B[maxn],r[maxn];
      int savlen;
      inline void getr(int len){
	    if(len==savlen)return;
	    int cnt=0;
	    for(register int t=1;t<len;t<<=1)++cnt;
	    for(register int t=1;t<len;++t)
		  r[t]=r[t>>1]>>1|(t&1)<<cnt>>1;
      }
      const int mod=1004535809;
      const int g=3;
      inline int ksm(int base,int p){
	    register int ret=1;
	    for(base%=mod;p;p>>=1,base=1ll*base*base%mod)
		  if(p&1) ret=1ll*ret*base%mod;
	    return ret;
      }
      const int gi=ksm(3,mod-2);
      
      inline void NTT(int*a,const int&len,const int&tag){
	    getr(len);
	    for(register int t=1;t<len;++t)
		  if(r[t]>t) swap(a[t],a[r[t]]);
	    int *a1,*a0,s=g;
	    if(tag!=1) s=gi;
	    for(register int t=1,wn;t<len;t<<=1){
		  wn=ksm(s,(mod-1)/(t<<1));
		  for(register int i=0;i<len;i+=t<<1){
			a1=(a0=a+i)+t;
			for(register int j=0,w=1,tm;j<t;++j,++a1,++a0,w=1ll*w*wn%mod){
			      tm=1ll**a1*w%mod;
			      *a1=(*a0-tm)%mod;
			      *a0=(*a0+tm)%mod;
			      if(*a1<0)*a1+=mod;
			}
		  }
	    }
	    if(tag!=1)
		  for(register int t=0,in=ksm(len,mod-2);t<len;++t)
			a[t]=1ll*a[t]*in%mod;
      }

      void INV(int*a,int*b,const int&len){
	    if(len==1){b[0]=ksm(a[0],mod-2);return;}
	    INV(a,b,len>>1);
	    for(register int t=0;t<len;++t) A[t]=a[t],B[t]=b[t];
	    NTT(A,len<<1,1);NTT(B,len<<1,1);
	    for(register int t=0,w=len<<1;t<w;++t) A[t]=1ll*A[t]*B[t]%mod*B[t]%mod;
	    NTT(A,len<<1,-1);
	    for(register int t=0;t<len;++t) b[t]=((b[t]+b[t])%mod-A[t]+mod)%mod;

      }

            
}
using poly::NTT;
using poly::mod;
using poly::ksm;
const int maxn=1e7+7;
int jc[maxn];
int inv[maxn];

inline void pre(){
      jc[0]=inv[0]=1;
      for(register int t=1;t<maxn;++t)
	    jc[t]=1ll*jc[t-1]*t%mod;
      inv[maxn-1]=ksm(jc[maxn-1],mod-2);
      for(register int t=maxn-2;t;--t)
	    inv[t]=1ll*inv[t+1]*(t+1)%mod;
      for(register int t=0;t<maxn;++t)
	    if(1ll*inv[t]*jc[t]%mod!=1) puts("WA");
}

inline int c(const int&n,const int&m){
      if(n<m)return 0;
      return 1ll*jc[n]*inv[m]%mod*inv[n-m]%mod;
}

int n,m,s,L;
int f[maxn];
int g[maxn];

int main(){
#ifndef ONLINE_JUDGE
      freopen("in.in","r",stdin);
      
#endif
      pre();
      n=qr();m=qr();s=qr();
      L=min(n/s,min(m,n));
      for(register int t=0;t<=L;++t){
	    f[t]=1ll*c(m,t)*jc[n]%mod*ksm(m-t,n-s*t)%mod*inv[n-s*t]%mod*ksm(ksm(jc[s],t),mod-2)%mod*jc[t]%mod;
      }
      for(register int t=0;t<=L;++t){
	    g[t]=inv[L-t];
	    if((L-t)&1)g[t]=mod-g[t];
      }
      int k=1;
      while(k<=L)k<<=1;
      NTT(f,k<<1,1);
      NTT(g,k<<1,1);
      for(register int t=0;t<k<<1;++t) g[t]=1ll*g[t]*f[t]%mod;
      NTT(g,k<<1,-1);
      int ans=0;
      for(register int t=0;t<=L;++t) ans=(ans+1ll*qr()*g[t+L]%mod*inv[t]%mod)%mod;
      printf("%d\n",ans);
      return 0;
}