概率论与数理统计

概率dp

校内赛 概率论与数理统计

题目描述：

这天空洞大师zzk在教gyx打隔膜，考虑到gyx的手残体质，大师决定教他刮痧流打法。

所谓刮痧流即打一次就跑，具体的,gyx每次攻击可以使选定的目标生命值-1，可以攻击m次。

教程关一共有n个怪物，第i个生命值为$a_i$，都没有攻击力，排成一排等着挨打。

但是由于不会操作,gyx还是只会乱砍，所以每次行动他会随机选一个当前生命值为正数的怪物攻击。

作为概率论大师的柴老师觉得这个过程很有趣，于是想了一道题：经过m次攻击后,gyx能击败的怪物期望数。

称一个怪物被击败，当且仅当某一时刻其生命值=0。

题目分析：

乍一看好像每一次让一个怪物扣血的概率都和场上存活的怪物数量有关系，又注意到n很小，容易想到使用状态压缩。记 $dp[s][i]$ 为攻击了 $i$ 次，有 $s$ 这个状态的怪物死了的概率。那么最后的答案就是攻击m次之后每个状态的popcount乘以概率乘以其余怪物活着的方案。记 $g[s][i]$ 为砍了 $i$ 次，有 $s$ 这些怪物活着的方案数。那么可以分别得到 $dp$ 和 $g$ 的状态转移方程：

如果第$i+1$刀砍死了一个怪物，那么 $dp[s][i]->dp[s|(1<<j)][i+1]$ ，其概率为 {\bionm{a[j]-1} {i - sum[s]}}/ popcount(s) 。因为要确保在i+1杀死这个怪物，那么第i刀已经杀了 $a[j]-1$ 刀了。如果没砍死，那么就直接转移到 $dp[s][i+1]$

考虑如何计算 $g[s][i]$ ,这个就根背包类似。

代码

//
// Created by mrx on 2022/10/4.
//
#include <vector>
#include <algorithm>
#include <iostream>
#include <array>
#include <cmath>
#include <iomanip>


int lowbit(int x) {
	return x & -x;
}

int main() {
	std::ios::sync_with_stdio(false);
	std::cin.tie(nullptr);
	std::cout.tie(nullptr);

	int n, m;
	std::cin >> n >> m;

	std::vector<std::vector<long double>> C(105, std::vector<long double>(105));
	C[0][0] = 1;
	for (int i = 1; i <= 100; ++i) {
		C[i][0] = 1;
		for (int j = 1; j <= 100; ++j) {
			C[i][j] = C[i - 1][j] + C[i - 1][j - 1];
		}
	}

	std::vector<int> a(n);
	for (int i = 0; i < n; ++i)std::cin >> a[i];

	std::vector<int> sum(1 << n);
	for (int s = 0; s < 1 << n; ++s) for (int i = 0; i < n; ++i) if (s >> i & 1)sum[s] += a[i];

	std::vector<std::vector<long double>> g(m + 1, std::vector<long double>(1 << n));
	//噶i刀，存活状态为s的方案数量。
	g[0][0] = 1;
	for (int s = 0; s < 1 << n; ++s) {
		for (int k = 0; k < n; ++k) {
			if ((s >> k) & 1) {
				int s0 = s ^ (1 << k);
				for (int i = 0; i <= m; ++i) {
					for (int j = 0; j <= std::min(a[k] - 1, i); ++j) {
						g[i][s] += g[i - j][s0] * C[i][j];
					}
				}
				break;
			}
		}
	}


	std::vector<std::vector<long double>> dp(m + 1, std::vector<long double>(1 << n));
	//噶i刀，状态为s的挂了的概率
	dp[0][0] = 1;
	for (int s = 0; s < 1 << n; ++s) {
		int cur = __builtin_popcount(s);
		int cnt = n - cur;
		for (int now = sum[s]; now < m; ++now) {
			int rem = ((1 << n) - 1) ^ s;
			for (int j = 0; j < n; ++j) {
				if (!((s >> j) & 1)) {
					if (a[j] - 1 <= now - sum[s])dp[now + 1][s | (1 << j)] += dp[now][s] * C[now - sum[s]][a[j] - 1] / cnt;
				}
			}
			dp[now + 1][s] += dp[now][s] / cnt;
		}
	}

	for (int s = 0; s < 1 << n; ++s) {
		for (int i = 0; i <= m; ++i) {
			if (i < sum[s])continue;
			int rem = ((1 << n) - 1) ^ s;
			dp[i][s] *= g[i - sum[s]][rem];
		}
	}

	long double ans = 0;
	for (int s = 0; s < 1 << n; ++s) {
		if (m < sum[s])continue;
		int cnt = __builtin_popcount(s);
		ans += dp[m][s] * cnt;
	}
	std::cout << std::fixed << std::setprecision(5) << ans << '\n';
	return 0;
}