linprog#
- scipy.optimize.linprog(c, A_ub=None, b_ub=None, A_eq=None, b_eq=None, bounds=(0, None), method='highs', callback=None, options=None, x0=None, integrality=None)[source]#
线性规划:在满足线性等式和不等式约束条件下,最小化线性目标函数。
线性规划解决以下形式的问题
\[\begin{split}\min_x \ & c^T x \\ \mbox{such that} \ & A_{ub} x \leq b_{ub},\\ & A_{eq} x = b_{eq},\\ & l \leq x \leq u ,\end{split}\]其中 \(x\) 是决策变量的向量;\(c\)、\(b_{ub}\)、\(b_{eq}\)、\(l\) 和 \(u\) 是向量;\(A_{ub}\) 和 \(A_{eq}\) 是矩阵。
或者,
最小化
c @ x
使得
A_ub @ x <= b_ub A_eq @ x == b_eq lb <= x <= ub
请注意,默认情况下
lb = 0
和ub = None
。可以通过bounds
指定其他边界。- 参数::
- c一维数组
要最小化的线性目标函数的系数。
- A_ub二维数组,可选
不等式约束矩阵。
A_ub
的每一行指定了对x
的线性不等式约束的系数。- b_ub一维数组,可选
不等式约束向量。每个元素代表对
A_ub @ x
对应值的约束上限。- A_eq二维数组,可选
等式约束矩阵。
A_eq
的每一行指定了对x
的线性等式约束的系数。- b_eq一维数组,可选
等式约束向量。
A_eq @ x
的每个元素必须等于b_eq
的对应元素。- bounds序列,可选
对
x
中每个元素的(min, max)
对的序列,定义该决策变量的最小值和最大值。如果提供一个单一的元组(min, max)
,那么min
和max
将作为所有决策变量的边界。使用None
表示没有边界。例如,默认边界(0, None)
表示所有决策变量是非负的,而对(None, None)
表示根本没有边界,即所有变量都允许为任何实数。- method字符串,可选
用于解决标准形式问题的算法。 ‘highs’ (默认)、‘highs-ds’、‘highs-ipm’、‘interior-point’ (遗留)、‘revised simplex’ (遗留) 和 ‘simplex’ (遗留) 都受支持。遗留方法已弃用,将在 SciPy 1.11.0 中删除。
- callback可调用对象,可选
如果提供了回调函数,它将至少在算法的每次迭代中被调用一次。回调函数必须接受单个
scipy.optimize.OptimizeResult
,其中包含以下字段- x一维数组
当前解向量。
- fun浮点数
目标函数
c @ x
的当前值。- success布尔值
当算法成功完成时为
True
。- slack一维数组
松弛的(名义上为正)值,
b_ub - A_ub @ x
。- con一维数组
等式约束的(名义上为零)残差,
b_eq - A_eq @ x
。- phase整数
正在执行的算法阶段。
- status整数
表示算法状态的整数。
0
: 优化过程正常进行。1
: 达到迭代次数限制。2
: 问题似乎不可行。3
: 问题似乎无界。4
: 遇到数值困难。- nit整数
当前迭代次数。
- message字符串
算法状态的字符串描述。
回调函数目前不受 HiGHS 方法支持。
- options字典,可选
求解器选项的字典。所有方法都接受以下选项
- maxiter整数
要执行的最大迭代次数。默认值:参见特定方法的文档。
- disp布尔值
设置为
True
以打印收敛消息。默认值:False
。- presolve布尔值
设置为
False
以禁用自动预处理。默认值:True
。
除了 HiGHS 求解器之外,所有方法还接受
- tol浮点数
确定何时将残差视为“足够接近”零以被认为完全为零的容差。
- autoscale布尔值
设置为
True
以自动执行均衡。如果约束中的数值相差几个数量级,请考虑使用此选项。默认值:False
。- rr布尔值
设置为
False
以禁用自动冗余消除。默认值:True
。- rr_method字符串
用于在预处理后识别并从等式约束矩阵中消除冗余行的算法。对于具有密集输入的问题,可用的冗余消除方法有
- “SVD”
对矩阵重复执行奇异值分解,根据对应于零奇异值的左奇异向量中的非零元素来检测冗余行。当矩阵接近满秩时,可能会很快。
- “pivot”
使用 [5] 中介绍的算法来识别冗余行。
- “ID”
使用随机插值分解。识别矩阵转置中未用于矩阵的满秩插值分解的列。
- None
如果矩阵接近满秩,即矩阵秩与行数之差小于 5,则使用“svd”。否则,使用“pivot”。此默认行为可能在未经事先通知的情况下发生变化。
默认值:None。对于具有稀疏输入的问题,此选项将被忽略,并将使用 [5] 中介绍的基于枢轴的算法。
有关特定方法的选项,请参见
show_options('linprog')
.- x0一维数组,可选
决策变量的猜测值,这些值将由优化算法进行细化。此参数目前仅由 “revised simplex” 方法使用,并且只能在 x0 表示一个基本可行解时使用。
- integrality一维数组或整数,可选
指示对每个决策变量的整数约束类型。
0
: 连续变量;没有整数约束。1
: 整数变量;决策变量必须是在 bounds 内的整数。2
: 半连续变量;决策变量必须在bounds内或取值0
。3
: 半整数变量;决策变量必须是bounds内的整数或取值0
。默认情况下,所有变量都是连续的。
对于混合整数约束,提供一个形状为c.shape的数组。要从较短的输入推断每个决策变量的约束,该参数将使用np.broadcast_to广播到c.shape。
此参数目前仅由
'highs'
方法使用,其他情况下会被忽略。
- 返回值:
- resOptimizeResult
一个
scipy.optimize.OptimizeResult
,包含以下字段。请注意,字段的返回类型可能取决于优化是否成功,因此建议在依赖其他字段之前检查OptimizeResult.status- x一维数组
决策变量的值,这些值在满足约束条件的同时使目标函数最小化。
- fun浮点数
目标函数
c @ x
的最优值。- slack一维数组
松弛变量的(名义上为正)值,
b_ub - A_ub @ x
。- con一维数组
等式约束的(名义上为零)残差,
b_eq - A_eq @ x
。- success布尔值
True
当算法成功找到最优解时。- status整数
一个整数,表示算法的退出状态。
0
: 优化成功终止。1
: 达到迭代次数限制。2
: 问题似乎不可行。3
: 问题似乎无界。4
: 遇到数值困难。- nit整数
所有阶段执行的迭代总数。
- message字符串
算法退出状态的字符串描述。
参见
show_options
求解器接受的附加选项。
备注
本节描述了可通过“method”参数选择的可用求解器。
‘highs-ds’ 和 ‘highs-ipm’ 是对 HiGHS 单纯形和内点法求解器 [13] 的接口,分别。‘highs’(默认)会在两者之间自动选择。这些是 SciPy 中最快的线性规划求解器,特别是对于大型稀疏问题;哪一个更快取决于问题。其他求解器(‘interior-point’、‘revised simplex’ 和 ‘simplex’)是遗留方法,将在 SciPy 1.11.0 中删除。
方法highs-ds 是 C++ 高性能对偶修正单纯形实现 (HSOL) [13]、[14] 的包装器。方法highs-ipm 是对 C++ 实现的 **i**nterior-**p**oint **m**ethod [13] 的包装器;它具有交叉例程,因此与单纯形求解器一样精确。方法highs 会在两者之间自动选择。对于涉及
linprog
的新代码,建议明确选择这三个方法值之一。在版本 1.6.0 中添加。
方法interior-point 使用如 [4] 中所述的原始对偶路径跟踪算法。该算法支持稀疏约束矩阵,通常比单纯形方法更快,尤其是在大型稀疏问题上。但是请注意,返回的解决方案可能比单纯形方法的解决方案略微不准确,并且通常不对应于由约束定义的多面体的顶点。
在版本 1.0.0 中添加。
方法revised simplex 使用如 [9] 中所述的修正单纯形方法,除了它有效地维护基矩阵的分解 [11] 而不是它的逆,并将其用于解决算法每次迭代的线性系统。
在版本 1.3.0 中添加。
方法simplex 使用 Dantzig 单纯形算法 [1]、[2] 的传统完整表实现(不是 Nelder-Mead 单纯形)。该算法是为了向后兼容和教育目的而包含的。
在版本 0.15.0 中添加。
在应用interior-point、revised simplex 或simplex 之前,一个基于 [8] 的预处理过程试图识别微不足道的不可行性、微不足道的无界性和潜在的问题简化。具体来说,它检查
A_eq
或A_ub
中的零行,表示微不足道的约束;A_eq
和A_ub
中的零列,表示无约束变量;A_eq
中的单列,表示固定变量;以及A_ub
中的单列,表示简单边界。
如果预处理表明问题是无界的(例如,无约束且无界变量具有负成本)或不可行的(例如,
A_eq
中的零行对应于b_eq
中的非零),则求解器会终止并显示相应的状态代码。请注意,预处理会尽快终止,只要检测到任何无界性的迹象;因此,问题可能会报告为无界,而实际上问题是不可行的(但不可行性尚未检测到)。因此,如果必须知道问题是否实际上不可行,请再次使用选项presolve=False
来解决问题。如果在预处理的一次扫描中未检测到不可行性或无界性,则尽可能收紧边界并从问题中删除固定变量。然后,
A_eq
矩阵的线性相关行将被删除(除非它们表示不可行性),以避免主要求解例程中的数值困难。请注意,几乎线性相关的行(在规定的容差范围内)也可能被删除,这在极少数情况下会改变最优解。如果这是一个问题,请从问题公式中消除冗余,并使用选项rr=False
或presolve=False
运行。这里可以进行一些潜在的改进:[8] 中概述的额外预处理检查应该被实现,预处理例程应该运行多次(直到无法进行进一步的简化),并且 [5] 中的更多效率改进应该在冗余删除例程中实现。
经过预处理后,问题将通过将(收紧的)简单边界转换为上界约束、为不等式约束引入非负松弛变量以及将无界变量表示为两个非负变量之差来转换为标准形式。可选地,问题可以通过平衡 [12] 自动缩放。选定的算法解决标准形式问题,后处理例程将结果转换为原始问题的解决方案。
参考文献
[1]Dantzig, George B., 线性规划及其扩展。Rand 公司研究报告 普林斯顿大学出版社,新泽西州普林斯顿,1963 年
[2]Hillier, S.H. and Lieberman, G.J. (1995), “数学规划导论”, McGraw-Hill, 第 4 章。
[3]Bland, Robert G. 单纯形方法的新有限枢轴规则。运筹学研究 (2), 1977: pp. 103-107。
[4]Andersen, Erling D., and Knud D. Andersen. “用于线性规划的 MOSEK 内点优化器:同类算法的实现。” 高性能优化。Springer US, 2000. 197-232。
[6]Freund, Robert M. “基于牛顿法的线性规划的原始对偶内点法。” 未出版课程笔记,2004 年 3 月。于 2017 年 2 月 25 日提供 https://ocw.mit.edu/courses/sloan-school-of-management/15-084j-nonlinear-programming-spring-2004/lecture-notes/lec14_int_pt_mthd.pdf
[7]Fourer, Robert. “通过内点法求解线性规划。” 未出版课程笔记,2005 年 8 月 26 日。于 2017 年 2 月 25 日提供 http://www.4er.org/CourseNotes/Book%20B/B-III.pdf
[9]Bertsimas, Dimitris, and J. Tsitsiklis. “线性规划导论。” Athena Scientific 1 (1997): 997。
[10]Andersen, Erling D., 等人。大型线性规划内点法实现。HEC/日内瓦大学,1996 年。
[11]Bartels, Richard H. “单纯形方法的稳定化。” Numerische Mathematik 杂志 16.5 (1971): 414-434。
[12]Tomlin, J. A. “关于缩放线性规划问题。” 数学规划研究 4 (1975): 146-166。
[13] (1,2,3)Huangfu, Q., Galabova, I., Feldmeier, M., and Hall, J. A. J. “HiGHS - 用于线性优化的高性能软件。” https://highs.dev/
[14]Huangfu, Q. and Hall, J. A. J. “对偶修正单纯形法的并行化。” 数学规划计算, 10 (1), 119-142, 2018. DOI: 10.1007/s12532-017-0130-5
示例
考虑以下问题
\[\begin{split}\min_{x_0, x_1} \ -x_0 + 4x_1 & \\ \mbox{such that} \ -3x_0 + x_1 & \leq 6,\\ -x_0 - 2x_1 & \geq -4,\\ x_1 & \geq -3.\end{split}\]该问题并未以
linprog
接受的形式呈现。这很容易通过将“大于”不等式约束转换为“小于”不等式约束来解决,方法是将两边乘以因子 \(-1\)。还请注意,最后一个约束实际上是简单边界 \(-3 \leq x_1 \leq \infty\)。最后,由于 \(x_0\) 没有边界,我们必须明确指定边界 \(-\infty \leq x_0 \leq \infty\),因为默认情况下变量是非负的。在将系数收集到数组和元组中之后,此问题的输入为>>> from scipy.optimize import linprog >>> c = [-1, 4] >>> A = [[-3, 1], [1, 2]] >>> b = [6, 4] >>> x0_bounds = (None, None) >>> x1_bounds = (-3, None) >>> res = linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]) >>> res.fun -22.0 >>> res.x array([10., -3.]) >>> res.message 'Optimization terminated successfully. (HiGHS Status 7: Optimal)'
边际(也称为对偶值/影子价格/拉格朗日乘子)和残差(松弛)也可以使用。
>>> res.ineqlin residual: [ 3.900e+01 0.000e+00] marginals: [-0.000e+00 -1.000e+00]
例如,因为与第二个不等式约束相关的边际为 -1,所以我们预计如果在第二个不等式约束的右侧添加一个小的量
eps
,目标函数的最优值将减少eps
>>> eps = 0.05 >>> b[1] += eps >>> linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]).fun -22.05
此外,因为第一个不等式约束的残差为 39,所以我们可以在不影响最优解的情况下将第一个约束的右侧减少 39。
>>> b = [6, 4] # reset to original values >>> b[0] -= 39 >>> linprog(c, A_ub=A, b_ub=b, bounds=[x0_bounds, x1_bounds]).fun -22.0