算法动态规划

分汤-概率DP与浮点数精度处理结合

字数统计: 1.8k字   |   阅读时长: 7分
,

摘要: 一道概率 DP 的问题，记忆化搜索

【对算法，数学，计算机感兴趣的同学，欢迎关注我哈，阅读更多原创文章】
我的网站：潮汐朝夕的生活实验室
我的公众号：算法题刷刷
我的知乎：潮汐朝夕
我的github：FennelDumplings
我的leetcode：FennelDumplings

各位好，前面我们详细学习过概率 DP 这类算法，参考文章 概率DP。【力扣808. 分汤】是力扣上比较少见的一个相关题目，本文我们就来看一下。

本题有三个要点，一个是含概率的动态规划的分析方法，一个是对浮点数精度的处理，一个是使用记忆化搜索避免无效状态的计算。

题目

808. 分汤

有 A 和 B 两种类型 的汤。一开始每种类型的汤有 n 毫升。有四种分配操作：

提供 100ml 的 汤A 和 0ml 的 汤B 。
提供 75ml 的 汤A 和 25ml 的 汤B 。
提供 50ml 的 汤A 和 50ml 的 汤B 。
提供 25ml 的 汤A 和 75ml 的 汤B 。

当我们把汤分配给某人之后，汤就没有了。每个回合，我们将从四种概率同为 0.25 的操作中进行分配选择。如果汤的剩余量不足以完成某次操作，我们将尽可能分配。当两种类型的汤都分配完时，停止操作。

注意 不存在先分配 100 ml 汤B 的操作。

需要返回的值： 汤A 先分配完的概率 + 汤A和汤B 同时分配完的概率 / 2。返回值在正确答案 1e-5 的范围内将被认为是正确的。

提示:

1
0 <= n <= 1e9

示例 1:
输入: n = 50
输出: 0.62500
解释:如果我们选择前两个操作，A 首先将变为空。
对于第三个操作，A 和 B 会同时变为空。
对于第四个操作，B 首先将变为空。
所以 A 变为空的总概率加上 A 和 B 同时变为空的概率的一半是 0.25 *(1 + 1 + 0.5 + 0)= 0.625。

示例 2:
输入: n = 100
输出: 0.71875

题解

算法1: 概率 DP

为了方便，我们将 25ml 视为一单位。这样每个回合操作便有 4 种操作：

  • 4 单位 A、0 单位 B
  • 3 单位 A、1 单位 B
  • 2 单位 A、2 单位 B
  • 1 单位 A、3 单位 B

相应地，一开始每种类型的汤有 n 毫升，相当于 N = ceil(n / 25.0) 单位。

在每一回合后，我们考察 A, B 两种汤还剩余的单位数，记为 i, j，有 4 种情况：

  • i = 0, j > 0: A 分完了，B 还没分完 (A 先分配完)
  • i > 0, j = 0: B 分完了，A 还没分完 (B 先分配完)
  • i = 0, j = 0: A, B 同时分完了 (A 和 B 同时分配完)
  • i > 0, j > 0: A, B 均未分完

如果是前三种则停止操作，第四种则进行下一回合操作。

如果将 (i, j) 视为一个状态，那么基于一回合可选的四种操作，当前回合后的状态 (i, j) 可以从前一回合的以下四种状态得到：

1
2
3
4
(i - 4, j)
(i - 3, j - 1)
(i - 2, j - 2)
(i - 1, j - 3)

(i, j) 状态来自以上四种状态的概率均为 0.25，记 dp[i][j] := A, B 两种汤剩余的单位数为 (i, j) 时，对应的最终所求的概率值，根据全概率公式有：

1
dp[i][j] = 0.25 * (dp[i - 4][j] + dp[i - 3][j - 1] + dp[i - 2][j - 2] + dp[i - 1][j - 3])

由此写出概率 DP 算法如下，其中每一步操作是一个阶段，取出的 A, B 两种汤的数量 i, j 是阶段的附加状态，也就是 i, j 隐含了阶段:

1
2
3
4
5
6
7
8
9
10
11
12
13
状态定义：
dp[i][j] := A, B 两种汤剩余的单位数为 (i, j) 时，对应的最终所求的概率值

答案：
dp[N][N]

初始化：
dp[0][0] = 0.5
dp[0][1..N] = 1.0
dp[1..N][0] = 0.0

状态转移：
dp[i][j] = 0.25 * (dp[i - 4][j] + dp[i - 3][j - 1] + dp[i - 2][j - 2] + dp[i - 1][j - 3])

这个动态规划算法的时间复杂度是 $O(N^{2})$，但是本题中 N 很大，因此还要进一步处理。

根据每一回合的四种可能的操作，我们看一下两种汤被分配的单位数的期望，汤 A 被分配的期望为 (1 + 2 + 3 + 4) * 0.25 = 2.5；汤 B 被分配的期望为 (0 + 1 + 2 + 3) * 0.25 = 1.5。

因此我们很自然地推测，在 n 很大的时候，汤 A 会有很大的概率比 B 先分配完，且这个概率非常接近 1，并且 n 越大，越接近 1。由于本题要求答案在 1e-5 的范围就是正确的，因此当 n 大到某个值，使得 dp[n][n] > 0.99999 时，直接返回 1.0 作为答案即可。

现在的问题是前面提到的某个值具体是多少，可以通过尝试 N 一些 取值，然后算出 dp[N][N]，看是否大于 0.99999，这个值不需要太准确，能解决问题就行。这里我们取 N = 500。

代码 (C++)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
class Solution {
public:
    double soupServings(int N) {
        N = ceil(N / 25.0);
        if(N >= 500)
            return 1.0;

        vector<vector<double> > dp(N + 1, vector<double>(N + 1, 0));
        dp[0][0] = 0.5;
        for(int j = 1; j <= N; ++j)
            dp[0][j] = 1.0;
        for(int i = 1; i <= N; ++i)
            for(int j = 1; j <= N; ++j)
            {
                dp[i][j] = 0.25 * (dp[relu(i - 4)][j] + dp[relu(i - 3)][j - 1]
                           + dp[relu(i - 2)][relu(j - 2)] + dp[relu(i - 1)][relu(j - 3)]);
            }
        return dp[N][N];
    }

private:
    int relu(int x)
    {
        return max(0, x);
    }
};

算法2: 记忆化搜索

观察状态转移方程，我们可以发现其中有很多无效状态。用记忆化搜索可以避免计算这些无效状态。如图：

(i, j) 可以转移到 (i - 4, j), (i - 3, j - 1), (i - 2, j - 2), (i - 1, j - 3)。

我们要求的是 (N, N)，我们只需要知道 (N - 4, N), (N - 3, N - 1), (N - 2, N - 2), (N - 1, N - 3) 即可。

但是按照前面的动态规划的实现方式，类似于 (N - 3, N), (N - 2, N) 这些状态也计算了，而这些状态对推出 (N, N) 没有帮助。记忆化搜索可以避免对这些状态的计算。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
class Solution {
public:
    double soupServings(int N) {
        N = ceil(N / 25.0);
        if(N >= 500)
            return 1.0;

        vector<vector<double> > dp(N + 1, vector<double>(N + 1, -1.0));
        dp[0][0] = 0.5;
        for(int j = 1; j <= N; ++j)
            dp[0][j] = 1.0;
        for(int i = 1; i <= N; ++i)
            dp[i][0] = 0.0;
        return solve(N, N, dp);
    }

private:
    int relu(int x)
    {
        return max(0, x);
    }

    double solve(int i, int j, vector<vector<double>>& dp)
    {
        if(dp[i][j] >= 0)
            return dp[i][j];

        double ans = 0.0;
        ans += solve(relu(i - 4), relu(j), dp);
        ans += solve(relu(i - 3), relu(j - 1), dp);
        ans += solve(relu(i - 2), relu(j - 2), dp);
        ans += solve(relu(i - 1), relu(j - 3), dp);

        return dp[i][j] = ans * 0.25;
    }
};
Share