[LuoguP4859]P4859 已经没有什么好害怕的了(二项式反演)

[LuoguP4859]P4859 已经没有什么好害怕的了(二项式反演)

题面

已经使 Madoka 有签订契约，和自己一起战斗的想法后，Mami 忽然感到自己不再是孤单一人了呢。

于是，之前的谨慎的战斗作风也消失了，在对 Charlotte 的傀儡使用终曲——Tiro Finale 后，Mami 面临着即将被 Charlotte 的本体吃掉的局面。

这时，已经多次面对过 Charlotte 的 Homura 告诉了学 OI 的你这样一个性质：Charlotte 的结界中有两种具有能量的元素，一种是“糖果”，另一种是“药片”，各有 n 个。在 Charlotte 发动进攻前，“糖果”和“药片”会两两配对，若恰好糖果比药片能量大的组数比“药片”比“糖果”能量大的组数多 k组，则在这种局面下，Charlotte 的攻击会丟失，从而 Mami 仍有消灭 Charlotte 的可能。

你必须根据 Homura 告诉你的“糖果”和“药片”的能量的信息迅速告诉 Homura 这种情况的个数.

分析

首先问题转化为\(a>b\)的个数有\(\frac{n+k}{2}\)

设\(f(i)\)表示保证(钦定)有\(i\)对\(a>b\),其余不一定。\(g(i)\)表示恰好\(i\)对a>b.显然\(f(i)=\sum_{j=i}^m C_j^i g(j)\),根据二项式反演定理\(g(i)=\sum_{j=i}^m C_{j}^i f(j)\)

考虑如何求\(f\).把\(a\)从小到大排序,记\(cnt(i)\)表示\(b\)中小于\(a_i\)的数的个数。设\(dp_{i,j}\)为\(a\)序列的前\(i\)个,保证有\(j\)对\(a>b\)的方案数,则\(f_i=dp_{n,i}(n-i)!\)

容易写出转移方程：
\(dp_{i,j}=dp_{i,j-1}+dp_{i-1,j-1}\cdot (cnt(i)-(j-1))\)

我们按\(a\)排序,因此\(b\)中小于\(a_i\)的数的集合一定包含\(b\)中小于\(a_{i-1}\)的数的集合.之前已经选了\(j-1\)个数,那只能选剩下的\(cnt(i)-(j-1)\)个数来使得\(a_i>b_x\),这样对数会增加1。

代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#define mod 1000000009
#define maxn 2000
using namespace std;
typedef long long ll;
inline ll fast_pow(ll x,ll k){
	ll ans=1;
	while(k){
		if(k&1) ans=ans*x%mod;
		x=x*x%mod;
		k>>=1;
	}
	return ans;
}
inline ll inv(ll x){
	return fast_pow(x,mod-2);
}
ll fact[maxn+5],invfact[maxn+5];
void ini(int n){
	fact[0]=1;
	for(int i=1;i<=n;i++) fact[i]=fact[i-1]*i%mod;
	invfact[n]=inv(fact[n]);
	for(int i=n-1;i>=0;i--) invfact[i]=invfact[i+1]*(i+1)%mod;
}
ll C(int n,int m){
	return fact[n]*invfact[n-m]%mod*invfact[m]%mod;
} 

int n,k;
int a[maxn+5],b[maxn+5];
int cnt[maxn+5];//a从小到大排序后,a[i]>b[j]的j的数量 
ll dp[maxn+5][maxn+5];
ll f[maxn+5],g[maxn+5];

void calc_inv(ll *f,ll *g,int n){
	for(int i=0;i<=n;i++){
		for(int j=i;j<=n;j++){
			if((j-i)%2==1) g[i]-=C(j,i)*f[j]%mod;
			else g[i]+=C(j,i)*f[j]%mod;
			g[i]=(g[i]+mod)%mod;
		}
	}
} 
int main(){
	scanf("%d %d",&n,&k);
	ini(n);
	for(int i=1;i<=n;i++) scanf("%d",&a[i]);
	for(int i=1;i<=n;i++) scanf("%d",&b[i]);
	sort(a+1,a+1+n);
	for(int i=1;i<=n;i++){
		for(int j=1;j<=n;j++) if(a[i]>b[j]) cnt[i]++;
	}
	for(int i=0;i<=n;i++) dp[i][0]=1;
	for(int i=1;i<=n;i++){
		for(int j=1;j<=n;j++){
			dp[i][j]=(dp[i-1][j-1]*max(cnt[i]-j+1,0)%mod+dp[i-1][j])%mod;
		}
	}
	for(int i=0;i<=n;i++) f[i]=dp[n][i]*fact[n-i]%mod;
	calc_inv(f,g,n);
	printf("%lld\n",g[(n+k)/2]);
}