文章目录
- 函数极值与规划模型
- 1. 线性代数和线性规划的联系
- 1.1 线性代数的基本概念
- 1.2 线性规划的基本概念
- 1.3 线性代数与线性规划的联系
- 矩阵和向量
- 线性方程组
- 单纯形法
- 内点法
- 凸优化
- 1.4 例子
- 2. Numpy有关矩阵运算示例
- 2.1 矩阵的创建
- 2.2 矩阵的基本运算
- 2.3 矩阵的合并
- 2.4 矩阵的分割
- 3. 线性规划工具包介绍
- 3.1 scipy的linprog函数
- 示例代码
- 3.2 PuLP工具包
- 安装PuLP
- 使用示例
- 参数说明
- 4. 线性规划的建模案例
- 4.1 油料的采购与加工计划
- 实现这个问题的代码如下:
- 代码详细解释
- 结果解释
- 4.2 农民承包土地问题
- 模型建立
- 决策变量
- 目标函数
- 约束条件
- 代码实现
- 参数解释
- 计算结果
- 补充pulp工具包的解题代码
- 5. 非线性规划工具包介绍
- 5.1 SciPy的`minimize`函数
- 安装SciPy
- 使用示例
- 参数说明
- 返回值
- 5.2 Pyomo
- 安装Pyomo
- 使用示例
- 参数说明
- 返回值
- 5.3 IPOPT
- 安装IPOPT
- 使用示例
- 5.4 CVXPY
- 安装CVXPY
- 使用示例
- 参数说明
- 返回值
- 5.5 各类工具包对比
- 6. 非线性规划的建模案例
- 6.1 工地选址
- 解题思路
- 数学建模公式
- 代码解析
- 6.2 职称晋级与评审规划
- 建立模型
- 1) 总工资预算约束
- 2) 编制人数约束
- 3) 晋升人数约束
- 目标函数
- 解题代码
- 理解松弛变量 d n − d_n^- dn− 和 d n + d_n^+ dn+
- 7. 整数规划与指派问题
- 7.1 整数规划的基本概念
- 7.2 分支定界法
- 7.3 指派问题与匈牙利法
- 8. 一些拓展
- 8.1 动态规划与贪心算法
- 8.2 博弈论与排队论
- 8.3 多目标规划
- 8.4 Monte Carlo模拟
- 9. 一些奇奇怪怪的知识点
- 9.1 半正定矩阵和正定矩阵分别是什么,两者有什么联系和区别
- 基本概念
- 联系与区别
- 总结
- 9.2 线性函数和非线性函数的区别
- 定义和基本形式
- 超位性和齐次性
- 求解方法
- 实际应用
- 9.3 分支定界法的时间复杂度
- 分支定界法的时间复杂度因素
- 指数级复杂度的原因
- 剪枝的影响
- 复杂度示例
- 10. 学习心得
函数极值与规划模型
1. 线性代数和线性规划的联系
1.1 线性代数的基本概念
线性代数是数学的一个分支,主要研究向量、向量空间(或线性空间)、线性变换以及线性方程组。它提供了用于描述和操作线性系统的工具。以下是一些线性代数的基本概念:
- 向量:具有方向和大小的量,可以表示为一维数组。
- 矩阵:矩阵是由向量组成的二维数组,用于表示线性变换和线性方程组。
- 线性方程组:由多个线性方程组成的方程组,可以用矩阵形式表示。
- 行列式:矩阵的一个标量值属性,反映了矩阵的一些几何性质,如是否可逆。
- 特征值和特征向量:矩阵的一些固有属性,反映了线性变换的一些重要性质。
1.2 线性规划的基本概念
线性规划是一种优化技术,旨在最大化或最小化一个线性目标函数,同时满足一组线性约束条件。以下是线性规划的基本概念:
- 目标函数:需要优化的线性函数,通常表示为 (c^T x),其中 (c) 是系数向量,(x) 是决策变量向量。
- 约束条件:限制决策变量的线性不等式或等式,通常表示为 (Ax \leq b) 或 (Ax = b),其中 (A) 是系数矩阵,(b) 是常数向量。
- 可行解空间:满足所有约束条件的决策变量的集合。
- 最优解:在可行解空间中使目标函数达到最大或最小值的解。
1.3 线性代数与线性规划的联系
矩阵和向量
在线性规划中,决策变量和约束条件通常用向量和矩阵表示。例如,目标函数 (c^T x) 和约束条件 (Ax \leq b) 都使用了线性代数中的矩阵和向量表示法。这使得线性规划问题可以形式化为矩阵运算问题,从而利用线性代数的方法来分析和求解。
线性方程组
线性规划中的约束条件通常包括线性方程组或线性不等式组。解决这些方程组和不等式组是线性规划的核心任务之一。线性代数提供了求解线性方程组的工具,如高斯消元法和矩阵分解方法。
单纯形法
单纯形法是一种用于求解线性规划问题的算法,基于线性代数的基本原理。它通过在约束条件定义的多面体(可行解空间)的顶点之间移动,逐步寻找最优解。这个过程涉及大量的矩阵运算和线性代数操作,如矩阵的行变换和基变换。
内点法
内点法是另一种解决线性规划问题的算法,它通过在可行解空间内部寻找路径来逼近最优解。内点法使用了线性代数中的矩阵分解和求解技术,以有效地处理大规模的线性规划问题。
凸优化
线性规划是凸优化的一种特例。线性规划问题的可行解空间是一个凸多面体,而线性代数提供了分析凸集和凸函数的工具。这些工具对于理解线性规划的几何性质和设计求解算法非常重要。
1.4 例子
以下是一个简单的线性规划问题,展示了如何使用线性代数方法来求解:
问题:
最大化 z = 3 x 1 + 2 x 2 z = 3x_1 + 2x_2 z=3x1+2x2
约束条件:
{ 2 x 1 + x 2 ≤ 20 4 x 1 − 5 x 2 ≥ − 10 x 1 + 2 x 2 ≤ 15 x 1 , x 2 ≥ 0 \begin{cases} 2x_1 + x_2 \leq 20 \\ 4x_1 - 5x_2 \geq -10 \\ x_1 + 2x_2 \leq 15 \\ x_1, x_2 \geq 0 \end{cases} ⎩ ⎨ ⎧2x1+x2≤204x1−5x2≥−10x1+2x2≤15x1,x2≥0
步骤:
- 将约束条件转换为标准形式:
{ 2 x 1 + x 2 ≤ 20 4 x 1 − 5 x 2 ≤ 10 x 1 + 2 x 2 ≤ 15 x 1 , x 2 ≥ 0 \begin{cases} 2x_1 + x_2 \leq 20 \\ 4x_1 - 5x_2 \leq 10 \\ x_1 + 2x_2 \leq 15 \\ x_1, x_2 \geq 0 \end{cases} ⎩ ⎨ ⎧2x1+x2≤204x1−5x2≤10x1+2x2≤15x1,x2≥0 - 使用线性代数方法(如单纯形法或内点法)求解。
2. Numpy有关矩阵运算示例
直接上代码:
2.1 矩阵的创建
# 创建数组a = np.array([1,2,3])# 指定数据类型a = np.array([1,2,3],dtype=np.int64)# 创建二维数组a = np.array([[1,2,3],[4,5,6]])# 创建多维数组a = np.array([[[1,2,3],[4,5,6]],[[1,2,3],[4,5,6]],[[1,2,3],[4,5,6]]]) # 3x2x3print(a.shape)# 创建数据全为0a = np.zeros((1,3))# 创建数据全为1a = np.ones((1,3))# 创建数据接近0a = np.empty((1,3))# 创建一个从10到20的数组,步长为1a = np.arange(10, 20, 1)print("数组 a:", a)# 创建一个从10到20的数组,步长为2b = np.arange(10, 20, 2)print("数组 b:", b)# 创建一个从10到20的数组,包含5个均匀分布的点c = np.linspace(10, 20, 5)print("数组 c:", c)# 创建一个从10到20的数组,包含10个均匀分布的点d = np.linspace(10, 20, 10)print("数组 d:", d)# 生成 3x2 的随机数数组# 指定范围low = 10high = 20random_array = np.random.uniform(low, high, (3, 2))print("3x2 的随机数数组:\n", random_array)# 生成 3x2 的随机整数数组random_int_array = np.random.randint(low, high, (3, 2))print("3x2 的随机整数数组:\n", random_int_array)
2.2 矩阵的基本运算
# 创建两个 3x2 的矩阵
A = np.array([[1, 2], [3, 4], [5, 6]])
B = np.array([[6, 5], [4, 3], [2, 1]])print("矩阵 A:\n", A)
print("矩阵 B:\n", B)# 矩阵加法
C = A + B
print("矩阵 A + B:\n", C)# 矩阵减法
D = A - B
print("矩阵 A - B:\n", D)# 矩阵乘法(元素乘法)
E = A * B
print("矩阵 A * B (元素乘法):\n", E)# 矩阵点乘
F = A.dot(B.T)
print("矩阵 A 和 B.T 的点乘:\n", F)# 矩阵转置
A_T = A.T
print("矩阵 A 的转置:\n", A_T)
2.3 矩阵的合并
# 创建一个 2x2 方阵
G = np.array([[1, 2], [3, 4]])
G_inv = np.linalg.inv(G)
print("矩阵 G:\n", G)
print("矩阵 G 的逆:\n", G_inv)# 垂直合并
H = np.vstack((A, B))
print("垂直合并 A 和 B:\n", H)# 水平合并
I = np.hstack((A, B))
print("水平合并 A 和 B:\n", I)
2.4 矩阵的分割
# 垂直分割
J1, J2, J3 = np.vsplit(A, 3)
print("垂直分割 A:\n", J1, "\n", J2, "\n", J3)# 水平分割
K1, K2 = np.hsplit(A, 2)
print("水平分割 A:\n", K1, "\n", K2)
3. 线性规划工具包介绍
scipy.optimize.linprog
是一个简单易用的函数,适用于基本的线性规划问题。- PuLP 提供了更强大的功能和灵活性,适用于复杂的线性规划问题。
下面是两者的详细介绍和使用示例
3.1 scipy的linprog函数
res=linprog(c=c, A_ub=A, b_ub=b, A_eq=aeq, b_eq=beq, bounds=bounds)
参数说明:
c
: 目标函数的系数向量。我们希望最小化目标函数 (c \cdot x)。A_ub
: 不等式约束矩阵(小于等于)。每一行代表一个不等式约束的系数。b_ub
: 不等式约束向量。每个元素代表对应不等式约束的右侧常数。A_eq
: 等式约束矩阵。每一行代表一个等式约束的系数。b_eq
: 等式约束向量。每个元素代表对应等式约束的右侧常数。bounds
: 变量的取值范围。每个元组表示一个变量的下限和上限。
返回值:
res
: 返回一个优化结果对象,其中包含以下属性:x
: 最优解。fun
: 目标函数在最优解处的值。success
: 优化是否成功的布尔值。status
: 优化的状态代码。message
: 描述优化状态的消息。
示例代码
问题描述:
最大化目标函数:
z = − x 1 + 4 x 2 z = -x_1 + 4x_2 z=−x1+4x2
约束条件:
{ − 3 x 1 + x 2 ≤ 6 x 1 + 2 x 2 ≤ 4 2 x 1 + x 2 = 1 x 1 , x 2 ≥ 0 \begin{cases} -3x_1 + x_2 \leq 6 \\ x_1 + 2x_2 \leq 4 \\ 2x_1 + x_2 = 1 \\ x_1, x_2 \geq 0 \end{cases} ⎩ ⎨ ⎧−3x1+x2≤6x1+2x2≤42x1+x2=1x1,x2≥0
转换为标准形式进行求解:
import numpy as np
from scipy.optimize import linprog# 定义目标函数的系数向量
c = [-1, 4] # 原目标函数是最大化 -x1 + 4x2,但 linprog 求解的是最小化问题# 定义不等式约束矩阵和向量
A_ub = [[-3, 1], [1, 2]]
b_ub = [6, 4]# 定义等式约束矩阵和向量
A_eq = [[2, 1]]
b_eq = [1]# 定义变量的取值范围
x0_bounds = (0, None)
x1_bounds = (0, None)# 调用 linprog 函数求解
res = linprog(c, A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq=b_eq, bounds=[x0_bounds, x1_bounds])# 打印结果
print('Optimal value:', res.fun, '\nX:', res.x)
3.2 PuLP工具包
PuLP 是一个用于线性规划的 Python 库,它提供了定义和求解线性规划问题的简单接口。与 scipy.optimize.linprog
相比,PuLP 更加灵活,适用于更复杂的线性规划问题。
安装PuLP
可以使用以下命令安装 PuLP:
pip install pulp
使用示例
问题描述:
最小化目标函数:
z = x 1 + 4 x 2 z = x_1 + 4x_2 z=x1+4x2
约束条件:
{ − 3 x 1 + x 2 ≤ 6 x 1 + 2 x 2 ≤ 4 2 x 1 + x 2 = 1 x 1 , x 2 ≥ 0 \begin{cases} -3x_1 + x_2 \leq 6 \\ x_1 + 2x_2 \leq 4 \\ 2x_1 + x_2 = 1 \\ x_1, x_2 \geq 0 \end{cases} ⎩ ⎨ ⎧−3x1+x2≤6x1+2x2≤42x1+x2=1x1,x2≥0
转换为标准形式进行求解:
import pulp# 创建一个最小化问题
model = pulp.LpProblem("Minimize Problem", pulp.LpMinimize)# 定义变量
x1 = pulp.LpVariable('x1', lowBound=0) # x1 >= 0
x2 = pulp.LpVariable('x2', lowBound=0) # x2 >= 0# 定义目标函数
model += x1 + 4 * x2, "Objective Function"# 定义约束条件
model += -3 * x1 + x2 <= 6, "Constraint 1"
model += x1 + 2 * x2 <= 4, "Constraint 2"
model += 2 * x1 + x2 == 1, "Constraint 3"# 求解问题
model.solve()# 输出结果
print("Status:", pulp.LpStatus[model.status])
print("Optimal value of x1:", x1.varValue)
print("Optimal value of x2:", x2.varValue)
print("Optimal value of the objective function:", pulp.value(model.objective))
参数说明
pulp.LpProblem
: 创建一个线性规划问题,可以指定问题名称和类型(最小化或最大化)。pulp.LpVariable
: 定义决策变量,可以指定变量名和范围。+=
: 添加目标函数或约束条件。model.solve()
: 求解线性规划问题。pulp.LpStatus
: 返回求解状态。varValue
: 返回变量的最优解值。pulp.value
: 返回目标函数在最优解处的值。
4. 线性规划的建模案例
4.1 油料的采购与加工计划
某加工厂加工一种油,原料为五种油(植物油1、植物油2、非植物油1、非植物油2、非植物油3),每种油的价格和硬度如图表所示。最终生产的成品将以150英镑/吨的价格卖出。每个月能够提炼的植物油不超过200吨,非植物油不超过250吨。假设提炼过程中油料没有损失,提炼费用忽略不计,并且最终的产品的硬度需要在(3-6)之间(假设硬度的混合时线性的)。根据以上信息,制定月采购和加工计划。
假设 x 1 , x 2 , x 3 , x 4 , x 5 x_1, x_2, x_3, x_4, x_5 x1,x2,x3,x4,x5 分别为每月需要采购的原料油吨数, x 6 x_6 x6 为每个月加工的成品油吨数,由于不考虑油料损失,存在关系:
x 6 = x 1 + x 2 + x 3 + x 4 + x 5 (3.4.2) x_6 = x_1 + x_2 + x_3 + x_4 + x_5 \tag{3.4.2} x6=x1+x2+x3+x4+x5(3.4.2)
平均的硬度为:
η = ∑ i = 1 5 w i x i x 6 (3.4.3) \eta = \frac{\sum_{i=1}^{5} w_i x_i}{x_6} \tag{3.4.3} η=x6∑i=15wixi(3.4.3)
加上植物油重量限制和非植物油重量限制,根据题意,可以列出规划模型如下:
minimize z = − 110 x 1 − 120 x 2 − 130 x 3 − 110 x 4 − 115 x 5 + 150 x 6 s.t. x 1 + x 2 ≤ 200 x 3 + x 4 + x 5 ≤ 250 8.8 x 1 + 6.1 x 2 + 2.0 x 3 + 4.2 x 4 + 5.0 x 5 ≤ 6 x 6 8.8 x 1 + 6.1 x 2 + 2.0 x 3 + 4.2 x 4 + 5.0 x 5 ≥ 3 x 6 x 1 + x 2 + x 3 + x 4 + x 5 = x 6 x i ≥ 0 , i = 1 , 2 , … , 6 (3.4.4) \begin{align} \text{minimize}~& z = -110x_1 - 120x_2 - 130x_3 - 110x_4 - 115x_5 + 150x_6 \\[0.5em] \text{s.t.}~ & x_1 + x_2 \leq 200 \\ & x_3 + x_4 + x_5 \leq 250 \\ & 8.8x_1 + 6.1x_2 + 2.0x_3 + 4.2x_4 + 5.0x_5 \leq 6x_6 \\ & 8.8x_1 + 6.1x_2 + 2.0x_3 + 4.2x_4 + 5.0x_5 \geq 3x_6 \\ & x_1 + x_2 + x_3 + x_4 + x_5 = x_6 \\ & x_i \geq 0, \quad i = 1, 2, \dots, 6 \end{align} \tag{3.4.4} minimize s.t. z=−110x1−120x2−130x3−110x4−115x5+150x6x1+x2≤200x3+x4+x5≤2508.8x1+6.1x2+2.0x3+4.2x4+5.0x5≤6x68.8x1+6.1x2+2.0x3+4.2x4+5.0x5≥3x6x1+x2+x3+x4+x5=x6xi≥0,i=1,2,…,6(3.4.4)
实现这个问题的代码如下:
from scipy.optimize import linprog# 定义目标函数系数
c = [110, 120, 130, 110, 115, -150]# 定义不等式约束矩阵
A = [[1, 1, 0, 0, 0, 0], # 植物油重量限制[0, 0, 1, 1, 1, 0], # 非植物油重量限制[8.8, 6.1, 2.0, 4.2, 5.0, -6], # 硬度上限[-8.8, -6.1, -2.0, -4.2, -5.0, 3] # 硬度下限
]# 定义不等式约束右侧常数
b = [200, 250, 0, 0]# 定义等式约束矩阵
aeq = [[1, 1, 1, 1, 1, -1]] # 总重量等于成品油重量# 定义等式约束右侧常数
beq = [0]# 定义变量的取值范围
bounds = [(0, None), (0, None), (0, None), (0, None), (0, None), (0, 450)]# 求解线性规划问题
res = linprog(c=c, A_ub=A, b_ub=b, A_eq=aeq, b_eq=beq, bounds=bounds)# 打印结果
print(f'Optimal solution: {res.x}')
print(f'Maximum profit: {-res.fun}') # 由于目标函数是取负,因此要取反得到最大利润
代码详细解释
-
定义目标函数系数:
c = [110, 120, 130, 110, 115, -150]
:目标函数的系数,其中最后一项 x 6 x_6 x6 的系数是 -150(因为我们要最大化利润,因此在求解时使用负号)。
-
定义不等式约束矩阵:
A = [...]
:包含四个不等式约束:[1, 1, 0, 0, 0, 0]
:植物油重量限制 x 1 + x 2 ≤ 200 x_1 + x_2 \leq 200 x1+x2≤200。[0, 0, 1, 1, 1, 0]
:非植物油重量限制 x 3 + x 4 + x 5 ≤ 250 x_3 + x_4 + x_5 \leq 250 x3+x4+x5≤250。[8.8, 6.1, 2.0, 4.2, 5.0, -6]
:硬度上限 8.8 x 1 + 6.1 x 2 + 2.0 x 3 + 4.2 x 4 + 5.0 x 5 ≤ 6 x 6 8.8x_1 + 6.1x_2 + 2.0x_3 + 4.2x_4 + 5.0x_5 \leq 6x_6 8.8x1+6.1x2+2.0x3+4.2x4+5.0x5≤6x6。[-8.8, -6.1, -2.0, -4.2, -5.0, 3]
:硬度下限 8.8 x 1 + 6.1 x 2 + 2.0 x 3 + 4.2 x 4 + 5.0 x 5 ≥ 3 x 6 8.8x_1 + 6.1x_2 + 2.0x_3 + 4.2x_4 + 5.0x_5 \geq 3x_6 8.8x1+6.1x2+2.0x3+4.2x4+5.0x5≥3x6。
-
定义不等式约束右侧常数:
b = [200, 250, 0, 0]
:对应四个不等式约束的右侧常数。
-
定义等式约束矩阵:
aeq = [[1, 1, 1, 1, 1, -1]]
:总重量等于成品油重量 x 1 + x 2 + x 3 + x 4 + x 5 = x 6 x_1 + x_2 + x_3 + x_4 + x_5 = x_6 x1+x2+x3+x4+x5=x6。
-
定义等式约束右侧常数:
beq = [0]
:等式约束右侧常数。
-
定义变量的取值范围:
bounds = [(0, None), (0, None), (0, None), (0, None), (0, None), (0, 450)]
:定义各变量的取值范围,所有 x i x_i xi 非负,其中 x 6 x_6 x6 上限为 450。
-
求解线性规划问题:
- 使用
linprog
函数求解线性规划问题,返回结果res
。
- 使用
-
打印结果:
print(f'Optimal solution: {res.x}')
:打印最优解,即各油料采购量。print(f'Maximum profit: {-res.fun}')
:打印最大利润(取负号)。
结果解释
代码执行后,五种原料油的采购量分别为 [159.25, 40.7407, 0, 250, 0]
吨,此时总利润可以达到最大,约为 17592 英镑/月。这展示了通过线性规划方法优化采购和加工计划,以实现利润最大化的过程。
4.2 农民承包土地问题
一个农民承包了6块耕地共300亩,准备播种小麦、玉米、水果和蔬菜四种农产品,各种农产品的计划播种面积、每块土地种植不同农产品的单产收益如下表:
问如何安排种植计划,可得到最大收益。
这是一个产销平衡的运输问题。可以建立下列的运输模型:
对上面问题建立模型,并且给出代码
要解决这个农民的种植计划问题,我们可以将其视为一个线性规划问题,目标是最大化收益。在这个问题中,我们有以下信息:
- 每块地的面积和各块地的面积总和。
- 每块地种植不同作物的收益。
- 每种作物的计划种植面积。
基于这些信息,我们可以建立线性规划模型。设 x i j x_{ij} xij 表示在地块 i i i 上种植作物 j j j 的面积,目标是最大化总收益。
模型建立
决策变量
设 x i j x_{ij} xij 表示在种植作物 i i i 在地块 j j j的种植面积,其中 i = 1 , 2 , … , 4 i = 1, 2, \ldots, 4 i=1,2,…,4 分别表示小麦、玉米、水果和蔬菜, j = 1 , 2 , … , 6 j = 1, 2, \ldots, 6 j=1,2,…,6 表示地块。
目标函数
最大化总收益:
max Z = ∑ i = 1 4 ∑ j = 1 6 p i j x i j \text{max } Z = \sum_{i=1}^{4} \sum_{j=1}^{6} p_{ij} x_{ij} max Z=i=1∑4j=1∑6pijxij
其中, p i j p_{ij} pij 表示种植作物 i i i 在地块 i i i 上的单产收益。
约束条件
- 每块地的面积约束:
∑ i = 1 4 x i j ≤ 地块 j 的面积 , j = 1 , 2 , … , 6 \sum_{i=1}^{4} x_{ij} \leq \text{地块 }j\text{ 的面积}, \quad j = 1, 2, \ldots, 6 i=1∑4xij≤地块 j 的面积,j=1,2,…,6
- 每种作物的种植面积约束:
∑ j = 1 6 x i j = 作物 i 的计划种植面积 , i = 1 , 2 , 3 , 4 \sum_{j=1}^{6} x_{ij} = \text{作物 }i\text{ 的计划种植面积}, \quad i = 1, 2, 3, 4 j=1∑6xij=作物 i 的计划种植面积,i=1,2,3,4
- 非负约束:
x i j ≥ 0 , ∀ i , j x_{ij} \geq 0, \quad \forall i, j xij≥0,∀i,j
代码实现
以下是使用 scipy.optimize.linprog
来求解该问题的代码:
from scipy.optimize import linprog
import numpy as np# 定义单产收益矩阵 p_ijp = np.array([[500, 800, 1000, 1200],[550, 700, 960, 1040],[630, 600, 840, 980],[1000, 950, 650, 860],[800, 900, 600, 880],[700, 930, 700, 780]]).T# 地块面积land_area = np.array([42, 56, 44, 39, 60, 59])# 作物计划种植面积crop_area = np.array([76, 88, 96, 40])# 将收益矩阵展平c = -p.flatten()# 不等式约束:每块地的总种植面积不超过其面积A_ub = np.zeros((6, 24))for j in range(6):A_ub[j,j::6] = 1b_ub = land_area# 等式约束:每种作物的总种植面积等于其计划面积A_eq = np.zeros((4, 24))for i in range(4):A_eq[i, i * 6:(i + 1) * 6] = 1b_eq = crop_area# 变量范围bounds = [(0, None) for _ in range(24)]# 求解线性规划问题res = linprog(c, A_ub=A_ub, b_ub=b_ub, A_eq=A_eq, b_eq=b_eq, bounds=bounds, method='highs')# 结果展示if res.success:x = res.x.reshape((4, 6))print("最优种植计划 (每块地种植的面积):")print(x)print(f"最大收益: {-res.fun}")else:print("未找到最优解")
参数解释
c
: 目标函数的系数向量,包含所有地块和作物组合的收益(取负号是因为 linprog 默认最小化)。A_ub
: 不等式约束矩阵,每行对应一个地块,每列对应一个组合。b_ub
: 不等式约束的右侧常数,每个元素对应地块的面积。A_eq
: 等式约束矩阵,每行对应一个作物,每列对应一个组合。b_eq
: 等式约束的右侧常数,每个元素对应作物的计划种植面积。bounds
: 决策变量的范围,每个变量的取值范围为 [0, ∞)。
该代码通过求解线性规划模型来确定每块地应种植的作物面积,以实现最大收益。
计算结果
最优种植计划 (每块地种植的面积):
[[ 0. 0. 6. 39. 31. 0.][ 0. 0. 0. 0. 29. 59.][ 2. 56. 38. 0. 0. 0.][40. 0. 0. 0. 0. 0.]]
最大收益: 284230.0
补充pulp工具包的解题代码
import pulp
import numpy as np
def transportation_problem(costs, x_max, y_max):row = len(costs)col = len(costs[0])prob = pulp.LpProblem('Transportation Proble',sense=pulp.LpMaximize)var = [[pulp.LpVariable(f'x{i}{j}',lowBound=0,cat=pulp.LpInteger) for j in range(col)] for i in range(row)]# 转为一维flatten = lambda x:[y for l in x for y in flatten(l)] if type(x) is list else [x]prob += pulp.lpDot(flatten(var),costs.flatten())for i in range(row):prob += (pulp.lpSum(var[i]) <= x_max[i])for j in range(col):prob += (pulp.lpSum([var[i][j] for i in range(row)]) <= y_max[j])prob.solve()return {'objective':pulp.value(prob.objective),'var':[[pulp.value(var[i][j]) for j in range(col)] for i in range(row)]}
costs = np.array([[500,550,630,1000,800,700],[800,700,600,950,900,930],[1000,960,840,650,600,700],[1200,1040,980,860,880,780]])
max_plant = [76,88,96,40]
max_cultivation = [42,56,44,39,60,59]
res = transportation_problem(costs, max_plant, max_cultivation)
print(f'最大值为{res["objective"]}')
print("各个变量的取值为:")
print(res['var'])
# 最大值为284230.0
# 各变量的取值为:
# [[0.0, 0.0, 6.0, 39.0, 31.0, 0.0],
# [0.0, 0.0, 0.0, 0.0, 29.0, 59.0],
# [2.0, 56.0, 38.0, 0.0, 0.0, 0.0],
# [40.0, 0.0, 0.0, 0.0, 0.0, 0.0]]
5. 非线性规划工具包介绍
非线性规划(Nonlinear Programming,NLP)涉及优化非线性目标函数,同时满足一组非线性约束。由于非线性问题的复杂性,相较于线性规划,求解非线性规划问题的方法更多样且更复杂。以下是一些常用的非线性规划工具包及其详细内容:
5.1 SciPy的minimize
函数
SciPy 库中的 minimize
函数是一个通用的优化器,适用于多种优化问题,包括无约束和有约束的非线性规划问题。
安装SciPy
pip install scipy
使用示例
问题描述:
最小化目标函数:
f ( x , y ) = ( x − 1 ) 2 + ( y − 2 ) 2 f(x, y) = (x - 1)^2 + (y - 2)^2 f(x,y)=(x−1)2+(y−2)2
约束条件:
{ x 2 + y 2 ≤ 1 x + y = 1 \begin{cases} x^2 + y^2 \leq 1 \\ x + y = 1 \end{cases} {x2+y2≤1x+y=1
import numpy as np
from scipy.optimize import minimize# 定义目标函数
def objective(x):return (x[0] - 1)**2 + (x[1] - 2)**2# 定义约束条件
def constraint1(x):return x[0]**2 + x[1]**2 - 1def constraint2(x):return x[0] + x[1] - 1# 初始猜测
x0 = [0.5, 0.5]# 定义约束
con1 = {'type': 'ineq', 'fun': constraint1}
con2 = {'type': 'eq', 'fun': constraint2}
cons = [con1, con2]# 定义变量的取值范围
bnds = [(0, None), (0, None)]# 使用 minimize 函数求解
solution = minimize(objective, x0, method='SLSQP', bounds=bnds, constraints=cons)# 打印结果
print('Optimal value:', solution.fun)
print('Optimal solution:', solution.x)
参数说明
objective
: 目标函数。x0
: 初始猜测。method
: 优化算法(例如,SLSQP
、trust-constr
等)。bounds
: 变量的取值范围。constraints
: 约束条件。
返回值
solution
: 包含优化结果的对象,属性包括fun
(目标函数值)、x
(最优解)等。
5.2 Pyomo
Pyomo 是一个强大的 Python 库,用于定义和求解各种优化问题,包括线性规划、非线性规划、混合整数规划等。
安装Pyomo
pip install pyomo
使用示例
问题描述:
最小化目标函数:
f ( x , y ) = ( x − 1 ) 2 + ( y − 2 ) 2 f(x, y) = (x - 1)^2 + (y - 2)^2 f(x,y)=(x−1)2+(y−2)2
约束条件:
{ x 2 + y 2 ≤ 1 x + y = 1 \begin{cases} x^2 + y^2 \leq 1 \\ x + y = 1 \end{cases} {x2+y2≤1x+y=1
from pyomo.environ import ConcreteModel, Var, Objective, Constraint, SolverFactory# 创建模型
model = ConcreteModel()# 定义变量
model.x = Var(bounds=(0, None))
model.y = Var(bounds=(0, None))# 定义目标函数
model.obj = Objective(expr=(model.x - 1)**2 + (model.y - 2)**2, sense=1)# 定义约束条件
model.con1 = Constraint(expr=model.x**2 + model.y**2 <= 1)
model.con2 = Constraint(expr=model.x + model.y == 1)# 求解模型
solver = SolverFactory('ipopt')
solver.solve(model)# 打印结果
print('Optimal value:', model.obj())
print('Optimal solution: x =', model.x(), ', y =', model.y())
参数说明
ConcreteModel
: 创建一个具体模型实例。Var
: 定义变量。Objective
: 定义目标函数。Constraint
: 定义约束条件。SolverFactory
: 指定求解器。
返回值
model.obj()
: 目标函数在最优解处的值。model.x()
,model.y()
: 变量的最优解值。
5.3 IPOPT
IPOPT(Interior Point OPTimizer)是一个用于大规模非线性优化的求解器,特别适合处理具有非线性约束的问题。
安装IPOPT
IPOPT 通常与 Pyomo 一起使用,安装命令如下:
pip install ipopt
然后,需要单独下载并安装 IPOPT 可执行文件,并将其路径添加到系统路径中。
使用示例
上面的 Pyomo 示例已经展示了如何使用 IPOPT 求解非线性规划问题。只需在求解器部分指定 ipopt
即可:
solver = SolverFactory('ipopt')
solver.solve(model)
5.4 CVXPY
CVXPY 是一个用于构建和求解凸优化问题的库,但它也可以处理某些非凸问题。CVXPY 提供了一种简单直观的方式来定义优化问题,并使用高效的求解器来找到解。
安装CVXPY
pip install cvxpy
使用示例
问题描述:
最小化目标函数:
f ( x , y ) = ( x − 1 ) 2 + ( y − 2 ) 2 f(x, y) = (x - 1)^2 + (y - 2)^2 f(x,y)=(x−1)2+(y−2)2
约束条件:
{ x 2 + y 2 ≤ 1 x + y = 1 \begin{cases} x^2 + y^2 \leq 1 \\ x + y = 1 \end{cases} {x2+y2≤1x+y=1
import cvxpy as cp# 定义变量
x = cp.Variable()
y = cp.Variable()# 定义目标函数
objective = cp.Minimize((x - 1)**2 + (y - 2)**2)# 定义约束条件
constraints = [x**2 + y**2 <= 1, x + y == 1]# 定义问题
problem = cp.Problem(objective, constraints)# 求解问题
problem.solve()# 打印结果
print('Optimal value:', problem.value)
print('Optimal solution: x =', x.value, ', y =', y.value)
参数说明
cp.Variable()
: 定义优化变量。cp.Minimize()
: 定义最小化目标函数。cp.Problem()
: 定义优化问题,包括目标函数和约束条件。problem.solve()
: 求解优化问题。
返回值
problem.value
: 目标函数在最优解处的值。x.value
,y.value
: 变量的最优解值。
5.5 各类工具包对比
scipy.optimize.minimize
是一个通用的优化函数,适用于各种非线性规划问题。- Pyomo 提供了定义和求解复杂优化问题的强大功能,适用于大规模和复杂的非线性规划问题。
- IPOPT 是一个高效的非线性求解器,特别适合处理大规模非线性约束问题。
- CVXPY 是一个用于构建和求解凸优化问题的库,但它也可以处理某些非凸问题。
选择哪个工具包取决于具体问题的复杂程度和求解需求。对于简单的非线性问题
scipy.optimize.minimize
可能足够;对于复杂和大规模的问题,Pyomo 和 IPOPT 提供了更强大的功能和灵活性。CVXPY 是一个用于凸优化问题的 Python 库,尽管其设计初衷是解决凸优化问题,但它也可以用于某些非线性规划问题,特别是那些可以表达为凸问题的非线性优化问题。
6. 非线性规划的建模案例
6.1 工地选址
某公司有 6 6 6个建筑工地要开工,每个工地的位置(用平面坐标系 a , b a,b a,b表示,距离单位:千米)及水泥日用 量 d d d(吨)由下表给出。规划设立两个料场位于 A A A, B B B,日储量各为20吨。假设从料场到工地之间均有 直线道路相连,试确定料场的位置,并制定每天的供应计划,即从 A A A, B B B两料场分别向各工地运送多少吨水泥,使总的吨千米数最小。
规划问题的核心有三样:决策变量,目标函数和约束条件。
决策变量包括哪些?首先两个料场的坐标未知吧,坐标有横纵坐标于是这就有了四个变量;两个料场到六个工地 12 12 12条线路上的运输量也未知吧,于是这就又来了12个变量,一共是 16 16 16个。距离可以用 Euclid 距离来计算。
解题思路
要解决这个问题,我们需要使用线性规划方法,通过设置目标函数和约束条件来优化料场位置和运输计划。具体步骤如下:
-
决策变量:
- 料场A的位置 (x1, y1)
- 料场B的位置 (x2, y2)
- 从料场A到工地i的运输量 (tA_i),i = 1,…,6
- 从料场B到工地i的运输量 (tB_i),i = 1,…,6
-
目标函数:
- 我们的目标是最小化总的吨千米数。即从料场A和B分别到各工地的运输量乘以距离之和:
minimize ∑ i = 1 6 ( t A i ⋅ dist ( A , 工地 i ) + t B i ⋅ dist ( B , 工地 i ) ) \text{minimize} \quad \sum_{i=1}^6 \left( tA_i \cdot \text{dist}(A, 工地i) + tB_i \cdot \text{dist}(B, 工地i) \right) minimizei=1∑6(tAi⋅dist(A,工地i)+tBi⋅dist(B,工地i))
- 我们的目标是最小化总的吨千米数。即从料场A和B分别到各工地的运输量乘以距离之和:
-
约束条件:
- 各工地的水泥需求必须得到满足:
t A i + t B i = d i , i = 1 , 2 , . . . , 6 tA_i + tB_i = d_i, \quad i = 1,2,...,6 tAi+tBi=di,i=1,2,...,6 - 每个料场的运输总量不能超过其日储量:
∑ i = 1 6 t A i ≤ 20 \sum_{i=1}^6 tA_i \leq 20 i=1∑6tAi≤20
∑ i = 1 6 t B i ≤ 20 \sum_{i=1}^6 tB_i \leq 20 i=1∑6tBi≤20 - 所有运输量非负:
t A i ≥ 0 , t B i ≥ 0 , i = 1 , 2 , . . . , 6 tA_i \geq 0, \quad tB_i \geq 0, \quad i = 1,2,...,6 tAi≥0,tBi≥0,i=1,2,...,6
- 各工地的水泥需求必须得到满足:
数学建模公式
我们用Python和优化库(如SciPy或PuLP)来解决这个优化问题。
import numpy as np
from scipy.optimize import minimize# 工地的位置和需求
a = np.array([1.25, 8.75, 0.5, 5.75, 3, 7.25])
b = np.array([1.25, 0.75, 4.75, 5, 6.5, 7.25])
d = np.array([3, 5, 4, 7, 6, 11])# 目标函数
def objective(x):s = 0for i in range(6):s += x[4+i] * np.sqrt((x[0] - a[i])**2 + (x[1] - b[i])**2)s += x[10+i] * np.sqrt((x[2] - a[i])**2 + (x[3] - b[i])**2)return s# 约束条件
constraints = []# 各工地的水泥需求必须得到满足
for i in range(6):constraints.append({'type': 'eq', 'fun': lambda x, i=i: x[4+i] + x[10+i] - d[i]})# 每个料场的运输总量不能超过其日储量
constraints.append({'type': 'ineq', 'fun': lambda x: 20 - sum(x[4:10])})
constraints.append({'type': 'ineq', 'fun': lambda x: 20 - sum(x[10:16])})# 初始猜测
x0 = np.array([0, 0, 0, 0] + [1]*12)# 求解
result = minimize(objective, x0, constraints=constraints)# 结果展示
print("料场A的位置: ({:.2f}, {:.2f})".format(result.x[0], result.x[1]))
print("料场B的位置: ({:.2f}, {:.2f})".format(result.x[2], result.x[3]))
for i in range(6):print("从料场A到工地{}的运输量: {:.2f}吨".format(i+1, result.x[4+i]))print("从料场B到工地{}的运输量: {:.2f}吨".format(i+1, result.x[10+i]))
print("总共运输的最少吨数:{:.2f}吨".format(np.sum(result.x[4:15])))
代码解析
- 目标函数
objective(x)
:计算总的吨千米数。 - 约束条件
constraints
:包括工地需求的等式约束和料场运输总量的不等式约束。 - 初始猜测
x0
:给定一个初始猜测值。 - 求解
minimize(objective, x0, constraints=constraints)
:使用SciPy的minimize
函数进行求解。
通过上述步骤,我们可以求得料场的最优位置和每天的运输计划。
料场A的位置: (5.75, 5.00)
料场B的位置: (5.75, 5.00)
从料场A到工地1的运输量: 1.50吨
从料场B到工地1的运输量: 1.50吨
从料场A到工地2的运输量: 2.51吨
从料场B到工地2的运输量: 2.49吨
从料场A到工地3的运输量: 2.00吨
从料场B到工地3的运输量: 2.00吨
从料场A到工地4的运输量: 3.51吨
从料场B到工地4的运输量: 3.49吨
从料场A到工地5的运输量: 3.00吨
从料场B到工地5的运输量: 3.00吨
从料场A到工地6的运输量: 5.49吨
从料场B到工地6的运输量: 5.51吨
总共运输的最少吨数:{}吨 30.49
6.2 职称晋级与评审规划
建立模型
为了设计符合要求的职工调薪方案,我们需要解决一个线性规划问题。具体而言,我们有以下几个目标和约束条件:
- 总工资预算不超过150万元。
- 每级职工人数不超过编制人数。
- II、III级职工的晋升人数尽量达到现有人数的20%。
我们将问题分解为三个子问题,并用线性规划的方法来求解。
1) 总工资预算约束
设:
- 由II级晋升为I级的人数为 x 1 x_1 x1
- 由III级晋升为II级的人数为 x 2 x_2 x2
- 新招聘的III级职工人数为 x 3 x_3 x3
总工资预算的约束为:
50000 ( 9 + x 1 ) + 30000 ( 12 − x 1 + x 2 ) + 20000 ( 15 − x 2 + x 3 ) ≤ 1500000 50000(9 + x_1) + 30000(12 - x_1+x_2) + 20000(15 - x_2 + x_3) \leq 1500000 50000(9+x1)+30000(12−x1+x2)+20000(15−x2+x3)≤1500000
用松弛变量 d 1 − d_1^- d1− 和 d 1 + d_1^+ d1+ 表示未满误差和过盈误差,我们有:
50000 ( 9 + x 1 ) + 30000 ( 12 − x 1 + x 2 ) + 20000 ( 15 − x 2 + x 3 ) + d 1 − − d 1 + = 1500000 50000(9 + x_1) + 30000(12 - x_1+x_2) + 20000(15 - x_2 + x_3) + d_1^- - d_1^+ = 1500000 50000(9+x1)+30000(12−x1+x2)+20000(15−x2+x3)+d1−−d1+=1500000
2) 编制人数约束
每级人数不超过编制规定人数的约束为:
9 + x 1 + d 2 − − d 2 + = 12 9 + x_1 + d_2^- - d_2^+ = 12 9+x1+d2−−d2+=12
12 − x 1 + x 2 + d 3 − − d 3 + = 15 12 - x_1 + x_2 + d_3^- - d_3^+ = 15 12−x1+x2+d3−−d3+=15
15 − x 2 + x 3 + d 4 − − d 4 + = 15 15 - x_2 + x_3 + d_4^- - d_4^+ = 15 15−x2+x3+d4−−d4+=15
3) 晋升人数约束
为了使II级、III级职工的晋升人数尽量达到现有人数的20%,我们有:
x 1 + d 5 − − d 5 + = 3 x_1 + d_5^- - d_5^+ = 3 x1+d5−−d5+=3
x 2 + d 6 − − d 6 + = 3 x_2 + d_6^- - d_6^+ = 3 x2+d6−−d6+=3
目标函数
目标函数是我们希望最小化的值,它在综合考虑所有约束的偏差后,给出一个最优解。目标函数如下:
minimize p 1 ⋅ d 1 + + p 2 ⋅ ( d 2 + + d 3 + + d 4 + ) + p 3 ⋅ ( d 5 − − d 5 + + d 6 − − d 6 + ) \text{minimize}~ p_1 \cdot d_1^+ + p_2 \cdot (d_2^+ + d_3^+ + d_4^+) + p_3 \cdot (d_5^- - d_5^+ + d_6^- - d_6^+) minimize p1⋅d1++p2⋅(d2++d3++d4+)+p3⋅(d5−−d5++d6−−d6+)
每一项的含义如下:
- p 1 ⋅ d 1 + p_1 \cdot d_1^+ p1⋅d1+:最小化预算超出的部分。即,尽量避免总工资超出预算。
- p 2 ⋅ ( d 2 + + d 3 + + d 4 + ) p_2 \cdot (d_2^+ + d_3^+ + d_4^+) p2⋅(d2++d3++d4+):最小化每一级职工人数超过编制人数的部分。即,尽量避免各级职工人数超过编制人数。
- p 3 ⋅ ( d 5 − − d 5 + + d 6 − − d 6 + ) p_3 \cdot (d_5^- - d_5^+ + d_6^- - d_6^+) p3⋅(d5−−d5++d6−−d6+):最小化晋升比例的偏差。即,尽量使得II级、III级职工的晋升人数达到预期的20%。
解题代码
from scipy.optimize import linprog
import numpy as np# 定义目标函数的系数 (c)
# c 中每个元素对应一个变量的系数,变量依次为:
# x1, x2, x3, d1-, d1+, d2-, d2+, d3-, d3+, d4-, d4+, d5-, d5+, d6-, d6+
# 其中 p1, p2, p3 是目标函数中不同项的权重
p1, p2, p3 = 1, 1, 1 # 可以根据实际情况调整权重
c = np.array([0, 0, 0, 0, p1, 0, p2, 0, p2, 0, p2, p3, p3, -p3, -p3])# 定义等式约束矩阵 (A_eq) 和等式约束向量 (b_eq)
# A_eq 中每一行对应一个等式约束,每一列对应一个变量
# b_eq 中每个元素对应等式约束的右侧常数值A_eq = np.array([# 预算约束:20000(x1) + 10000(x2) + 20000(x3) + d1- - d1+ = 390000[20000, 10000, 20000, -1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],# 编制约束:9 + x1 + d2- - d2+ = 12[1, 0, 0, 0, 0, -1, 1, 0, 0, 0, 0, 0, 0, 0, 0],# 编制约束:12 - x1 + x2 + d3- - d3+ = 15[0, -1, 1, 0, 0, 0, 0, -1, 1, 0, 0, 0, 0, 0, 0],# 编制约束:15 - x2 + x3 + d4- - d4+ = 15[0, 0, -1, 1, 0, 0, 0, 0, 0, -1, 1, 0, 0, 0, 0],# 晋升人数约束:x1 + d5- - d5+ = 3[1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 1, 0, 0],# 晋升人数约束:x2 + d6- - d6+ = 3[0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 1]
])# 定义等式约束右侧的常数值
b_eq = np.array([390000, 3, 15, 15, 3, 3])# 定义变量的边界
# x1 的范围是 [0, 12],因为现有II级职工人数最多可以升为I级
# x2 的范围是 [0, 15],因为现有III级职工人数最多可以升为II级
# x3 的范围是 [0, 15],因为现有III级职工最多可以招聘至编制上限
# d1-, d1+, d2-, d2+, d3-, d3+, d4-, d4+, d5-, d5+, d6-, d6+ 的范围是 [0, 10],误差范围设置为 0 到 10
bounds = [(0, 12), # x1(0, 15), # x2(0, 15), # x3(0, 10), # d1-(0, 10), # d1+(0, 10), # d2-(0, 10), # d2+(0, 10), # d3-(0, 10), # d3+(0, 10), # d4-(0, 10), # d4+(0, 10), # d5-(0, 10), # d5+(0, 10), # d6-(0, 10) # d6+
]
# 使用 scipy.optimize.linprog 求解线性规划问题
# 参数说明:
# - c: 目标函数的系数
# - A_eq: 等式约束矩阵
# - b_eq: 等式约束右侧的常数值
# - bounds: 每个变量的取值范围
# - method: 使用的方法,这里选择 'highs',是目前 scipy 中较新的线性规划求解方法
res = linprog(c, A_eq=A_eq, b_eq=b_eq, bounds=bounds, method='highs')# 打印结果
print(res)
理解松弛变量 d n − d_n^- dn− 和 d n + d_n^+ dn+
对于约束条件:
50000 ( 9 + x 1 ) + 30000 ( 12 − x 1 + x 2 ) + 20000 ( 15 − x 2 + x 3 ) ≤ 1500000 50000(9 + x_1) + 30000(12 - x_1+x_2) + 20000(15 - x_2 + x_3) \leq 1500000 50000(9+x1)+30000(12−x1+x2)+20000(15−x2+x3)≤1500000
我们引入松弛变量 d 1 − d_1^- d1− 和 d 1 + d_1^+ d1+ 将其转化为等式:
50000 ( 9 + x 1 ) + 30000 ( 12 − x 1 + x 2 ) + 20000 ( 15 − x 2 + x 3 ) + d 1 − − d 1 + = 1500000 50000(9 + x_1) + 30000(12 - x_1+x_2) + 20000(15 - x_2 + x_3) + d_1^- - d_1^+ = 1500000 50000(9+x1)+30000(12−x1+x2)+20000(15−x2+x3)+d1−−d1+=1500000
- d 1 + d_1^+ d1+ 表示“过盈误差”,即超出预算的部分。它是一个正数或零。
- d 1 − d_1^- d1− 表示“未满误差”,即预算不足的部分。它也是一个正数或零。
在此等式中,如果总工资 50000 ( 9 + x 1 ) + 30000 ( 12 − x 1 + x 2 ) + 20000 ( 15 − x 2 + x 3 ) 50000(9 + x_1) + 30000(12 - x_1+x_2) + 20000(15 - x_2 + x_3) 50000(9+x1)+30000(12−x1+x2)+20000(15−x2+x3) 刚好等于 1500000,那么 d 1 + d_1^+ d1+ 和 d 1 − d_1^- d1− 都为零。如果总工资超出1500000元,那么 d 1 + d_1^+ d1+ 为正,表示超出的金额, d 1 − d_1^- d1− 为零。如果总工资不足1500000元,那么 d 1 − d_1^- d1− 为正,表示不足的金额, d 1 + d_1^+ d1+ 为零。
7. 整数规划与指派问题
7.1 整数规划的基本概念
整数规划(Integer Programming, IP) 是一种特殊的线性规划问题,其中一些或所有决策变量必须取整数值。整数规划可以分为以下几类:
- 纯整数规划:所有决策变量必须取整数值。
- 混合整数规划(Mixed Integer Programming, MIP):部分决策变量取整数值,其他决策变量可以取连续值。
- 二进制整数规划(Binary Integer Programming, BIP):决策变量仅能取0或1的值。
整数规划在实际应用中有广泛的应用,例如在调度问题、资源分配问题、物流问题和资本预算问题中,整数约束往往是必需的。
整数规划的标准形式可以表示为:
minimize c T x subject to A x ≤ b x ∈ Z n 或 x i ∈ { 0 , 1 } , ∀ i \begin{aligned} & \text{minimize} & & c^T x \\ & \text{subject to} & & Ax \leq b \\ & & & x \in \mathbb{Z}^n \quad \text{或} \quad x_i \in \{0, 1\}, \forall i \end{aligned} minimizesubject tocTxAx≤bx∈Zn或xi∈{0,1},∀i
其中, c c c 和 b b b 是已知的向量, A A A 是已知的矩阵, x x x 是待求解的整数向量。
7.2 分支定界法
分支定界法(Branch and Bound, B&B) 是求解整数规划问题的一种常用方法。它是一种系统的搜索方法,通过构建解空间树来探索所有可能的解,但通过剪枝技术避免了对每一个可能的解进行穷举。
基本步骤:
-
初始问题求解:首先忽略整数约束,求解相应的线性规划问题,得到松弛问题的最优解。
-
分支:如果松弛问题的最优解不是整数解,选择一个非整数变量,将问题分解为两个子问题,其中一个子问题约束该变量向下取整,另一个子问题约束该变量向上取整。
-
界定:对于每个子问题,计算其松弛问题的最优解。如果某个子问题的松弛问题解不可行,或其最优值不优于已知的整数解,则舍弃该子问题。
-
剪枝:如果子问题的松弛问题解是整数解,并且优于已知最优整数解,则更新当前最优解。继续搜索其他子问题。
-
终止条件:当所有子问题都被处理或舍弃时,算法终止,当前最优整数解即为最优解。
7.3 指派问题与匈牙利法
指派问题(Assignment Problem) 是一种特殊的整数规划问题,旨在将任务分配给工人、机器或其他资源,以最小化总成本或总时间。
指派问题的标准形式可以表示为:
minimize ∑ i = 1 n ∑ j = 1 n c i j x i j subject to ∑ j = 1 n x i j = 1 , i = 1 , … , n ∑ i = 1 n x j = 1 , j = 1 , … , n x i j ∈ { 0 , 1 } \begin{aligned} & \text{minimize} & & \sum_{i=1}^n \sum_{j=1}^n c_{ij} x_{ij} \\ & \text{subject to} & & \sum_{j=1}^n x_{ij} = 1, \quad i=1, \ldots, n \\ & & & \sum_{i=1}^n x_{j} = 1, \quad j=1, \ldots, n \\ & & & x_{ij} \in \{0, 1\} \end{aligned} minimizesubject toi=1∑nj=1∑ncijxijj=1∑nxij=1,i=1,…,ni=1∑nxj=1,j=1,…,nxij∈{0,1}
其中, c i j c_{ij} cij 是将任务 i i i 分配给工人 j j j 的成本, x i j x_{ij} xij 是二进制变量,当任务 i i i 被分配给工人 j j j 时, x i j = 1 x_{ij} = 1 xij=1,否则 x i j = 0 x_{ij} = 0 xij=0。
匈牙利法(Hungarian Method) 是一种有效求解指派问题的算法,能够在多项式时间内找到最优解。
匈牙利法的步骤:
-
构造初始矩阵:从原始成本矩阵中减去每行的最小值,然后减去每列的最小值。
-
覆盖零:用尽可能少的水平线和垂直线覆盖所有的零。
-
调整矩阵:如果覆盖零的线数等于矩阵的阶数,则找到最优指派;否则,对未覆盖的元素进行调整,重复覆盖零的过程,直到线数等于矩阵的阶数。
-
确定最优解:根据调整后的矩阵确定最优指派。
匈牙利法通过一系列矩阵变换和调整,能够有效找到指派问题的最优解,并广泛应用于资源分配、任务调度等领域。
8. 一些拓展
8.1 动态规划与贪心算法
8.2 博弈论与排队论
8.3 多目标规划
8.4 Monte Carlo模拟
9. 一些奇奇怪怪的知识点
9.1 半正定矩阵和正定矩阵分别是什么,两者有什么联系和区别
基本概念
半正定矩阵和正定矩阵是线性代数中的重要概念,常用于二次型、优化理论等领域。它们都与矩阵的特征值和二次型有关。
- 半正定矩阵 (Positive Semidefinite Matrix)
一个对称矩阵 A ∈ R n × n A \in \mathbb{R}^{n \times n} A∈Rn×n 被称为半正定矩阵,如果对于所有的非零向量 x ∈ R n x \in \mathbb{R}^n x∈Rn,都有:
x T A x ≥ 0 x^T A x \geq 0 xTAx≥0
这意味着该矩阵的所有特征值都是非负的。
- 正定矩阵 (Positive Definite Matrix)
一个对称矩阵 A ∈ R n × n A \in \mathbb{R}^{n \times n} A∈Rn×n 被称为正定矩阵,如果对于所有的非零向量 x ∈ R n x \in \mathbb{R}^n x∈Rn,都有:
x T A x > 0 x^T A x > 0 xTAx>0
这意味着该矩阵的所有特征值都是正的。
联系与区别
-
特征值:
- 半正定矩阵:所有特征值都是非负的,即 λ i ≥ 0 \lambda_i \geq 0 λi≥0。
- 正定矩阵:所有特征值都是正的,即 λ i > 0 \lambda_i > 0 λi>0。
-
二次型:
- 半正定矩阵:对于所有的 x ≠ 0 x \neq 0 x=0, x T A x ≥ 0 x^T A x \geq 0 xTAx≥0。
- 正定矩阵:对于所有的 x ≠ 0 x \neq 0 x=0, x T A x > 0 x^T A x > 0 xTAx>0。
-
定义的强度:
- 正定矩阵是半正定矩阵的一个特例。所有正定矩阵都是半正定的,但并非所有半正定矩阵都是正定的。
-
判定方法:
- 半正定矩阵:检查矩阵的所有特征值是否为非负;或者使用主子式判定法(所有顺序主子式均非负)。
- 正定矩阵:检查矩阵的所有特征值是否为正;或者使用主子式判定法(所有顺序主子式均正)。
- 半正定矩阵示例
A = ( 1 0 0 0 ) A = \begin{pmatrix} 1 & 0 \\ 0 & 0 \end{pmatrix} A=(1000)
特征值为 λ 1 = 1 \lambda_1 = 1 λ1=1 和 λ 2 = 0 \lambda_2 = 0 λ2=0,因此 A A A 是半正定的。
- 正定矩阵示例
B = ( 2 1 1 2 ) B = \begin{pmatrix} 2 & 1 \\ 1 & 2 \end{pmatrix} B=(2112)
特征值为 λ 1 = 3 \lambda_1 = 3 λ1=3 和 λ 2 = 1 \lambda_2 = 1 λ2=1,因此 B B B 是正定的。
总结
半正定矩阵和正定矩阵都涉及矩阵的特征值和二次型。半正定矩阵的特征值非负,正定矩阵的特征值正。正定矩阵是半正定矩阵的一个特例。理解它们之间的联系和区别对深入掌握矩阵理论和应用非常重要。
9.2 线性函数和非线性函数的区别
定义和基本形式
-
线性函数:
- 一个线性函数具有以下基本形式:
f ( x ) = a x + b f(x) = ax + b f(x)=ax+b
其中 a a a 和 b b b 是常数, x x x 是自变量。 - 线性函数的特征是图像是一条直线,其斜率由系数 a a a 决定,截距由 b b b 决定。
- 一般形式的线性函数(多元线性函数)可以表示为:
f ( x 1 , x 2 , … , x n ) = a 1 x 1 + a 2 x 2 + … + a n x n + b f(x_1, x_2, \ldots, x_n) = a_1 x_1 + a_2 x_2 + \ldots + a_n x_n + b f(x1,x2,…,xn)=a1x1+a2x2+…+anxn+b
其中 a 1 , a 2 , … , a n a_1, a_2, \ldots, a_n a1,a2,…,an 是系数, b b b 是常数。
- 一个线性函数具有以下基本形式:
-
非线性函数:
- 一个非线性函数不满足线性函数的特性,其形式可以非常多样。例如:
f ( x ) = a x 2 + b x + c f(x) = ax^2 + bx + c f(x)=ax2+bx+c
或者:
f ( x ) = e x , f ( x ) = log ( x ) , f ( x ) = sin ( x ) f(x) = e^x, \quad f(x) = \log(x), \quad f(x) = \sin(x) f(x)=ex,f(x)=log(x),f(x)=sin(x) - 非线性函数的图像通常不是直线,可以是曲线、指数函数、对数函数、三角函数等。
- 非线性函数在多元情况下的形式也更加复杂,例如:
f ( x 1 , x 2 ) = x 1 2 + x 2 2 f(x_1, x_2) = x_1^2 + x_2^2 f(x1,x2)=x12+x22
或:
f ( x 1 , x 2 ) = e x 1 x 2 f(x_1, x_2) = e^{x_1 x_2} f(x1,x2)=ex1x2
- 一个非线性函数不满足线性函数的特性,其形式可以非常多样。例如:
超位性和齐次性
-
线性函数:
- 满足叠加性(超位性):
f ( x 1 + x 2 ) = f ( x 1 ) + f ( x 2 ) f(x_1 + x_2) = f(x_1) + f(x_2) f(x1+x2)=f(x1)+f(x2) - 满足齐次性:
f ( k x ) = k f ( x ) f(kx) = kf(x) f(kx)=kf(x)
其中 k k k 是常数。 - 因此,线性函数的解可以通过线性组合来求得。
- 满足叠加性(超位性):
-
非线性函数:
- 不满足叠加性和齐次性。例如,对于 f ( x ) = x 2 f(x) = x^2 f(x)=x2:
f ( x 1 + x 2 ) ≠ f ( x 1 ) + f ( x 2 ) f(x_1 + x_2) \neq f(x_1) + f(x_2) f(x1+x2)=f(x1)+f(x2)
f ( k x ) ≠ k f ( x ) f(kx) \neq kf(x) f(kx)=kf(x)
- 不满足叠加性和齐次性。例如,对于 f ( x ) = x 2 f(x) = x^2 f(x)=x2:
求解方法
-
线性函数:
- 可以使用线性代数和线性规划等方法求解。
- 线性问题通常可以通过矩阵运算、线性方程组的解法(如高斯消元法)等方法解决。
-
非线性函数:
- 求解方法更加复杂,通常需要数值方法,例如牛顿法、梯度下降法等。
- 非线性规划问题的求解需要专门的优化算法,如拉格朗日乘数法、约束优化中的KKT条件、遗传算法等。
实际应用
-
线性函数:
- 应用于许多简单和基础的问题,如线性回归、资源分配、网络流量优化等。
- 线性模型在工程、经济、物理等领域有广泛应用。
-
非线性函数:
- 用于描述更复杂的现象,如机器学习中的非线性模型、复杂系统的模拟、生物过程、化学反应、经济学中的非线性动态模型等。
- 非线性模型能够捕捉复杂关系和行为,因此在现代科学和工程中占据重要地位。
线性函数由于其简单性和易解性,适用于许多基础和简单的问题。而非线性函数虽然复杂,但能够描述和解决更多实际中的复杂现象和问题。因此,在实际应用中,选择线性还是非线性函数取决于问题的性质和要求。
9.3 分支定界法的时间复杂度
分支定界法的时间复杂度因素
- 问题规模:变量的数量 N N N。
- 分支数量:每个分支点可能会分解成多个子问题。
- 搜索树的深度:树的深度等于变量的数量 N N N,因为每个变量都可能需要分支。
- 剪枝效果:剪枝能够有效减少需要探索的节点数。
指数级复杂度的原因
分支定界法通过递归地分解问题,构建一棵搜索树。每个节点对应一个子问题,并且每个子问题的解空间都是原问题的一个子集。
- 树的构建:假设每次分支将一个问题分成两个子问题,那么在最坏的情况下,树的每一层的节点数会翻倍。因此,对于 N N N 个变量,最坏情况下树的叶子节点数为 2 N 2^N 2N。
- 分支和搜索:在最坏的情况下,需要对每个可能的解进行探索。虽然剪枝可以减少实际探索的节点数,但在没有有效剪枝的情况下,复杂度仍然是指数级的。
剪枝的影响
剪枝技术通过排除一些不必要的分支来减少计算量。有效的剪枝可以显著降低实际需要探索的节点数,从而提高算法效率。然而,剪枝的效果很大程度上取决于具体问题和剪枝策略的有效性。
复杂度示例
考虑一个简单的二进制整数规划问题,其中每个变量 x i x_i xi 可以取 0 或 1:
minimize c 1 x 1 + c 2 x 2 + … + c N x N subject to ∑ i = 1 N a i j x i ≤ b j , ∀ j ∈ { 1 , 2 , … , M } x i ∈ { 0 , 1 } , ∀ i ∈ { 1 , 2 , … , N } \begin{aligned} & \text{minimize} & & c_1 x_1 + c_2 x_2 + \ldots + c_N x_N \\ & \text{subject to} & & \sum_{i=1}^N a_{ij} x_i \leq b_j, \quad \forall j \in \{1, 2, \ldots, M\} \\ & & & x_i \in \{0, 1\}, \quad \forall i \in \{1, 2, \ldots, N\} \end{aligned} minimizesubject toc1x1+c2x2+…+cNxNi=1∑Naijxi≤bj,∀j∈{1,2,…,M}xi∈{0,1},∀i∈{1,2,…,N}
对于这个问题,分支定界法在最坏情况下需要探索所有可能的 2 N 2^N 2N 个解,即复杂度为 O ( 2 N ) O(2^N) O(2N)。
10. 学习心得
本章学习了线性规划的基本概念、掌握使用Numpy对矩阵进行基本操作(创建、运算、合并、分解)。
处理线性规划问题时,需要对问题进行建模,将线性规划整理为标准形式(程序如matlab底层会将标准形式转换为规范形式进行求解),进行coding得到最优解。
本章的示例中列举两种编程求解方法:scipy的linprog函数 和 pulp工具包
scipy.optimize.linprog
是一个简单易用的函数,适用于基本的线性规划问题。- PuLP 提供了更强大的功能和灵活性,适用于复杂的线性规划问题。
对于非线性规划问题,可以使用SciPy的minimize函数、CVXPY等工具包进行编程求解。
求解4.2 农民承包土地问题时,分别使用scipy的linprog函数 和 pulp工具包对模型进行求解;
在求解问题6.2 职称晋级与评审规划时,掌握了一个新技巧:对于未知范围的参数,可以引入松弛变量 d n − d_n^- dn− 和 d n + d_n^+ dn+建立模型。
本次还学习了一些拓展知识:动态规划与贪心算法、博弈论与排队论、多目标规划和Monte Carlo模拟。
超级喜欢这句话:我们解决问题的思路是将不熟悉的问题,转换成熟悉的问题。