【题解】[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 floor) ，然后限制一下(Lle n ,Lle m)
先有关系式
[f(x)=sum_{i=x}^L {ichoose x}g(i)
]
(恰好出现多次，会有({ichoose x})被钦定的可能。)
再容斥(或者叫二项式反演?)得到
[g(x)=sum_{i=x}^L (-1)^{i-x}{ichoose x}f(i)
]
解释一波这个容斥的含义
由于对于每个有(u)种颜色出现(S)的方案，有三种状态：
-
当(u)少于目标(x)时：在(f(x))中被计算了(0)次。
-
当(u)等于目标(x)时：在(f(x))中只被计算了(1)次。
-
当(u)大于目标(x)时：
考虑这个方案，在(f(x))到(f(u))之中分别计算了
[{uchoose x},{uchoose x+1},dots{uchoose u} ]次，很像一个二项式定理展开。
然后我们考虑一下({x+ichoose x}f(x+i))的情况
[{uchoose x}{xchoose x},{uchoose x+1}{x+1choose x},dots{uchoose u}{uchoose x} ]晓得一个公式，那就是
[{achoose b}{bchoose c}={achoose c}{a-cchoose b-c} ]用这个公式化简一下上面那个式子
[sum_{i=0}^{u-x}{uchoose x}{u-xchoose i} ]好家伙，右边就是一个二项式定理啊！
继续化简
[={uchoose x}sum_{i=0}^{u-x}{u-xchoose i}={uchoose x}(1+1)^{u-x} ]现在我们要让等式右边变成
[0={uchoose x} imes 0={uchoose x}(1-1)^{u-x} ]我们在(f(i))前面乘上一个容斥系数即可。此时，(ule x)的情况不受影响。
故
[sum_{i=x}^L (-1)^{i-x}{ichoose x}f(i)=g(x) ]
然后我们把(g(x))的式子化简一下
[g(x)=sum_{i=x}^L (-1)^{i-x}{ichoose x}f(i)=sum_{i=x}^L (-1)^{i-x}dfrac {i!} {(i-x)!x!}f(i)
\
herefore x!g(x)=sum_{i=x}^Ldfrac {(-1)^{i-x}}{(i-x)!} imes i!f(i)
]
拿出来
[x!g(x)=sum_{i=x}^Ldfrac {(-1)^{i-x}}{(i-x)!} imes 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
",ans);
return 0;
}