P4491 [HAOI2018]染色
反思
二项式反演的板子,但是还要加上一个NTT,就变得恶心人了
我NTT写挂调了好长时间。。。
注意
NTT应该这么写
int tmp=pow((opt)?invG:G,(MOD-1)/i);
而不这样写(因为单位元是\(G^{\frac{mod-1}{i}}\))
int tmp=pow((opt)?invG:G,n/i);
这里应该这样写
for(int i=2;i<=n;i<<=1){
...
}
而不应该这样写(每次长度要倍增)
for(int i=2;i<=n;i++){
...
}
多项式还是写的少了
思路
\(\{x_1*y_1,x_2*y_2,x_3*y_3,\dots,x_k*y_k\}\)的一个多重集合,设大小为n,它的排列数是
\[\frac{n!}{x_1!x_2!x_3!\dots x_4!} \]
题目让求的是有恰好出现s次的颜色个数恰好是i个的方案数
考虑容斥一发
设\(F_i\)表示出现s次的颜色个数最少是i个的方案数
考虑钦定i个
则有
\[F_i=\left(\begin{matrix}m\\i\end{matrix}\right)(m-i)^{n-i\times s}\frac{n!}{(s!)^i(n-s*i)!}
\]
从前到后分别是选哪i个,后面的\(m-i\)种颜色随便选\(n-i\times s\)个位置,剩下的是一个多重集合的排列
设\(G_i\)表示出现s次的颜色个数恰好是i个的方案
则有
\[F_i=\sum_{j=i}^m\left(\begin{matrix}j\\i\end{matrix}\right)G_j
\]
二项式反演一下,有
\[G_i=\sum_{j=i}^n(-1)^{j-i}\left(\begin{matrix}j\\i\end{matrix}\right)F_j
\]
把组合数拆开
\[G_i=\sum_{j=i}^n(-1)^{j-i}\frac{j!}{(j-i)!(i)!}F_j
\]
可以得到
\[i!G_i=\sum_{j=i}^n\frac{(-1)^{j-i}F_jj!}{(j-i)!}
\]
设\(f_i=F_ii!\),\(g_i=\frac{(-1)^{j-i}}{(j-i)!}\)
则
\[i!G_i=\sum_{j=i}^nf_ig_{j-i}
\]
是一个卷积的形式,按照力的做法,把\(f\)反向成\(f'\),卷积之后的第\(lim-i\)项就对应原来的第\(i\)项
上NTT就好了
代码
#include <cstdio>
#include <algorithm>
#include <cstring>
#define int long long
using namespace std;
const int G = 3,MOD = 1004535809 , invG = 334845270;
int pow(int a,int b){
int ans=1;
while(b){
if(b&1)
ans=(ans*a)%MOD;
a=(a*a)%MOD;
b>>=1;
}
return ans;
}
void NTT(int *a,int n,int opt){
int lim=0;
while((1<<lim)<n)
lim++;
n=(1<<lim);
for(int i=0;i<n;i++){
int t=0;
for(int j=0;j<lim;j++)
if((i>>j)&1)
t|=(1<<(lim-j-1));
if(i<t)
swap(a[i],a[t]);
}
for(int i=2;i<=n;i<<=1){
int len=i/2;
int tmp=pow((opt)?invG:G,(MOD-1)/i);
for(int j=0;j<n;j+=i){
int mid=1;
for(int k=j;k<len+j;k++){
int t=mid*a[k+len]%MOD;
a[k+len]=((a[k]-t)%MOD+MOD)%MOD;
a[k]=(a[k]+t)%MOD;
mid=(mid*tmp)%MOD;
}
}
}
if(opt){
int invn=pow(n,MOD-2);
for(int i=0;i<n;i++)
a[i]=(a[i]*invn)%MOD;
}
}
int a[20001000],b[20001000],Gx[20001000],f[20001000],jc[20001000],jc_inv[20001000],n,m,s,w[20001000];
int C(int n,int m){
if(n<m)
return 0;
return jc[n]*jc_inv[n-m]%MOD*jc_inv[m]%MOD;
}
signed main(){
// printf("INVG=%lld\n",pow(G,MOD-2));
scanf("%lld %lld %lld",&n,&m,&s);
for(int i=0;i<=m;i++)
scanf("%lld",&w[i]);
int lim=min(n/s,m),maxx=max(n,m);
jc[0]=1;
jc_inv[0]=1;
for(int i=1;i<=maxx;i++)
jc[i]=(jc[i-1]*i)%MOD;
jc_inv[maxx]=pow(jc[maxx],MOD-2);
for(int i=maxx-1;i>=1;i--)
jc_inv[i]=(jc_inv[i+1]*(i+1))%MOD;
for(int i=0;i<=min(n/s,m);i++){
f[i]=(C(m,i)*jc[n]%MOD*pow(pow(jc[s],i),MOD-2)%MOD*(pow(jc[n-i*s],MOD-2))%MOD*pow(m-i,n-i*s)%MOD)%MOD;
// printf("f[%lld]=%lld\n",i,f[i]);
}
for(int i=0;i<=lim;i++){
a[i]=f[i]*jc[i]%MOD;
b[i]=(((i&1)?-1:1)*jc_inv[i]%MOD+MOD)%MOD;
// b[i]=jc_inv[lim-i];
// if((lim-i)&1)
// b[i]=MOD-b[i];
// printf("A[%lld]=%lld B[%lld]=%lld\n",i,a[i],i,b[i]);
}
for(int i=0,j=lim;i<j;i++,j--)
swap(a[i],a[j]);
int ttt=0;
while((1<<ttt)<(2*lim+4))
ttt++;
int lent=(1<<ttt);
NTT(a,lent,0);
NTT(b,lent,0);
for(int i=0;i<lent;i++)
Gx[i]=(a[i]*b[i])%MOD;
// for(int i=0;i<lent;i++)
// printf("A[%lld]=%lld B[%lld]=%lld\n",i,a[i],i,b[i]);
NTT(Gx,lent,1);
int ans=0;
for(int i=0;i<=lim;i++)
ans=(ans+Gx[lim-i]*jc_inv[i]%MOD*w[i]%MOD)%MOD;
printf("%lld\n",ans);
return 0;
}