Python实现单纯形法

单纯形法是一种用于线性规划的算法，它可以在多项式时间内求解线性规划问题。下面是一个简单的 Python 实现。

假设我们要最小化目标函数 $z = 2x_1 + 3x_2$，其中 $x_1$ 和 $x_2$ 分别是变量，满足以下约束条件：

$$\begin{aligned} 3x_1 + 2x_2 &\geq 6 \\ 2x_1 + x_2 &\geq 4 \\ x_1, x_2 &\geq 0 \end{aligned}$$

我们可以将上述约束条件转化为标准形式：

$$\begin{aligned} -3x_1 - 2x_2 + s_1 &= -6 \\ -2x_1 - x_2 + s_2 &= -4 \\ x_1, x_2, s_1, s_2 &\geq 0 \end{aligned}$$

其中 $s_1$ 和 $s_2$ 分别是松弛变量。现在我们可以写出单纯形法的 Python 代码：

import numpy as np

def simplex(c, A, b):
    n, m = A.shape
    c = np.concatenate([c, np.zeros(m)])
    A = np.hstack([A, np.eye(m)])
    b = b.reshape(-1, 1)

    basis = np.arange(n+m, n+m-m-1, -1)
    while True:
        # compute reduced cost coefficients
        c_B = c[basis]
        c_N = c[:-m]
        c_r = c_N - A.T @ c_B

        # optimal solution found
        if np.all(c_r >= 0):
            x_star = np.zeros(n+m)
            x_star[basis] = np.linalg.inv(A[:, basis]) @ b
            return c @ x_star

        # choose entering variable
        j = np.argmin(c_r)

        # compute direction
        d = np.linalg.inv(A[:, basis]) @ A[:, j]

        # unbounded problem
        if np.all(d <= 0):
            return np.inf

        # choose leaving variable
        ratios = b / A[:, j]
        ratios[ratios < 0] = np.inf
        i = np.argmin(ratios)

        # update basis
        basis[i] = j

        # update tableau
        pivot = A[i, j]
        A[i, :] /= pivot
        b[i] /= pivot
        for k in range(n+m):
            if k != i:
                A[k, :] -= A[k, j] * A[i, :]
                b[k] -= A[k, j] * b[i]

现在我们可以使用上述代码求解线性规划问题：

c = np.array([2, 3])
A = np.array([
    [3, 2],
    [2, 1]
])
b = np.array([6, 4])
z_star = simplex(c, A, b)
print(z_star)  # output: 6.0

因此，原线性规划问题的最优解为 $z = 6$，此时 $x_1 = 2$，$x_2 = 0$。