纳什均衡求解器

1 实验目的

• 掌握求解纳什均衡的相关算法
• 锻炼数学基础能力及编程求解问题的能力

2 实验内容

本次实验要求使用Python语言在给定代码框架下编程求解纳什均衡 (Nash Equilibrium, NE), 包括纯策略 NE 与混合策略 NE, 并提交相应源码、输出文件以及实验报告。

3 实现过程

本次大作业要求我们用Python语言写一个纳什均衡求解器，除 Python 标准库外, 仅允许额外引入Numpy和Scipy两个包。Numpy是一个用于科学计算的基础包，使用它可以提高我们代码运行的效率；Scipy是一个用于执行数学、科学和工程计算的python开源代码，我们会用到它的optimize模块求解线性规划问题，用到它的spatial模块求解空间几何问题。下面将简要介绍纳什均衡求解器的实现过程，并主要聚焦于文件读取、纯策略均衡求解和混合策略均衡求解。

3.1 文件读取

每个博弈数据存储在以.nfg为后缀的文件中，我们需要读取这个文件并解析它。经观察发现，文件的最后一行总是payoff列表，倒数第三行总是包含了每个玩家可选的策略数目，根据这两行内容，我们就可以获取关于这个博弈的所有信息。

根据Gambit format的规定，玩家可选的策略数目总是包含在一对花括号里，并用空格分隔，比如{ 3 2 5 }表示玩家1有3种策略，玩家2有2种策略，玩家3有5种策略。知道了它的表示形式后，我们可以使用正则表达式把每个玩家可选的策略数目提取出来，同时我们也相当于知道了玩家的数量。

根据Gambit format的规定，如果总共有$N$个玩家，那么payoff列表中的前$N$个数表示这$N$个玩家都选第一种策略时，每个玩家的收益分别是多少。第一个数对应于玩家1的收益，第二个数对应于玩家2的收益，以此类推。第$N+1$个数到第$2N$个数表示玩家1选第二种策略，其余玩家选第一种策略时，每个玩家的收益。根据约定，首先是第一个玩家的策略递增，当它递增完自己所有可选择的策略之后，它的策略会“回滚”到1，然后第二个玩家的策略依次递增。在我的实现中，我选择了将每个玩家的收益从总的payoff列表中分离出来。分离的方法也很简单，对于payoff列表中第$i$个数（从0开始编号），它对应于玩家$i\%N$（玩家从0开始编号）的收益。现在分离出来的每个玩家的收益仍然是列表的形式，我们需要将其表示为矩阵的形式。

对于有$N$个玩家的博弈，每个玩家的收益都应该是一个$N$维的矩阵。按照前面说过的收益列表与策略的对应关系，每个玩家收益列表中的第一个数，对应于所有玩家都选第一种策略的情况，在收益矩阵中的坐标应该是(1, 1, $\cdots$, 1)。第二个数对应于玩家1选第二种策略，其余玩家选第一种策略的情况，在收益矩阵中的坐标应该是(2, 1, $\cdots$, 1)。如果玩家1有$s_1$种策略，那么收益列表中的第$s_1+1$个数对应于玩家2选策略2，其余玩家选策略1的情况，在收益矩阵中的坐标应该是(1, 2, $\cdots$, 1)。收益列表中第$s_1+2$个数对应于玩家1选策略2，玩家2选策略2，其余玩家选策略1的情况，在收益矩阵中的坐标应该是(2, 2, $\cdots$, 1)，第$s_1+3$个数对应于收益矩阵中的坐标应该是(3, 2, $\cdots$, 1)。可以看出，如果对收益矩阵按照高维优先访问，那么就可以得到收益列表。同理，从收益列表转换到收益矩阵也需要按照高维优先的原则。在Numpy的reshape函数中，提供了一个名为order的参数，当指定其为C时，它会按照C语言的方式，即低维优先进行转换，当指定其为F时，它会按照Fortran的方式，即高维优先进行转换。因此，我们只需要指定order为F就可以很方便的将收益列表转换为收益矩阵。请注意，任意一个收益矩阵中的第$k$维对应于玩家$k$的策略选择。

至此，我们已经从.nfg文件中读取并解析出了每个玩家可选择的策略数，以及每个玩家的收益矩阵。

3.2 纯策略求解

当有三个及以上的玩家时，我们仅需要求解它所有的纯策略纳什均衡点。在上一小节中我们已经介绍过，每个玩家都有自己的收益矩阵，并且收益矩阵中的第$i$维对应于第$i$个玩家的策略。假如对于玩家$i$，我们想要计算在其余玩家选定策略后，它选择哪个策略可以达到最大收益，只需计算玩家$i$的收益矩阵中第$i$维的最大值即可。Numpy自带了求最大值位置的函数argmax，要求在第$i$维上的最大值，只需要指定参数axis=$i$即可，但是此函数却并不适用于本次任务。我们的目标是求所有最大值的位置，当有多个最大值时，要全部计算出来，而argmax只会返回第一个最大值的位置，这显然不符合我们的要求。经过查阅资料，我最终采用了将max函数和where函数结合起来的方案。首先用max函数求第$i$维上的最大值，然后再用where函数找出所有与最大值相同的元素的位置。

当对所有的玩家$i\in [N]$，其中$N$表示玩家数量，都找出了它在第$i$维上的最大值的坐标集合后，将所有玩家的坐标集合求交集即为纯策略纳什均衡点。

3.3 混合策略求解

当只有两个人进行博弈时，我们需要求解它所有的纯策略和混合策略纳什均衡点（无穷多个时求极点）。在我的求解中，参考了Algorithmic Game Theory中3.3节的Labeled polytopes算法。下面将先简要介绍此算法，然后强调几个实现过程中的重点。

3.3.1 Labeled polytopes算法理论

在求解混合策略纳什均衡时，我们需要知道均衡策略中的哪些策略的概率不为0，而最优反应多面体（best response polyhedron）就是一种非常直接的求解方法。记$\mathbf{x}$和$\mathbf{y}$分别表示玩家1和玩家2的混合策略，$\mathbf{A}$和$\mathbf{B}$分别表示玩家1和玩家2的收益矩阵，$u$和$v$分别表示玩家1和玩家2期望收益的上包络（upper envelope），$M$和$N$分别表示玩家1和玩家2可选择的策略数量。玩家1和玩家2的最优反应多面体$\bar{P}$和$\bar{Q}$分别被定义为
P ˉ = { ( x , v ) ∈ R M × R ∣ x ≥ 0 ,   1 ⊤ x = 1 ,   B ⊤ x ≤ 1 v } Q ˉ = { ( y , u ) ∈ R N × R ∣ A y ≤ 1 u ,   y ≥ 0 ,   1 ⊤ y = 1 }

​Pˉ={(x,v)∈RM×R∣x≥0, 1⊤x=1, B⊤x≤1v}Qˉ​={(y,u)∈RN×R∣Ay≤1u, y≥0, 1⊤y=1}​
可以看出，$\bar{P}$中是一些形如$(x_1,x_2,\cdots ,x_M,v)$的元组，$\bar{Q}$中是一些形如$(y_{M+1},y_{M+2},\cdots ,y_{M+N},u)$的元组，其中$x_i$和$y_j$分别表示玩家1选择策略$i$和玩家2选择策略$j$的概率。下面根据作业指导手册上的例子，对$\bar{P}$和$\bar{Q}$中的不等式约束做进一步的说明。玩家1和玩家2的收益矩阵分别如表1和表2所示。

玩家1有3种策略，即$\mathbf{x}=(x_1,x_2,x_3)$，玩家2有两种策略，即$\mathbf{y}=(y_4,y_5)$。对于$\bar{P}$来说，$\mathbf{x}\ge 0$可以产生三个不等式，${\mathbf{B}}^{\top }\mathbf{x}\le 1v$可以产生两个不等式，因此$\bar{P}$的所有不等式约束为
x 1 ≥ 0 x 2 ≥ 0 x 3 ≥ 0 1 ⋅ x 1 + 2 ⋅ x 2 + 2 ⋅ x 3 ≤ v 1 ⋅ x 1 + 3 ⋅ x 2 + 0 ⋅ x 3 ≤ v

x1​x2​x3​1⋅x1​+2⋅x2​+2⋅x3​1⋅x1​+3⋅x2​+0⋅x3​​≥0≥0≥0≤v≤v​

对于$\bar{Q}$来说，$\mathbf{A}\mathbf{y}\le 1u$可以产生三个不等式，$\mathbf{y}\ge 0$可以产生两个不等式，因此$\bar{Q}$的所有不等式约束为
1 ⋅ y 4 + 1 ⋅ y 5 ≤ u 0 ⋅ y 4 + 0 ⋅ y 5 ≤ u 0 ⋅ y 4 + 2 ⋅ y 5 ≤ u y 4 ≥ 0 y 5 ≥ 0

1⋅y4​+1⋅y5​0⋅y4​+0⋅y5​0⋅y4​+2⋅y5​y4​y5​​≤u≤u≤u≥0≥0​

下面引入标记的概念。总的来说，如果某个点被标记为$l$，那么这个点能使第$l$个不等式中的等号成立。因为$\bar{P}$和$\bar{Q}$中不等式的顺序不同，具体来说，在$\bar{P}$中，是先$\mathbf{x}\ge 0$，后${\mathbf{B}}^{\top }\mathbf{x}\le 1v$，而在$\bar{Q}$中，是先$\mathbf{A}\mathbf{y}\le 1u$，后$\mathbf{y}\ge 0$，所以标记的含义也有所区别。令 M = { 1 , 2 , 3 } M = \{1, 2, 3\} M={1,2,3}表示玩家1的标记集合， N = { 4 , 5 } N = \{4, 5\} N={4,5}表示玩家2的标记集合。对于某个点$(\mathbf{y},u)\in \bar{Q}$，它的标记$k\in M\cup N$，如果$k=i\in M$，那么说明第$i$个不等式等号成立，即${\sum }_{j\in N}{a}_{ij}{y}_j=u$；如果$k=j\in N$，那么说明第$j$个不等式等号成立，即$y_j=0$。对于某个点$(\mathbf{x},v)\in \bar{P}$，它的标记$k\in M\cup N$，如果$k=i\in M$，那么$x_i=0$，如果$k=j\in N$，那么${\sum }_{i\in M}{b}_{ij}{x}_i=v$。

对于$((\mathbf{x},v),(\mathbf{y},u))\in \bar{P}\times \bar{Q}$而言，当它是被完全标记时，$(\mathbf{x},\mathbf{y})$是纳什均衡点。完全标记意味着对于任意一个标记$k\in M\cup N$，$k$要么出现在$(\mathbf{x},v)$中，要么出现在$(\mathbf{y},u)$中（可以同时出现）。

从上面的算法描述中可以看出，计算$\bar{P}$和$\bar{Q}$需要知道$u$和$v$，而当收益矩阵$\mathbf{A}$和${\mathbf{B}}^{\top }$非负时，我们可以简化掉$u$和$v$。得到的新的最优反应多面体为
P = { x ∈ R M ∣ x ≥ 0 , B ⊤ x ≤ 1 } Q = { y ∈ R N ∣ A y ≤ 1 , y ≥ 0 }

​P={x∈RM∣x≥0,B⊤x≤1}Q={y∈RN∣Ay≤1,y≥0}​

和原先的相比，去掉了等值约束，并将上包络正则化为1，此时集合中的元素为相应玩家的混合策略（未正则化）。在使用$P$和$Q$时，我们需要去掉其中的0向量，对于非零的向量$\mathbf{x}\in P$或$\mathbf{y}\in Q$，我们可以通过除以$1^{\top }\mathbf{x}$或$1^{\top }\mathbf{y}$将其正则化为概率向量。

上面的$P$和$Q$要求收益矩阵$\mathbf{A}$和${\mathbf{B}}^{\top }$非负，而当不满足这个条件时，可以通过平移收益矩阵，即给矩阵中的每一个元素加上一个大于0的常数，使其满足非负性的条件，并且这样做不会改变纳什均衡的求解结果。因此，对于任意的收益矩阵，都可以用于计算简化后的最优反应多面体。

至此，Labeled polytopes的主要算法就介绍完了。当计算出了$P$和$Q$之后，我们可以枚举这两个集合中的元素（$\mathbf{x}=0$或$\mathbf{y}=0$的情况除外）来找到均衡点。对于 x ∈ P − { 0 } \boldsymbol{x} \in P - \{\boldsymbol{0}\} x∈P−{0}和 y ∈ Q − { 0 } \boldsymbol{y} \in Q - \{\boldsymbol{0}\} y∈Q−{0}，如果$(\mathbf{x},\mathbf{y})$是被完全标记的，那么$\left(\mathbf{x}/(1^{\top }\mathbf{x}),\mathbf{y}/(1^{\top }\mathbf{y})\right)$即为混合策略纳什均衡。

3.3.2 Labeled polytopes算法实现

对于一个二人博弈，我们编写代码的重点是如何求解$P$和$Q$，以及这两个集合里每个元素的标记。显然，对于$P$或$Q$的求解是一个线性规划问题。为了能方便的使用现有的库函数，我们对$P$和$Q$的定义做一个等价转化，即都转换成小于等于的形似：
P = { x ∈ R M ∣ − x ≤ 0 , B ⊤ x ≤ 1 } Q = { y ∈ R N ∣ A y ≤ 1 , − y ≤ 0 }

​P={x∈RM∣−x≤0,B⊤x≤1}Q={y∈RN∣Ay≤1,−y≤0}​

仍然以表1和表2表示的博弈为例，此时$P$中的不等式约束变为
− x 1 ≤ 0 − x 2 ≤ 0 − x 3 ≤ 0 1 ⋅ x 1 + 2 ⋅ x 2 + 2 ⋅ x 3 ≤ 1 1 ⋅ x 1 + 3 ⋅ x 2 + 0 ⋅ x 3 ≤ 1 (1)

\tag{1} −x1​−x2​−x3​1⋅x1​+2⋅x2​+2⋅x3​1⋅x1​+3⋅x2​+0⋅x3​​≤0≤0≤0≤1≤1​(1)
$Q$中的不等式约束变为
1 ⋅ y 4 + 1 ⋅ y 5 ≤ 1 0 ⋅ y 4 + 0 ⋅ y 5 ≤ 1 0 ⋅ y 4 + 2 ⋅ y 5 ≤ 1 − y 4 ≤ 0 − y 5 ≤ 0 (2)

1⋅y4+1⋅y5≤10⋅y4+0⋅y5≤10⋅y4+2⋅y5≤1−y4≤0−y5≤0" role="presentation">1⋅y4+1⋅y50⋅y4+0⋅y50⋅y4+2⋅y5−y4−y5≤1≤1≤1≤0≤0$$\begin{array}{rl}1\cdot y_4+1\cdot y_5& \le 1\\ 0\cdot y_4+0\cdot y_5& \le 1\\ 0\cdot y_4+2\cdot y_5& \le 1\\ -y_4& \le 0\\ -y_5& \le 0\end{array}$$

\tag{2} 1⋅y4​+1⋅y5​0⋅y4​+0⋅y5​0⋅y4​+2⋅y5​−y4​−y5​​≤1≤1≤1≤0≤0​(2)
为了表述方便，下面以$P$为例来介绍求解过程。$P$里面的不等式约束构成了一个空间，现在我们要计算此空间的极点。在我的实现中，极点的计算是由scipy.spatial模块下的HalfspaceIntersection函数完成的。这个函数有两个必要的参数，第一个是halfspaces，它是一个二维的矩阵，形状为(nineq, ndim+1)，它是$\mathbf{M}\mathbf{x}+\mathbf{b}\le 0$这种不等式以$[\mathbf{M};\mathbf{b}]$这种形式表示的系数，其中nineq表示不等式的数目，它等于$\mathbf{M}$的第一维的大小，ndim表示$\mathbf{x}$的长度，它等于$\mathbf{M}$的第二维的大小。比如对于$P$中的不等式约束(1)，它的halfspaces为
[ − 1 0 0 0 0 − 1 0 0 0 0 − 1 0 1 2 2 − 1 1 3 0 − 1 ]

[−10000−10000−10122−1130−1]" role="presentation">⎡⎣⎢⎢⎢⎢⎢⎢−100110−102300−120000−1−1⎤⎦⎥⎥⎥⎥⎥⎥$$\left[\begin{array}{cccc}-1& 0& 0& 0\\ 0& -1& 0& 0\\ 0& 0& -1& 0\\ 1& 2& 2& -1\\ 1& 3& 0& -1\end{array}\right]$$

​−10011​0−1023​00−120​000−1−1​ ​

在由收益矩阵生成halfspaces时需要注意，$P$和$Q$的生成方式不同。$P$是混合策略（即$\mathbf{x}$）小于等于0的约束在前，而$Q$是混合策略（即$\mathbf{y}$)小于等于0的约束在后。因此，不能使用完全相同的代码来生成$P$和$Q$的halfspaces，而是要根据情况来生成。

另一个参数是interior_point，它表示已经确定在halfspaces所定义的区域内部的点，也被称为可行点，可以通过线性规划求得。在HalfspaceIntersection的官方文档中，提供了一种求interior_point的{方法}。考虑具有$\mathbf{M}\mathbf{x}+\mathbf{b}\le 0$形式的半空间，解决线性规划问题
$$\underset{y}{\max }s.t.\mathbf{M}\mathbf{x}+y||{\mathbf{M}}_i||\le -\mathbf{b}$$

其中${\mathbf{M}}_i$表示$\mathbf{M}$的第$i$行，即每一行的范式。通过这种方法，可以产生一个凸多面体最内部的点$\mathbf{x}$。

上面的线性规划求极值问题由scipy.optimize模块的linprog函数完成。这个函数用于在等式和不等式约束下求一个线性函数的最小值。形式化的来说，对于如下形式的问题，
min ⁡ x   c T x s . t . A u b x ≤ b u b A e q x = b e q l ≤ x ≤ u

xmin​s.t.​cTxAub​x≤bub​Aeq​x=beq​l≤x≤u​

其中$x$是决策变量，$c,b_{ub},b_{eq},l$和$u$都是向量，$A_{ub}$和$A_{eq}$都是矩阵，它可以求$c^{\top }x$的最小值。如果想要求一个线性函数的最大值，只需要将其加个负号就可以转化成极小值的问题。在HalfspaceIntersection的官方文档中，提供了求解interior_point的代码，只需要稍作修改就可以应用于纳什均衡求解器。

求得了约束空间的极点之后，接下来需要计算每一个非0极点的标记。对于halfspaces为$[\mathbf{M},\mathbf{b}]$的非0极点$\mathbf{z}$，$\mathbf{M}\mathbf{z}-\mathbf{b}$会得到一个向量，这个向量中值为0的元素所在的位置就是$\mathbf{z}$的标记。在实际编写代码时发现，由于浮点数计算精度的问题，有时候理论上为0的元素实际上并不为0，而是一个非常小的数，比如$e^{-17}$，因此，在我的代码中使用了np.isclose($\mathbf{M}\mathbf{z}$, $-\mathbf{b}$)取代$\mathbf{M}\mathbf{z}$ == $-\mathbf{b}$。

现在我们已经计算出了$P$和$Q$中的极点及其标记，接下来就是枚举这两个集合，找到完全标记的组合，此部分实现较为简单，不再详述。

4 实验结果

以examples文件夹下的.nfg文件为样例输入，.ne文件为样例输出，对代码进行了测试。测试结果显示，程序的输出结果和样例输出完全一致，这说明了代码的正确性。

为了测试代码的运行效率，以input文件夹下的.nfg文件为输入，重复运行程序10次，取运行时间的平均值来作为度量。对于32个纳什均衡求解问题，最终测得的平均运行时间为25.17秒。考虑到写文件操作可能会比较耗时，无法真正的展现求解纳什均衡的运行时间，因此上述结果不包含写文件的耗时。

5 源代码

import re
import numpy as np
from scipy.optimize import linprog
from scipy.spatial import HalfspaceIntersection
import os


def write_rst(equilibria, out_path):
    """
    将求得的纳什均衡点写入文件中
    :param equilibria: 所有纳什均衡点
    :param out_path: 要写入的文件路径
    :return:
    """
    with open(out_path, 'wt') as f:
        for eq in equilibria:
            write_content = ''
            for idx, value in enumerate(eq):
                if type(value) is int:
                    value = str(abs(value))
                else:
                    value = str(abs(value))

                if write_content == '':
                    write_content = value
                else:
                    write_content += ',' + value
            f.write(write_content + '\n')


def read_payoff(in_path):
    """
    读取.nfg文件，并将列表形式的收益解析成矩阵形式
    :param in_path: .nfg文件的路径
    :return: 每个玩家矩阵形式的收益，每个玩家可选择的策略数量
    """
    with open(in_path) as fin:
        all_lines = fin.readlines()
        payoff = all_lines[-1].strip().split()
        player_info = all_lines[-3]
    res = re.search('{ (\d+ )*}', player_info)
    strategy_space = [int(num) for num in res.group()[1: -1].split()]
    player_num = len(strategy_space)
    all_payoff = []
    for player in range(player_num):
        payoff_matrix = [int(pay) for idx, pay in enumerate(payoff) if idx % player_num == player]
        payoff_matrix = np.array(payoff_matrix).reshape(strategy_space, order='F')
        all_payoff.append(payoff_matrix)
    return all_payoff, strategy_space


def create_half_spaces(payoff_matrix, is_p=False):
    """
    根据收益矩阵创建half_spaces
    :param payoff_matrix: 收益矩阵，正是根据它来创建half_spaces
    :param is_p: 根据参考书上的算法，polyhedron P和Q中half_spaces的顺序不一样，此变量用于确定是哪一种情况
    :return: 创建好的half_spaces，加入payoff_matrix的维度维d1 * d2，那么half_spaces的维度就是(d1+d2) * (d2 + 1)
    """
    d1, d2 = payoff_matrix.shape  # 收益矩阵的第一维和第二维的大小
    if is_p:
        b = np.append(np.zeros(d2), -np.ones(d1))
        payoff_matrix = np.append(-np.eye(d2), payoff_matrix, axis=0)
    else:
        b = np.append(-np.ones(d1), np.zeros(d2))
        payoff_matrix = np.append(payoff_matrix, -np.eye(d2), axis=0)

    half_spaces = np.column_stack((payoff_matrix, b.transpose()))
    return half_spaces


def find_feasible_point(half_spaces):
    """
    求线性规划空间内部的点，参考了HalfspaceIntersection官方文档，
    地址：https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.HalfspaceIntersection.html?highlight=halfspaceintersection#scipy.spatial.HalfspaceIntersection
    :param half_spaces:
    :return: half_spaces所表示的线性空间的内部的点
    """
    norm_vector = np.reshape(
        np.linalg.norm(half_spaces[:, :-1], axis=1), (half_spaces.shape[0], 1)
    )
    c = np.zeros((half_spaces.shape[1],))
    c[-1] = -1
    A_ub = np.hstack((half_spaces[:, :-1], norm_vector))
    b_ub = -half_spaces[:, -1:]
    res = linprog(c, A_ub=A_ub, b_ub=b_ub)
    return res.x[:-1]


def compute_vertices_and_labels(half_spaces):
    """
    根据half_spaces来计算vertices，需要注意结果中不能含有0向量
    :param half_spaces:
    :return:
    """
    feasible_point = find_feasible_point(half_spaces)
    hs = HalfspaceIntersection(half_spaces, feasible_point)
    hs.close()
    rst = []
    for v in hs.intersections:
        if not np.all(np.isclose(v, 0)):  # 点v不能是(0, 0)，也不能是无穷大（小）
            # 根据Av - b = 0来寻找点v的标记
            b = half_spaces[:, -1]
            A = half_spaces[:, :-1]
            labeled_v = np.where(np.isclose(np.dot(A, v), -b))[0]
            rst.append((v, labeled_v))
    return rst


def mixed_nash_equilibrium(all_payoff):
    """
    计算两个玩家的混合策略纳什均衡（含纯策略）
    :param all_payoff: 两个玩家的收益矩阵
    :return: 求得的所有纯策略纳什均衡点
    """
    payoff_mat1, payoff_mat2 = all_payoff
    # 对收益矩阵平移不会影响纳什均衡点，当收益矩阵非负时，可以简化掉u和v
    if np.min(payoff_mat1) < 0:
        payoff_mat1 = payoff_mat1 + abs(np.min(payoff_mat1))
    if np.min(payoff_mat2) < 0:
        payoff_mat2 = payoff_mat2 + abs(np.min(payoff_mat2))

    p1_strategy_num, p2_strategy_num = payoff_mat1.shape  # 玩家1的策略数量和玩家2的策略数量
    all_strategy_num = p1_strategy_num + p2_strategy_num

    p1_half_spaces = create_half_spaces(payoff_mat2.transpose(), is_p=True)
    p2_half_spaces = create_half_spaces(payoff_mat1, is_p=False)
    rst = []
    p_vertices = compute_vertices_and_labels(p1_half_spaces)  # 玩家1的所有顶点，这里用p是想与参考书统一起来，因为参考书用P表示的polyhedron，下同
    q_vertices = compute_vertices_and_labels(p2_half_spaces)
    for vec_x, x_l in p_vertices:  # 点x（以向量的形式表示）和x的标记，分别表示交点的坐标，以及是哪几个不等式相交
        all_labels = np.zeros(all_strategy_num, dtype=int)
        all_labels[x_l] = 1
        for vec_y, y_l in q_vertices:
            all_labels[y_l] += 1
            if np.all(all_labels):
                rst.append(np.append(vec_x / sum(vec_x), vec_y / sum(vec_y)))  # 做归一化
            all_labels[list(y_l)] -= 1
    return rst


def compare_eq(computed_NE, expected_NEs):
    '''
    用于测试求解结果是否正确
    :param computed_NE: 通过程序计算出来的均衡点
    :param expected_NEs: 正确的纳什均衡点
    :return: 如果计算出来的NE是正确的，则返回True，否则返回False
    '''
    for one_NE in expected_NEs:
        diff = np.abs(computed_NE - one_NE)
        if np.sum(diff) <= 1e-4:
            return True
    return False


def pure_nash_equilibrium(all_payoff):
    """
    计算三个及以上的玩家的混合策略纳什均衡
    :param all_payoff: 是一个list，第i个元素是第i个玩家的收益矩阵
    :return: 求得的所有混合策略纳什均衡点
    """
    player_num = len(all_payoff)
    strategy_space = all_payoff[0].shape
    max_payoff_coordinate = set()  # 玩家的最大收益在收益矩阵中的坐标
    for player, payoff_matrix in enumerate(all_payoff):
        # 第i个玩家在第i维上选取动作
        temp_max_payoff_coordinate = set()
        max_pos = np.where(payoff_matrix == np.max(payoff_matrix, axis=player, keepdims=True))
        max_pos = list(max_pos)
        # 因为np.where返回的是tuple，所以需要对其进一步解析才能得到坐标
        for i in range(len(max_pos[0])):
            temp_list = []
            for j in range(player_num):
                temp_list.append(max_pos[j][i])
            temp_max_payoff_coordinate.add(tuple(temp_list))
        if player == 0:
            max_payoff_coordinate = temp_max_payoff_coordinate
        else:
            max_payoff_coordinate = max_payoff_coordinate.intersection(temp_max_payoff_coordinate)

    # 将求得的纯策略纳什均衡转换为指定的输出格式
    equilibria = []
    for equilibrium_coord in max_payoff_coordinate:
        equilibrium = [0] * (sum(strategy_space))
        pre_len = 0
        for pos, strategy_num in zip(equilibrium_coord, strategy_space):
            equilibrium[pre_len + pos] = 1
            pre_len += strategy_num
        equilibria.append(equilibrium)
    return equilibria


def nash(in_path, out_path):
    """
    求解纳什均衡，对于两个玩家的博弈，需要求解所有的纯策略均衡点和混合策略均衡点，对于三个及以上玩家的博弈，只需求解纯策略均衡点
    :param in_path: .nfg文件的路径，此文件包含了与博弈有关的所有信息
    :param out_path: 均衡点输出文件的路径
    :return:
    """
    all_payoff, strategy_space = read_payoff(in_path)
    if len(strategy_space) == 2:
        # 二人博弈，求解混合策略纳什均衡
        equilibria = mixed_nash_equilibrium(all_payoff)
    else:
        # 三人以上博弈，求解纯策略纳什均衡
        equilibria = pure_nash_equilibrium(all_payoff)

    write_rst(equilibria, out_path)


def main():
    for f in os.listdir('input'):
        if f.endswith('.nfg') and not f.startswith('._'):
            nash('input/' + f, 'output/' + f.replace('nfg', 'ne'))


if __name__ == '__main__':
    main()
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
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217

实验中用到的测试文件：Github