单纯形法是线性规划中最经典且广泛应用的求解方法,通过在可行解的边界上移动,逐步逼近最优解。它从一个初始基本可行解开始,不断优化目标函数值,直到找到最优解。对偶单纯形法则是单纯形法的一种变形,尤其适用于特定类型的线性规划问题。不同于标准的单纯形法,对偶单纯形法从一个对偶可行但原始不可行的初始解出发,通过逐步改善解的可行性和最优性,最终找到最优解。对偶单纯形法在快速调整和重新优化方面表现出色,特别是在处理新的约束或变量被添加到现有模型中的情况时,这使得它在实际应用中具有独特的优势。
单纯形法比较 | 对偶单纯形流程图 |
---|---|
一、对偶单纯形算法
对偶单纯形法是一种用于求解线性规划问题的数学方法,由美国数学家C.莱姆基在1954年提出。这种方法基于对偶理论,通过迭代过程逐步搜索原始问题的最优解。与单纯形法不同,对偶单纯形法从满足对偶可行性条件出发,通过迭代逐步逼近原始问题的最优解。在迭代过程中,始终保持基解的对偶可行性,使不可行性逐步消失。
对偶单纯形法的基本思想是:从原规划的一个对偶可行解出发;然后检验原规划的基本解是否可行,即是否有负的分量,如果有小于零的分量,则进行迭代,求另一个基本解,此基本解对应着另一个对偶可行解(检验数非正)。如果得到的基本解的分量皆非负,则该基本解为最优解。也就是说,对偶单纯形法在迭代过程中始终保持对偶解的可行性(即检验数非正),使原规划的基本解由不可行逐步变为可行,当同时得到对偶规划与原规划的可行解时,便得到原规划的最优解。
1.1 对偶问题算法步骤
a)根据线性规划典式形式,建立初始对偶单纯形表。此表对应原规划的一个基本解。要求检验数行各元素一定非正,原规划的基本解允许有小于零的分量;
\[\theta = \min \left\{ \frac{\sigma_j}{a_{ij}} \middle| a_{ij} < 0 \right\} = \frac{\sigma_k}{a_{ik}} \]
b)若基本解的所有分量皆非负,则得到原规划的最优解,停止计算;若基本解中有小于零的分量\(b_j < 0\),并且\(b_j\)所在行各系数\(a_{ij} \geq 0\),则原规划没有可行解,停止计算;若\(b_j < 0\),并且存在\(a_{ij} < 0\),则确定\(x_i\)为出基变量,并计算确定\(x_k\)为进基变量。若有多个\(b_j < 0\),则选择最小的进行分析计算;
c)以\(a_{ik}\)为中心元素,按照与单纯形类似的方法,在表中进行迭代计算,返回步骤 b)。
1.2 练习
例1:用对偶单纯形法求解下面线性规划
\[\begin{aligned} \text{min}\quad w= 2x_1 + 3x_2 + 4x_3 \\ \text{s.t. } \begin{cases} x_1 + 2x_2 + x_3 \geq 3 \\ 2x_1 - x_2 + 3x_3 \geq 4 \\ x_1, x_2, x_3 \geq 0 \end{cases} \end{aligned} \]\[\begin{aligned} \text{max } \quad w= -2x_1 - 3x_2 - 4x_3 \\ \text{s.t. } \begin{cases} -x_1 - 2x_2 - x_3 + x_4 = -3 \\ -2x_1 + x_2 - 3x_3 + x_5 = -4 \\ x_1, \ldots, x_5 \geq 0 \end{cases} \end{aligned} \]得原问题的最优解为$x_1 = 11/5,x_2 = 2/5,x_3 = 3$;最优值为$w = 28/5$。
二、例题
利用对偶单纯形法求解线性规划问题:
\[\begin{cases} \min f(x) = 4x_{1} + 12x_{2} + 18x_{3} \\ \text{s.t. } x_{1} + 3x_{3} \geq 3 \\ \quad \quad 2x_{2} + 2x_{3} \geq 5 \\ \quad \quad x_{1}, x_{2}, x_{3} \geq 0 \end{cases} \]解:先标准化:
\[\begin{cases} \max f(x) =- 4x_{1} - 12x_{2} - 18x_{3} \\ \text{s.t.}\quad x_{1} + 3x_{3} - x_{4} = 3 \\ \quad \quad 2x_{2} + 2x_{3} - x_{5} = 5 \\ \quad \quad x_{1}, x_{2}, x_{3}, x_{4}, x_{5} \geq 0 \end{cases} \]为了将 $$ B = (p_{4}, p_{5}) $$ 作为初始基,对约束条件两边同乘 \(-1\):
\[\begin{cases} -x_{1} - 3x_{3} + x_{4} = -3 \\ -2x_{2} - 2x_{3} + x_{5} = -5 \end{cases} \]\[C^{T}-c_{B}^{T}BA = (-4, -12, -18, 0, 0)^{T} \]对应的初始单纯形表如下:
变量 | $$ x_{1} $$ | $$ x_{2} $$ | $$ x_{3} $$ | $$ x_{4} $$ | $$ x_{5} $$ | |
---|---|---|---|---|---|---|
$$ f $$ | 0 | -4 | -12 | -18 | 0 | 0 |
$$ x_{4} $$ | -3 | -1 | 0 | -3 | 1 | 0 |
$$ x_{5} $$ | -5 | 0 | -2 | -2 | 0 | 1 |
当前解 $$ x = (0, 0, 0, -3, -5)^{T} $$ 为非可行解,但是它是对偶可行解。
因为:
\[\min\left\{ -3, -5 \right\} = -5, \quad r = 2, \quad x_{j2} = x_{5} \text{ 为离基变量} \]\[\frac{b_{0s}}{b_{rs}} = \min\left\{ \frac{b_{0j}}{b_{rj}} \mid b_{rj} < 0 \right\} = \min\left\{ \frac{-12}{-2}, \frac{-18}{-2} \right\} = 6 = \frac{b_{02}}{b_{22}}, \quad s = 2 \]因此,进基变量 $$ x_{s} = x_{2} $$,且 $$ b_{22} = -2 $$ 为主元。
经旋转变换后,有:
变量 | $$ x_{1} $$ | $$ x_{2} $$ | $$ x_{3} $$ | $$ x_{4} $$ | $$ x_{5} $$ | |
---|---|---|---|---|---|---|
$$ f $$ | 30 | -4 | 0 | -6 | 0 | -6 |
$$ x_{4} $$ | -3 | -1 | 0 | -3 | 1 | 0 |
$$ x_{2} $$ | 5/2 | 0 | 1 | 1 | 0 | -1/2 |
此时 $$ b_{i0} \geq 0 , (i=1,...,m) $$ 不成立,仍为不可行解。
因为:
\[b_{r0} = \min\left\{b_{i0} \mid b_{i0} < 0, i=1,...,m \right\} = \min\left\{ -3\right\} = -3, \quad x_{4} \text{ 为离基变量}, \quad r = 1 \]\[\frac{b_{0s}}{b_{rs}} = \min\left\{ \frac{b_{0j}}{b_{rj}} \mid b_{rj} < 0 \right\} = \min\left\{ \frac{-4}{-1}, \frac{-6}{-3} \right\} = 2 = \frac{b_{03}}{b_{r3}}, \quad s = 3 \]因此,进基变量 $$ x_{s} = x_{3} $$,且 $$ b_{13} = -3 $$ 为主元。
经旋转变换后,有:
变量 | $$ x_{1} $$ | $$ x_{2} $$ | $$ x_{3} $$ | $$ x_{4} $$ | $$ x_{5} $$ | |
---|---|---|---|---|---|---|
$$ f $$ | 36 | -2 | 0 | 0 | -2 | -6 |
$$ x_{3} $$ | 1 | 1/3 | 0 | 1 | -1/3 | 0 |
$$ x_{2} $$ | 3/2 | -1/3 | 1 | 0 | 1/3 | -1/2 |
此时 $$ b_{i0} \geq 0, (i=1,...,m) $$ 成立,是可行解,算法终止。
最优解和最优值:
\[\bar{x} = \left(0, \frac{3}{2}, 1\right)^{T} \quad f(\bar{x}) = 36 \]三、对偶单纯形的Python求解
假设食物的各种营养成分、价格如下表:
Food | Energy (能量) | Protein (蛋白质) | Calcium (钙) | Price |
---|---|---|---|---|
Oatmeal (燕麦) | 110 | 4 | 2 | 3 |
Whole milk (全奶) | 160 | 8 | 285 | 9 |
Cherry pie (樱桃派) | 420 | 4 | 22 | 20 |
Pork with beans (猪肉) | 260 | 14 | 80 | 19 |
要求我们买的食物中,至少要有2000的能量,55的蛋白质,800的钙,怎样买最省钱?
设买燕麦、全奶、樱桃派、猪肉为\(x_1, x_2, x_3, x_4\),于是可以写出如下的不等式组:
\[\begin{aligned} & \min \quad 3x_1 + 9x_2 + 20x_3 + 19x_4 \quad \text{money} \\ & \text{s.t.} \\ & 110x_1 + 160x_2 + 420x_3 + 260x_4 \geq 2000 \quad \text{energy} \\ & 4x_1 + 8x_2 + 4x_3 + 14x_4 \geq 55 \quad \text{protein} \\ & 2x_1 + 285x_2 + 22x_3 + 80x_4 \geq 800 \quad \text{calcium} \\ & x_1, x_2, x_3, x_4 \geq 0 \end{aligned}\]"""
对偶单纯形法
min z=15*x1+24*x2+5*x3
6*x2+x3>=2
5*x1+2*x2+x3>=1
x1,x2>=0
标准化后格式
max z=-15*x1-24*x2-5*x3
-6*x2-x3+x4=-2
-5*x1-2*x2-x3+x5=-1
x1,x2,x3,x4,x5>=0
"""
import numpy as np
# 根据线性规划要改变的数据
Cj = [-3, -9, -20, -19, 0, 0, 0] # Coefficients of the objective function
constraints_matrix = [
[-2000, -110, -160, -420, -260, 1, 0, 0],
[-55, -4, -8, -4, -14, 0, 1, 0],
[-800, -2, -285, -22, -80, 0, 0, 1]
] # Constraints matrix (with RHS)
z = [0, -3, -9, -20, -19, 0, 0, 0] # Z-row
Xb = [5, 6, 7] # Basis variables (indices of slack variables)
Cb = [0, 0, 0] # Coefficients of basis variables
# 为便于表示而生成
C = [0] + Cj
# Output simplex table
def show():
print("--" * 50)
print("| Cj ", end='|')
for i in Cj:
print(f"{i:.2f}".center(10), end='|')
print()
print("--" * 50)
print("| Cb | Xb | b ", end='|')
for i in range(len(Cj)):
print(("X" + str(i + 1)).center(10), end='|')
print()
print("--" * 50)
for i in range(len(constraints_matrix)):
print("|", end="")
print(f"{Cb[i]:.2f}".center(10), end='|')
print(("X" + str(Xb[i])).center(10), end='|')
for j in range(len(constraints_matrix[0])):
print(f"{constraints_matrix[i][j]:.2f}".center(10), end='|')
print()
print("--" * 50)
print("| Z ", end='|')
for i in range(len(z)):
print(f"{z[i]:.2f}".center(10), end='|')
print()
print("--" * 50)
print("**" * 50)
def fu():
global constraints_matrix, z, Xb, Cb # Declare global variables
# Extracting the first column
first_column = [row[0] for row in constraints_matrix]
# Finding the index of the minimum value in the first column
min_index = first_column.index(min(first_column))
# Extract the row corresponding to min_index
min_row = constraints_matrix[min_index]
# Calculate the ratios z / min_row
ratios = [z[i] / min_row[i] if min_row[i] != 0 else float('inf') for i in range(len(z))]
# Filter out the non-positive values from z and corresponding min_row
filtered_ratios = [ratios[i] for i in range(len(ratios)) if z[i] < 0 and min_row[i] < 0]
# Find the minimum ratio and its index in z
if filtered_ratios:
min_ratio = min(filtered_ratios)
min_ratio_index = ratios.index(min_ratio)
else:
print("No valid ratio found")
return
# Updating Xb
Xb[min_index] = min_ratio_index
print("Updated basis variables Xb:")
print(Xb)
# Updating Cb based on Xb
Cb = [C[Xb[i]] for i in range(len(Xb))]
print("Updated coefficients of basis variables Cb:")
print(Cb)
# 根据 Xb 提取列以形成矩阵 B
B = [[row[xb] for xb in Xb] for row in constraints_matrix]
try:
# Calculate the inverse of matrix B
B_matrix = np.array(B)
B_inv = np.linalg.inv(B_matrix)
except np.linalg.LinAlgError:
print("Matrix B is singular, skipping this iteration.")
return
# Update z
z_update = (np.array(C[1:]) - np.dot(np.array(Cb), np.dot(B_inv, np.array(constraints_matrix)[:, 1:]))).tolist()
z = [z[0]] + z_update # Add the first element back to z
# Update constraints_matrix with the inverse of B
constraints_matrix = np.dot(B_inv, np.array(constraints_matrix)).tolist()
# 输出最优解和最优值
def print_optimal_solution():
# 最优解:基变量的值
optimal_solution = [0] * len(Cj) # 初始所有变量的值为0
for i, xb in enumerate(Xb):
optimal_solution[xb - 1] = constraints_matrix[i][0] # 设置基变量的值
# 最优值:计算目标函数值
optimal_value = sum(Cj[i] * optimal_solution[i] for i in range(len(Cj)))
optimal_value = -optimal_value # 取相反数
# 打印最优解和最优值,保留两位小数
print(f"最优解(变量的取值):{[round(val, 2) for val in optimal_solution]}")
print(f"最优值(目标函数值的相反数):{optimal_value:.2f}")
# 迭代执行对偶单纯形法
while any(row[0] < 0 for row in constraints_matrix):
fu()
show()
# 输出最后的最优解和最优值
print_optimal_solution()
#最终的单纯形表
Updated basis variables Xb:
[1, 6, 2]
Updated coefficients of basis variables Cb:
[-3, 0, -9]
----------------------------------------------------------------------------------------------------
| Cj | -3.00 | -9.00 | -20.00 | -19.00 | 0.00 | 0.00 | 0.00 |
----------------------------------------------------------------------------------------------------
| Cb | Xb | b | X1 | X2 | X3 | X4 | X5 | X6 | X7 |
----------------------------------------------------------------------------------------------------
| -3.00 | X1 | 14.24 | 1.00 | 0.00 | 3.74 | 1.98 | -0.01 | 0.00 | 0.01 |
| 0.00 | X6 | 23.63 | 0.00 | 0.00 | 11.38 | -3.96 | -0.04 | 1.00 | -0.01 |
| -9.00 | X2 | 2.71 | 0.00 | 1.00 | 0.05 | 0.27 | 0.00 | 0.00 | -0.00 |
----------------------------------------------------------------------------------------------------
| Z | 0.00 | 0.00 | 0.00 | -8.31 | -10.67 | -0.03 | 0.00 | -0.02 |
----------------------------------------------------------------------------------------------------
最优解(变量的取值):[14.24, 2.71, 0, 0, 0, 23.63, 0]
最优值(目标函数值的相反数):67.10
总结
单纯形法和对偶单纯形法是线性规划中最常用的两种算法,它们在解决优化问题方面发挥着关键作用。单纯形法通过从一个初始的基本可行解出发,沿着多面体的边逐步移动,寻找使目标函数值不断改善的邻近顶点,直到找到最优解。相较之下,对偶单纯形法则采用了不同的思路。它从一个对偶可行解开始,重点在于原问题的对偶问题的可行性,而非原问题的直接可行性。通过逐步调整,使得对偶条件和原始问题的可行性条件同时得到满足,从而找到最优解。对偶单纯形法特别适用于那些初始解不满足原问题约束但满足对偶约束的情况。这在需要快速适应变化的约束条件或在求解过程中频繁调整约束的场景中具有明显优势。
两种方法的有效结合和运用,可以大幅提升线性规划问题的求解效率和灵活性。理解和掌握这两种方法对于优化问题的解决至关重要,它不仅帮助优化计算速度和准确性,还能在面对不同类型的线性规划问题时,选择最适合的算法。无论是在理论研究还是实际应用中,这两种方法都提供了强大的工具,帮助我们应对各种复杂的决策问题。