scipy.optimize.

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)[源码]#

线性规划:在满足线性和非线性约束的条件下,最小化线性目标函数。

线性规划解决以下形式的问题:

\[\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 = 0ub = None。其他边界可以通过 bounds 指定。

参数:
c1-D 数组

要最小化的线性目标函数的系数。

A_ub2-D 数组,可选

不等式约束矩阵。 A_ub 的每一行指定了 x 上的线性不等式约束的系数。

b_ub1-D 数组,可选

不等式约束向量。每个元素代表 A_ub @ x 对应值的上限。

A_eq2-D 数组,可选

等式约束矩阵。 A_eq 的每一行指定了 x 上的线性等式约束的系数。

b_eq1-D 数组,可选

等式约束向量。 A_eq @ x 的每个元素必须等于 b_eq 的对应元素。

bounds序列,可选

一个序列,包含 (min, max) 中每个元素的 x 对,定义了该决策变量的最小值和最大值。如果提供单个元组 (min, max),则 minmax 将作为所有决策变量的边界。使用 None 表示没有边界。例如,默认边界 (0, None) 意味着所有决策变量都是非负的,而 (None, None) 表示完全没有边界,即所有变量都可以是任意实数。

method字符串,可选

用于解决标准形式问题的算法。支持以下算法:

旧版方法已弃用,将在 SciPy 1.11.0 中移除。

callback可调用对象,可选

如果提供了回调函数,它将在算法每次迭代时至少被调用一次。回调函数必须接受一个包含以下字段的 scipy.optimize.OptimizeResult 对象:

x1-D 数组

当前解向量。

fun浮点数

目标函数 c @ x 的当前值。

success布尔值

当算法成功完成时为 True

slack1-D 数组

松弛变量(名义上为正)的值,即 b_ub - A_ub @ x

con1-D 数组

等式约束的(名义上为零的)残差,即 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:

如果矩阵接近满秩(即矩阵秩与行数之差小于五),则使用 svd。否则,使用 pivot。此默认行为可能会在不事先通知的情况下发生变化。

默认值:None。对于稀疏输入问题,此选项将被忽略,并使用 [5] 中介绍的基于轴的算法。

有关特定方法的选项,请参阅 show_options('linprog')

x01-D 数组,可选

决策变量的猜测值,将由优化算法进行细化。此参数目前仅用于 ‘revised simplex’ 方法,并且仅当 x0 代表一个基本可行解时才能使用。

integrality1-D 数组或整数,可选

指示每个决策变量的整数约束类型。

0:连续变量;无整数约束。

1:整数变量;决策变量必须在 bounds 内为整数。

2:半连续变量;决策变量必须在 bounds 内或取值 0

3:半整数变量;决策变量必须在 bounds 内为整数或取值 0

默认情况下,所有变量都是连续的。

对于混合整数约束,提供一个形状为 c.shape 的数组。为了从较短的输入中推断出每个决策变量的约束,参数将使用 numpy.broadcast_to 广播到 c.shape

此参数目前仅由 ‘highs’ 方法使用,否则将被忽略。

返回:
resOptimizeResult

一个 scipy.optimize.OptimizeResult 对象,包含以下字段。请注意,字段的返回类型可能取决于优化是否成功,因此建议在依赖其他字段之前检查 OptimizeResult.status

x1-D 数组

最小化目标函数同时满足约束条件的决策变量值。

fun浮点数

目标函数 c @ x 的最优值。

slack1-D 数组

松弛变量(名义上为正)的值,即 b_ub - A_ub @ x

con1-D 数组

等式约束的(名义上为零的)残差,即 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 中最快的线性规划求解器,尤其适用于大型稀疏问题;这两种方法中哪一个更快取决于具体问题。callback 参数被 HiGHS 方法支持后,其他求解器将被移除,因为它们是旧版方法。

方法 ‘highs-ds’ 是 C++ 高性能双修订单纯形法实现 (HSOL) [13], [14] 的封装。方法 ‘highs-ipm’ 是 C++ 内点法实现 [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_eqA_ub 中的零行,表示简单约束;

  • A_eq A_ub 中的零列,表示无约束变量;

  • A_eq 中的单列,表示固定变量;以及

  • A_ub 中的单列,表示简单边界。

如果预处理显示问题是无界的(例如,无约束且无界变量的成本为负)或不可行的(例如,A_eq 中的零行对应于 b_eq 中的非零值),则求解器将以适当的状态码终止。请注意,预处理会在检测到任何无界迹象时立即终止;因此,一个问题可能被报告为无界,而实际上该问题是不可行的(但尚未检测到不可行性)。因此,如果了解问题是否确实不可行很重要,请使用选项 presolve=False 再次求解问题。

如果在预处理的单次通过中既未检测到不可行性也未检测到无界性,则会在可能的情况下收紧边界,并从问题中移除固定变量。然后,移除 A_eq 矩阵中线性相关的行(除非它们表示不可行性),以避免主求解例程中的数值困难。请注意,接近线性相关的行(在规定容差范围内)也可能被移除,这在极少数情况下可能会改变最优解。如果这引起担忧,请从问题公式中消除冗余,并使用选项 rr=Falsepresolve=False 运行。

这里可以进行一些潜在的改进:应实现 [8] 中概述的额外预处理检查,预处理例程应多次运行(直到无法再进行简化),并且应在冗余删除例程中实现 [5] 中更多的效率改进。

预处理后,问题通过将(收紧的)简单边界转换为上限约束、为不等式约束引入非负松弛变量以及将无界变量表示为两个非负变量之间的差来转换为标准形式。此外,问题可选择通过均衡自动缩放 [12]。所选算法求解标准形式问题,后处理例程将结果转换为原始问题的解。

参考文献

[1]

Dantzig, George B., Linear programming and extensions. Rand Corporation Research Study Princeton Univ. Press, Princeton, NJ, 1963

[2]

Hillier, S.H. and Lieberman, G.J. (1995), “Introduction to Mathematical Programming”, McGraw-Hill, Chapter 4.

[3]

Bland, Robert G. New finite pivoting rules for the simplex method. Mathematics of Operations Research (2), 1977: pp. 103-107.

[4]

Andersen, Erling D., and Knud D. Andersen. “The MOSEK interior point optimizer for linear programming: an implementation of the homogeneous algorithm.” High performance optimization. Springer US, 2000. 197-232.

[5] (1,2,3)

Andersen, Erling D. “Finding all linearly dependent rows in large-scale linear programming.” Optimization Methods and Software 6.3 (1995): 219-227.

[6]

Freund, Robert M. “Primal-Dual Interior-Point Methods for Linear Programming based on Newton’s Method.” Unpublished Course Notes, March 2004. Available 2/25/2017 at 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. “Solving Linear Programs by Interior-Point Methods.” Unpublished Course Notes, August 26, 2005. Available 2/25/2017 at http://www.4er.org/CourseNotes/Book%20B/B-III.pdf

[8] (1,2)

Andersen, Erling D., and Knud D. Andersen. “Presolving in linear programming.” Mathematical Programming 71.2 (1995): 221-245.

[9]

Bertsimas, Dimitris, and J. Tsitsiklis. “Introduction to linear programming.” Athena Scientific 1 (1997): 997.

[10]

Andersen, Erling D., et al. Implementation of interior point methods for large scale linear programming. HEC/Universite de Geneve, 1996.

[11]

Bartels, Richard H. “A stabilization of the simplex method.” Journal in Numerische Mathematik 16.5 (1971): 414-434.

[12]

Tomlin, J. A. “On scaling linear programming problems.” Mathematical Programming Study 4 (1975): 146-166.

[13] (1,2,3)

Huangfu, Q., Galabova, I., Feldmeier, M., and Hall, J. A. J. “HiGHS - high performance software for linear optimization.” https://highs.dev/

[14]

Huangfu, Q. and Hall, J. A. J. “Parallelizing the dual revised simplex method.” Mathematical Programming Computation, 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