运筹学练习Python精解——线性规划
练习1
某企业有三个车间生产同一种产品。每件产品由四个零件 1 和三个零件 2 组成。两个零件需耗用两种原材料 A 和 B。已知这两种原材料的供应量分别为 300kg 和 500kg。由于三个车间拥有的设备及工艺条件不同,每个工班原材料耗用量和零件产量也不同。见下表(三个车间每班用料和生产情况),问三个车间就各开多少个班,才能使该产品的配套数达到最大?
| 车间 | A 材料 | B 材料 | 零件1 | 零件2 |
|---|---|---|---|---|
| 一车间 | 8 | 6 | 7 | 5 |
| 二车间 | 5 | 9 | 6 | 9 |
| 三车间 | 3 | 8 | 8 | 4 |
1.1问题描述:
我们有三个车间,每个车间的每班原材料耗用量和零件产量不同。我们的目标是确定每个车间开多少个班,以使得产品的配套数(即最终的产品数)最大。每个产品由4个零件1和3个零件2组成。我们有两个原材料A和B,供应量分别为300kg和500kg。
- 决策变量:
- \(x_1\) 表示一车间的班数
- \(x_2\) 表示二车间的班数
- \(x_3\) 表示三车间的班数
- 目标函数: 我们要最大化产品的配套数。每个产品需要4个零件1和3个零件2。假设最终产品的配套数为P,那么我们可以定义一个辅助变量来表示该目标:\[P = \max \left( \frac{7x_1 + 6x_2 + 8x_3}{4}, \frac{5x_1 + 9x_2 + 4x_3}{3} \right) \]
- 约束条件:
- 原材料A的供应量不能超过300kg: \(8x_1 + 5x_2 + 3x_3 \leq 300\)
- 原材料B的供应量不能超过500kg: \(6x_1 + 9x_2 + 8x_3 \leq 500\)
- 零件1和零件2必须满足产品配套需求: \(4P \leq 7x_1 + 6x_2 + 8x_3 \quad 3P \leq 5x_1 + 9x_2 + 4x_3\)
- 非负约束: \(x_1 \geq 0, \quad x_2 \geq 0, \quad x_3 \geq 0\)
- 求解方法: 使用线性规划求解工具,如Python的SciPy库或其他优化工具,求解这个线性规划问题。
1.2 Python求解
我们使用Python的SciPy库来求解这个线性规划问题。
import pulp
# 创建一个线性规划问题
lp = pulp.LpProblem("Maximize_P", pulp.LpMaximize)
# 定义变量,x1, x2, x3为非负连续变量,P为非负整数变量
x1 = pulp.LpVariable('x1', lowBound=0, cat='Continuous')
x2 = pulp.LpVariable('x2', lowBound=0, cat='Continuous')
x3 = pulp.LpVariable('x3', lowBound=0, cat='Continuous')
P = pulp.LpVariable('P', lowBound=0, cat='Integer')
# 添加目标函数
lp += P
# 添加约束条件
lp += 7*x1 + 6*x2 + 8*x3 - 4*P >= 0
lp += 5*x1 + 9*x2 + 4*x3 - 3*P >= 0
lp += 8*x1 + 5*x2 + 3*x3 <= 300
lp += 6*x1 + 9*x2 + 8*x3 <= 500
# 求解问题
status = lp.solve()
# 输出结果
if status == pulp.LpStatusOptimal:
print("Optimal integer value of P:", pulp.value(P))
print("x1:", pulp.value(x1))
print("x2:", pulp.value(x2))
print("x3:", pulp.value(x3))
else:
print("Failed to find an optimal solution.")
运行上述代码将会得到各车间的最优班次安排和最大产品配套数。
Optimal integer value of P: 116.0
x1: 12.571429
x2: 16.190476
x3: 34.857143
练习2
2.1 问题分析
某染化料厂要用 C、P、H 三种原料混合配制出 A、B、D 三种不同规格的产品。原料 C、P、H 每天的最大供应量分别为 100kg、100kg、60kg,每千克单价分别为 65 元、25 元、35 元。产品 A 要求原料 C 含量不少于 50%,含原料 P 不超过 25%;产品 B 含 C 不得少于 25%,P 不超过 50%;产品 D 的原料配比没有限制。产品 A、B、D 每千克的单价分别为 50 元、35 元、25 元。问应如何安排生产,使得利润为最大?
2.2 数学建模
设:
- \(x_A\) 为生产产品 A 的重量 (kg)
- \(x_B\) 为生产产品 B 的重量 (kg)
- \(x_D\) 为生产产品 D 的重量 (kg)
我们需要最大化利润:
其中 \(C_A\), \(P_A\) \(H_A\) 分别表示生产产品 A 所需的原料 C, P, H 的数量,同理,\(C_B\), \(P_B\), \(H_B\),\(C_D\), \(P_D\), \(H_D\) 分别表示生产产品 B 和 D 所需的原料 C, P, H 的数量。
根据题意,有如下约束条件:
-
原料限制:
\[\begin{aligned} C_A + C_B + C_D &\leq 100 \\ P_A + P_B + P_D &\leq 100 \\ H_A + H_B + H_D &\leq 60 \end{aligned} \] -
产品 A 的原料配比:
\[\begin{aligned} C_A &\geq 0.5x_A \\ P_A &\leq 0.25x_A \end{aligned} \] -
产品 B 的原料配比:
\[\begin{aligned} C_B &\geq 0.25x_B \\ P_B &\leq 0.5x_B \end{aligned} \] -
产品 D 的原料配比没有限制。
2.3 Python 求解程序
我们可以使用 scipy.optimize 库中的 linprog 函数来求解这个线性规划问题。
import numpy as np
from scipy.optimize import linprog
# 目标函数的系数(求最大化,需要取相反数)
c = np.array([-50, -35, -25, 65, 65, 65, 25, 25, 25, 35, 35, 35])
# 定义不等式约束的系数矩阵和常数项向量(包括原材料供应限制和产品成分要求)
A_ub = np.array([
[0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0], # 原料C的消耗约束
[0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0], # 原料P的消耗约束
[0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1], # 原料H的消耗约束
[0.5, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0], # 产品A的C含量要求
[-0.25, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0], # 产品A的C含量要求
[0, 0.25, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0], # 产品B的C含量要求
[0, -0.5, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0] # 产品B的C含量要求
])
b_ub = np.array([100, 100, 60, 0, 0, 0, 0])
# 定义变量的非负约束
bounds = [(0, None) for _ in range(12)]
# 求解线性规划问题
res = linprog(c, A_ub=A_ub, b_ub=b_ub, bounds=bounds, method='highs')
# 输出结果
if res.success:
print(f"最大利润: {-res.fun}")
for i in range(12):
print(f"x[{i+1}]: {res.x[i]}")
else:
print("优化未成功")
print(res.message)
最大利润: 92.14285714285714
产品A的生产量: 1.0
产品B的生产量: 1.0
产品D的生产量: 0.2857142857142857
练习3
3.1 问题分析
纽约市在平时的一天里警察巡逻的需求如下:
| 时段 | 时间 | 巡逻警察数(最低限) |
|---|---|---|
| 1 | 凌晨 12 点 - 凌晨 4 点 | 6 |
| 2 | 凌晨 4 点 - 早上 8 点 | 4 |
| 3 | 早上 8 点 - 正午 12 点 | 14 |
| 4 | 正午 12 点 - 下午 4 点 | 8 |
| 5 | 下午 4 点 - 晚上 8 点 | 12 |
| 6 | 晚上 8 点 - 凌晨 12 点 | 16 |
每个警察的班次为连续工作8小时,我们需要确定这6个班次各安排多少巡警可以使得满足人员需求的条件下使用的警力最少。
3.2 数学建模
设每个时段起始时的巡逻警数量分别为 \(x_1, x_2, x_3, x_4, x_5, x_6\),其中:
- \(x_1\) 表示午夜12点到早上8点的巡逻警数量
- \(x_2\) 表示凌晨4点到中午12点的巡逻警数量
- \(x_3\) 表示早上8点到下午4点的巡逻警数量
- \(x_4\) 表示中午12点到晚上8点的巡逻警数量
- \(x_5\) 表示下午4点到凌晨12点的巡逻警数量
- \(x_6\) 表示晚上8点到凌晨4点的巡逻警数量
根据需求,我们得到以下约束条件:
\(x_1 + x_6 \geq 6\)
\(x_1 + x_2 \geq 4\)
\(x_2 + x_3 \geq 14\)
\(x_3 + x_4 \geq 8\)
\(x_4 + x_5 \geq 12\)
\(x_5 + x_6 \geq 16\)
目标是最小化巡逻警总数,即最小化 \(x_1 + x_2 + x_3 + x_4 + x_5 + x_6\)
3.3 Python 求解程序
我们可以使用 scipy.optimize 库中的 linprog 函数来求解这个线性规划问题。
import numpy as np
from scipy.optimize import linprog
# 系数矩阵
A = np.array([
[-1, 0, 0, 0, 0, -1],
[-1, -1, 0, 0, 0, 0],
[ 0, -1, -1, 0, 0, 0],
[ 0, 0, -1, -1, 0, 0],
[ 0, 0, 0, -1, -1, 0],
[ 0, 0, 0, 0, -1, -1]
])
# 需求矩阵
b = np.array([-6, -4, -14, -8, -12, -16])
# 目标函数系数
c = np.array([1, 1, 1, 1, 1, 1])
# 求解线性规划问题
res = linprog(c, A_ub=A, b_ub=b, bounds=(0, None), method='highs')
# 输出结果
print(f"最小巡逻警数: {res.fun}")
print(f"各时段巡逻警数: {res.x}")
最小巡逻警数: 32.0
各时段巡逻警数: [ 2. 2. 12. 0. 12. 4.]
练习4
某工厂在今后4个月内需租用仓库堆放物资。已知各月份所需仓库面积如下表;仓库租用费用随和同期定,期限越长折扣越大,具体数据见下表。租用仓库的合同每月初都可办理,每份合同将具体规定租用面积数和期限。因此该厂可根据需要在任何一个月初办理租用合同,每次可签一份或多份和合同。工厂需决定如何租用可使总的所付租用费用最小。请建立此问题的线性规划模型。
| 月份 | 1 | 2 | 3 | 4 | 合同租用期限(个月) | 1 | 2 | 3 | 4 |
|---|---|---|---|---|---|---|---|---|---|
| 所需仓库面积(百平米) | 15 | 10 | 20 | 12 | 租用费用(元/百平米) | 2800 | 4500 | 6000 | 7300 |
4.1 建模三要素
-
决策变量
我们定义决策变量 \(x_{i,j}\),表示在第\(i\)个月租用\(j\)个月的仓库面积数(百平米)。 -
目标函数
最小化总租用费用。 -
约束条件
满足各月份的仓库面积需求;租用面积必须非负。
4.2 建立线性规划模型
我们将模型用数学公式表达如下:
- 目标函数:
其中,\(c_j\) 是租用j个月的费用。
- 约束条件:
对于每个月的需求约束:
- 非负约束:
4.3 Python代码
import pulp
# 定义模型
model = pulp.LpProblem("Minimize_Warehouse_Rental_Cost", pulp.LpMinimize)
# 定义决策变量
x = pulp.LpVariable.dicts("x", [(i, j) for i in range(1, 5) for j in range(1, 5)], lowBound=0, cat='Continuous')
# 定义租用费用(元/百平米)
cost = {1: 2800, 2: 4500, 3: 6000, 4: 7300}
# 添加目标函数
model += pulp.lpSum(cost[j] * x[i, j] for i in range(1, 5) for j in range(1, 5))
# 添加约束条件
# 月份需求
demand = {1: 15, 2: 10, 3: 20, 4: 12}
for i in range(1, 5):
model += pulp.lpSum(x[k, j] for j in range(1, 5) for k in range(max(1, i - j + 1), i + 1)) >= demand[i]
# 求解模型
model.solve()
# 输出结果
if model.status == pulp.LpStatusOptimal:
print("Optimal solution found.")
for i in range(1, 5):
for j in range(1, 5):
print(f"x[{i},{j}] = {pulp.value(x[i, j])}")
print(f"Total cost: {pulp.value(model.objective)}")
else:
print("Failed to find an optimal solution.")
Optimal solution found.
x[1,1] = 3.0 x[1,2] = 0.0 x[1,3] = 0.0
x[1,4] = 12.0 x[2,1] = 0.0 x[2,2] = 0.0
x[2,3] = 0.0 x[2,4] = 0.0 x[3,1] = 8.0
x[3,2] = 0.0 x[3,3] = 0.0 x[3,4] = 0.0
x[4,1] = 0.0 x[4,2] = 0.0 x[4,3] = 0.0
x[4,4] = 0.0
Total cost: 118400.0