【UOJ #390】【UNR #3】百鸽笼(指数型生成函数,二项式定理)

博客探讨了一个概率与组合优化问题,涉及指数型生成函数(EGF)和动态规划(DP)的应用。内容包括将问题转化为求特定序列的概率,使用EGF计算合法序列,并通过背包DP处理每个位置成为放满人的最后一列的概率。最终,通过广义二项式展开简化计算,实现算法的优化,避免超时。

摘要生成于 C知道 ,由 DeepSeek-R1 满血版支持, 前往体验 >

题面

🔗
在这里插入图片描述

限制

在这里插入图片描述

题解

参考

我们将题目等价为:每次随机任意一个位置放人,当人数 ≥ a i \geq a_i ai 时认为 i i i 放满,一直放,询问 i i i 成为最后放满人的一列的概率。

我们限定 I I I 为最后一列,那么其他列就必须放满,合法的放人序列(到 I I I 填满为止)可以用指数型生成函数(EGF)计算。(放心,不会很难理解的)

我们为了方便表示,令
G k ( x ) = ∑ i = 0 a k − 1 x i i ! G_k(x)=\sum_{i=0}^{a_k-1} \frac{x^i}{i!} Gk(x)=i=0ak1i!xi

其他列的 EGF 就是 e x − G k ( x ) e^x-G_k(x) exGk(x) ,第 I I I 列则是 x a k − 1 ( a k − 1 ) ! \frac{x^{a_k-1}}{(a_k-1)!} (ak1)!xak1 ,那么 I I I除去最后一项的合法序列 EGF 就是
F ( x ) = x a k − 1 ( a k − 1 ) ! ∏ i ≠ I ( e x − G i ( x ) ) F(x)=\frac{x^{a_k-1}}{(a_k-1)!}\prod_{i\not=I}(e^x-G_i(x)) F(x)=(ak1)!xak1i=I(exGi(x))

I I I 列的答案就是所有合法序列的概率之和(由于除去了最后一项, 1 n \frac{1}{n} n1 的幂为长度+1):
∑ i = 0 ∞ i ! ⋅ [ x i ] F ( x ) n i + 1 \sum_{i=0}^{\infty} \frac{i!\cdot[x^i]F(x)}{n^{i+1}} i=0ni+1i![xi]F(x)

由于 e x e^x ex 长度是无穷的,我们就把 F ( x ) F(x) F(x) 当成 x , e x x,e^x x,ex 的二元生成函数,处理出每个 x i ( e x ) j x^i(e^x)^j xi(ex)j 的系数 f i , j f_{i,j} fi,j ,这个可以背包DP做。然后我们单独算每个 f i , j x i ( e x ) j f_{i,j}x^i(e^x)^j fi,jxi(ex)j 的贡献:
f i , j n i + 1 ∑ d = 0 ∞ ( d + i ) ! j d d ! n d = f i , j ⋅ i ! n i + 1 ∑ d = 0 ∞ ( j n ) d ( d + i i ) \frac{f_{i,j}}{n^{i+1}}\sum_{d=0}^{\infty}(d+i)!\frac{j^d}{d!n^{d}}\\ =\frac{f_{i,j}\cdot i!}{n^{i+1}}\sum_{d=0}^{\infty}(\frac{j}{n})^d{d+i\choose i} ni+1fi,jd=0(d+i)!d!ndjd=ni+1fi,ji!d=0(nj)d(id+i)

右边那个形式感觉很有规律,其实它是广义二项式展开
( x + y ) − n = ∑ i = 0 ∞ ( − 1 ) i ( n − 1 + i n − 1 ) x i y − n − i (x+y)^{-n}=\sum_{i=0}^{\infty}(-1)^i{n-1+i\choose n-1}x^{i}y^{-n-i} (x+y)n=i=0(1)i(n1n1+i)xiyni

所以上式可以转化为
f i , j ⋅ i ! n i + 1 ∑ d = 0 ∞ ( − j n ) d ( − 1 ) d ( d + i i ) = f i , j ⋅ i ! n i + 1 ( 1 − j n ) − ( i + 1 ) = f i , j ⋅ i ! n i + 1 ( n n − j ) i + 1 = f i , j ⋅ i ! ( n − j ) i + 1 \frac{f_{i,j}\cdot i!}{n^{i+1}}\sum_{d=0}^{\infty}(-\frac{j}{n})^d(-1)^d{d+i\choose i}\\ =\frac{f_{i,j}\cdot i!}{n^{i+1}}(1-\frac{j}{n})^{-(i+1)}\\ =\frac{f_{i,j}\cdot i!}{n^{i+1}}(\frac{n}{n-j})^{i+1}\\ =\frac{f_{i,j}\cdot i!}{(n-j)^{i+1}}\\ ni+1fi,ji!d=0(nj)d(1)d(id+i)=ni+1fi,ji!(1nj)(i+1)=ni+1fi,ji!(njn)i+1=(nj)i+1fi,ji!

如此简洁。

但是,我们对于每个 I I I 都求一次背包DP刚好要超时。由于DP过程只有加减法,我们可以先处理出所有列的背包DP,然后再去掉 I I I (假想把 I I I 放到最后,按照DP的转移反着做一次)。这样复杂度就对了, O ( n 3 ∑ a i ) O(n^3\sum a_i) O(n3ai)

CODE

#include<map>
#include<set>
#include<cmath>
#include<ctime>
#include<queue>
#include<stack>
#include<random>
#include<bitset>
#include<vector>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
#include<unordered_map>
// #pragma GCC optimize(2)
using namespace std;
#define MAXN 905
#define LL long long
#define ULL unsigned long long
#define ENDL putchar('\n')
#define DB double
#define lowbit(x) (-(x) & (x))
#define FI first
#define SE second
#define PR pair<int,int>
#define UIN unsigned int
int xchar() {
	static const int maxn = 1000000;
	static char b[maxn];
	static int pos = 0,len = 0;
	if(pos == len) pos = 0,len = fread(b,1,maxn,stdin);
	if(pos == len) return -1;
	return b[pos ++];
}
// #define getchar() xchar()
LL read() {
	LL f = 1,x = 0;int s = getchar();
	while(s < '0' || s > '9') {if(s<0)return -1;if(s=='-')f=-f;s = getchar();}
	while(s >= '0' && s <= '9') {x = (x<<1) + (x<<3) + (s^48);s = getchar();}
	return f*x;
}
void putpos(LL x) {if(!x)return ;putpos(x/10);putchar((x%10)^48);}
void putnum(LL x) {
	if(!x) {putchar('0');return ;}
	if(x<0) putchar('-'),x = -x;
	return putpos(x);
}
void AIput(LL x,int c) {putnum(x);putchar(c);}

const int MOD = 998244353;
int n,m,s,o,k;
int qkpow(int a,int b) {
	int res = 1;
	while(b > 0) {
		if(b & 1) res = res *1ll* a % MOD;
		a = a *1ll* a % MOD; b >>= 1;
	} return res;
}
int fac[MAXN],inv[MAXN],invf[MAXN];
int A[MAXN];
int dp[35][MAXN],dp2[35][MAXN];
int main() {
	n = read();
	fac[0]=fac[1]=inv[0]=inv[1]=invf[0]=invf[1]=1;
	for(int i = 2;i <= 900;i ++) {
		fac[i] = fac[i-1] *1ll* i % MOD;
		inv[i] = (MOD - inv[MOD%i]) *1ll* (MOD/i) % MOD;
		invf[i] = invf[i-1] *1ll* inv[i] % MOD;
	}
	int sm = 0;
	dp[0][0] = 1;
	for(int i = 1;i <= n;i ++) {
		A[i] = read();
		sm += A[i]-1;
		for(int j = i;j > 0;j --) {
			for(int k = sm;k >= 0;k --) {
				for(int s = 0;s <= k && s < A[i];s ++)
				(dp[j][k] += MOD - dp[j-1][k-s]*1ll*invf[s] % MOD) %= MOD;
			}
		}
	}
	for(int i = 1;i <= n;i ++) {
		memcpy(dp2,dp,sizeof(dp));
		for(int j = 1;j <= n;j ++) {
			for(int k = 0;k <= sm;k ++) {
				for(int s = 0;s <= k && s < A[i];s ++)
				(dp2[j][k] += dp2[j-1][k-s]*1ll*invf[s] % MOD) %= MOD;
			}
		}
		for(int j = n;j > 0;j --) {
			for(int k = A[i]-1;k <= sm;k ++) {
				(dp2[j][k] = dp2[j-1][k-A[i]+1]*1ll*invf[A[i]-1] % MOD) %= MOD;
			}
		}
		int as = 0;
		for(int j = 0;j < n;j ++) {
			int pw = qkpow(inv[n-j],A[i]-1);
			for(int k = A[i]-1;k <= sm;k ++) {
				pw = pw *1ll* inv[n-j] % MOD;
				(as += pw *1ll* dp2[n-j][k] % MOD * fac[k] % MOD) %= MOD;
			}
		}
		AIput(as,i==n ? '\n':' ');
	}
	return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值