HAOI2018 染色

传送门

一道非常好的容斥+NTT，对我这样的菜鸡难度稍高。
符合要求的颜色最多有(lim = min(m,lfloorfrac{n}{S} floor))种。
首先，考虑恰好出现S次不是很容易，那么我们换一种想法，先考虑至少出现S次。设(f[i])表示至少有i种颜色出现S次的情况数。
首先要从m种颜色选i种，这是(C_m^i)种情况。
之后我们考虑取出的位置，i种颜色占用的位置是iS个，这是(C_n^{iS})种情况。
再考虑取出iS个位置中颜色的 摆放情况。因为颜色是相同的，所以应该是无标号的全排列，也就是(frac{(iS!)}{(S!)^i})种情况。
最后剩下的位置可以随便乱染，就是((m-i)^{n-iS})种情况。

用乘法原理乘起来，有(f[i] = C_m^iC_n^{iS}frac{(iS!)}{(S!)^i}(m-i)^{n-iS})

之后我们考虑(ans[i])表示恰好有i种颜色出现S次的情况。
那么根据容斥原理可以想得到，(ans[i] = sum_{j=i}^{lim}(-1)^{j-i}C_j^if[j])
把组合数拆开，把(i!)移到等式左边，(ans[i]*i! = sum_{j=i}^{lim}frac{(-1)^{j-i}}{(j-i)!}j!f[j])
我们设(F(x) = frac{-1^i}{i!}x^i,G(x) = f[i]i!x^i)
然后这个就已经非常像一个卷积的形式了……我们考虑把G反转就可以卷积……
之后再把卷积出来得到的答案数组反转过来对应的就是答案了。
最后对应的乘起来累加即为答案。
注意这道题预处理阶乘逆元的时候，要处理到(max(n,m))，否则可能会导致访问到没有处理的位置。
出题人比较良心的给了NTT模数，原根是3.

#include<bits/stdc++.h>
#define rep(i,a,n) for(int i = a;i <= n;i++)
#define per(i,n,a) for(int i = n;i >= a;i--)
#define enter putchar('
')
#define fr friend inline
#define y1 poj
#define mp make_pair
#define pr pair<int,int>
#define fi first
#define sc second
#define pb push_back

using namespace std;
typedef long long ll;
const int M = 2e5+5;
const int N = 2e7+5;
const int mod = 1004535809;
const int G = 3;
const int invG = 334845270;
const double eps = 1e-7;

int read()
{
   int ans = 0,op = 1;char ch = getchar();
   while(ch < '0' || ch > '9') {if(ch == '-') op = -1;ch = getchar();}
   while(ch >= '0' && ch <= '9') ans = ans * 10 + ch - '0',ch = getchar();
   return ans * op;
}

int n,m,s,w[M],fac[N],inv[N],rev[M<<2],F[M<<2],g[M<<2],ans[M<<2],f[M<<2],lim,tot;

int inc(int a,int b){return (a+b) % mod;}
int mul(int a,int b){return 1ll * a * b % mod;}
int qpow(int a,int b)
{
   int p = 1;
   while(b)
   {
      if(b & 1) p = mul(p,a);
      a = mul(a,a),b >>= 1;
   }
   return p;
}

void NTT(int *a,int l,int f)
{
   rep(i,0,l-1) if(i < rev[i]) swap(a[i],a[rev[i]]);
   for(int i = 1;i < l;i <<= 1)
   {
      int w1 = qpow(f ? G : invG,(mod-1) / (i<<1));
      for(int j = 0;j < l;j += (i<<1))
      {
	 int w = 1;
	 rep(k,0,i-1)
	 {
	    int kx = a[k+j],ky = mul(a[k+j+i],w);
	    a[k+j] = inc(kx,ky),a[k+j+i] = inc(kx,mod-ky),w = mul(w,w1);
	 }
      }
   }
   if(!f)
   {
      int inv = qpow(l,mod-2);
      rep(i,0,l-1) a[i] = mul(a[i],inv);
   }
}

void init()
{
   int k = max(n,m);
   fac[0] = inv[0] = inv[1] = 1;
   rep(i,1,k) fac[i] = mul(fac[i-1],i);
   inv[k] = qpow(fac[k],mod-2);
   per(i,k-1,1) inv[i] = mul(inv[i+1],i+1);
}

int C(int n,int m)
{
   if(n < m) return 0;
   return mul(mul(fac[n],inv[m]),inv[n-m]);
}

int main()
{
   n = read(),m = read(),s = read(),init(),lim = min(m,n/s);
   rep(i,0,m) w[i] = read();
   rep(i,0,lim) f[i] = mul(mul(mul(C(m,i),C(n,i*s)),mul(fac[i*s],qpow(inv[s],i))),qpow(m-i,n-i*s));
   g[0] = 1,F[0] = f[0];
   rep(i,1,lim)
   {
      F[i] = mul(f[i],fac[i]);
      (i & 1) ? g[i] = mod - inv[i] : g[i] = inv[i];
   }
   int l = 1,L = 0;
   while(l <= lim<<1) l <<= 1,L++;
   rep(i,0,l-1) rev[i] = (rev[i>>1] >> 1) | ((i&1) << (L-1));
   reverse(F,F+lim+1);
   NTT(F,l,1),NTT(g,l,1);
   rep(i,0,l-1) F[i] = mul(F[i],g[i]);
   NTT(F,l,0);
   reverse(F,F+lim+1);
   rep(i,0,lim) tot = inc(tot,mul(mul(F[i],inv[i]),w[i]));
   printf("%d
",tot);
   return 0;
}