混合整数线性规划 (MILP) 算法

混合整数线性规划定义

混合整数线性规划 (MILP) 问题具有以下要素：

线性目标函数 f^Tx，其中 f 是由常数组成的列向量，x 是由未知数组成的列向量
边界和线性约束，但没有非线性约束（有关定义，请参阅编写约束）
对 x 的某些分量的限制，使其必须具有整数值

以数学语言表达，即根据向量 f、lb 和 ub，矩阵 A 和 Aeq，对应的向量 b 和 beq，以及索引集 intcon，求解向量 x 使下式成立

$\min_{x} f^{T} x subject to {\begin{array}{l} x (intcon) are integers \\ A \cdot x \leq b \\ A e q \cdot x = b e q \\ l b \leq x \leq u b . \end{array}$

算法概述

intlinprog 使用此基本策略来求解混合整数线性规划。intlinprog 可以在任一阶段完成问题的求解。如果它在某个阶段成功求解了问题，intlinprog 不会执行后面的阶段。

使用线性规划预处理缩减问题的规模。
使用线性规划求解初始松弛（非整数）问题。
执行混合整数规划预处理以收紧混合整数问题的 LP 松弛。
尝试切割生成以进一步收紧混合整数问题的 LP 松弛。
尝试使用启发式方法求得整数可行解。
使用分支定界算法系统地搜索最优解。此算法通过限制整数变量的可能值范围来求解 LP 松弛问题。它尝试在最优目标函数值上生成一系列更新边界。

线性规划预处理

根据混合整数线性规划定义，矩阵 A 和 Aeq 以及对应的向量 b 和 beq 以如下形式编写一组线性不等式和线性等式

$\begin{matrix} A \cdot x \leq b \\ A e q \cdot x = b e q . \end{matrix}$

这些线性约束对解 x 施加约束。

通常，我们可以减少问题中的变量数量（x 的分量数量），并减少线性约束的数量。虽然执行这些缩减可能会花费求解器的时间，但它们通常会缩短求解的总时间，并使更大型的问题可求解。这些算法可以使解在数值上更加稳定。此外，这些算法有时可以检测不可行问题。

预处理步骤旨在消除冗余变量和约束，改善模型的尺度和约束矩阵的稀疏性，加强变量的边界，检测模型的原始和对偶不可行性。

有关详细信息，请参阅 Andersen 和 Andersen 的论著 [2] 以及 Mészáros 和 Suhl 的论著 [8]。

线性规划

初始松弛问题是线性规划问题，其目标和约束与混合整数线性规划定义相同，但没有整数约束。x_LP 称为松弛问题的解，x 称为具有整数约束的原始问题的解。显然，

f^Tx_LP ≤ f^Tx，

因为 x_LP 在较少的约束下最小化相同的函数。

此初始松弛 LP（根节点 LP）以及分支定界算法期间生成的所有 LP 松弛都使用线性规划方法来求解。

混合整数规划预处理

在混合整数规划预处理期间，intlinprog 分析线性不等式 A*x ≤ b 以及完整性约束，以确定是否有以下情形：

此问题不可行。
有些边界可以收紧。
有些不等式是冗余的，因此可以忽略或删除。
一些不等式可以得到加强。
一些整数变量可以固定。

通过 IntegerPreprocess 选项，您可以选择 intlinprog 采取某些步骤、采取所有步骤还是几乎不采取任何步骤。如果包含 x0 参数，intlinprog 会在预处理中使用该值。

混合整数规划预处理的主要目标是简化后续的分支定界计算。预处理包括快速预检查和消除一些无用的子问题候选项，以免分支定界算法对其进行分析。

有关整数预处理的详细信息，请参阅 Savelsbergh 的 [10]。

切割生成

切割是 intlinprog 为问题增加的额外线性不等式约束。这些不等式尝试限制 LP 松弛的可行域，使其解更接近整数。您可以通过 CutGeneration 选项控制 intlinprog 使用的切割类型。

'basic' 切割包括：

混合整数舍入切割
Gomory 切割
Clique 切割
覆盖切割
流覆盖切割

此外，如果问题为纯整数型（所有变量均为整数值），则 intlinprog 还使用以下切割：

强 Chvatal-Gomory 切割
零-半切割

'intermediate' 切割包括所有 'basic' 切割，还包括：

简单的提升投影切割
简单的主元缩减切割
缩减拆分切割

'advanced' 切割包括缩减拆分切割以外的所有 'intermediate' 切割，还包括：

强 Chvatal-Gomory 切割
零-半切割

对于纯整数型问题，'intermediate' 使用的切割类型最多，因为它使用缩减拆分切割，而 'advanced' 不使用。

另一个选项 CutMaxIterations 指定 intlinprog 迭代切割生成的次数上界。

有关切割生成算法（也称为割平面方法）的详细信息，请参阅 Cornuéjols·的论著 [5]；有关 clique 切割的详细信息，请参阅 Atamtürk、Nemhauser 和 Savelsbergh·的论著 [3]。

使用启发式方法求出可行解

为了得到目标函数的上界，分支定界过程必须求得可行点。分支定界过程中，LP 松弛问题的解可能是整数可行的，可用于改进原始 MILP 的上界。某些方法能够在分支定界之前或期间更快地找到可行点。intlinprog 在根节点上以及一些分支定界迭代过程中使用这些方法。这些方法是启发式的，意味着它们是可能成功但也可能失败的算法。

启发式方法可以是初始点启发式方法，它们帮助求解器求得初始或新的整数可行解。启发式方法也可以是改进点启发式方法，它们从整数可行点开始，尝试找到更好的整数可行点，即目标函数值更低的点。intlinprog 改进点启发式方法有 'rins'、'rss'、1-opt、2-opt 和引导潜水。

使用 'Heuristics' 选项设置 intlinprog 使用的启发式方法。选项包括：

选项	说明
`'basic'`（默认值）	求解器使用不同参数运行两次舍入启发式方法，再使用不同参数运行两次潜水启发式方法，然后运行 `'rss'`。如果前期的启发式方法能得到足够好的整数可行解，则求解器不再运行其后的启发式方法。
`'intermediate'`	求解器使用不同参数运行两次舍入启发式方法，然后使用不同参数运行两次潜水启发式方法。如果有整数可行解，则求解器运行 `'rins'`，然后运行 `'rss'`。如果 `'rss'` 求得新解，求解器将再次运行 `'rins'`。如果前期的启发式方法能得到足够好的整数可行解，则求解器不再运行其后的启发式方法。
`'advanced'`	求解器使用不同参数运行两次舍入启发式方法，然后使用不同参数运行两次潜水启发式方法。如果有整数可行解，则求解器运行 `'rins'`，然后运行 `'rss'`。如果 `'rss'` 求得新解，求解器将再次运行 `'rins'`。如果前期的启发式方法能得到足够好的整数可行解，则求解器不再运行其后的启发式方法。
`'rins'` 或等效的 `'rins-diving'`	`intlinprog` 搜索当前最佳整数可行解点（如有）的邻域，以求得更好的新解。请参阅 Danna、Rothberg 和 Le Pape·的论著 [6]。当您选择 `'rins'` 时，求解器使用不同参数运行两次舍入启发式方法，再使用不同参数运行两次潜水启发式方法，然后运行 `'rins'`。
`'rss'` 或等效的 `'rss-diving'`	`intlinprog` 将 `'rins'` 理念和局部分支理念整合，以混合方式搜索整数可行解。当您选择 `'rss'` 时，求解器使用不同参数运行两次舍入启发式方法，再使用不同参数运行两次潜水启发式方法，然后运行 `'rss'`。如果前期的启发式方法能得到足够好的整数可行解，则求解器不再运行其后的启发式方法。这些设置执行与 `'basic'` 相同的启发式方法。
`'round'`	`intlinprog` 取节点处松弛问题的 LP 解，并以尝试保持可行性的方式舍入整数分量。当您选择 `'round'` 时，求解器在根节点上使用不同参数运行两次舍入启发式方法，然后使用不同参数运行两次潜水启发式方法。此后，求解器仅在一些分支定界节点上运行舍入启发式方法。
`'round-diving'`	求解器的工作方式与 `'round'` 相似，但求解器还在一些分支定界节点上运行潜水启发式方法（在舍入启发式方法的基础上），而不仅限于根节点。
`'diving'`	`intlinprog` 使用类似于分支定界步骤的启发式方法，但只沿树的一个分支向下，而不创建其他分支。此单分支使得算法能够快速下“潜”到树的片段，因此命名为“潜水”。当前，`intlinprog` 按照以下顺序使用六种潜水启发式方法：向量长度潜水系数潜水小数潜水伪代价潜水线搜索潜水引导式潜水（当求解器已找到至少一个整数可行点时适用）潜水启发式方法通常选择一个应该取整数值、但当前解为小数的变量。然后，启发式方法引入一个强制变量取整数值的边界，并再次求解相关联的松弛 LP。各种潜水启发式方法之间的主要区别在于如何选择要定界的变量。请参阅 Berthold 的论著·[4]，第 3.1 节。
`'none'`	`intlinprog` 不搜索可行点。求解器只选取它在分支定界搜索中遇到的任何可行点。

'intermediate' 和 'advanced' 之间的主要区别在于 'advanced' 在分支定界迭代期间更频繁地运行启发式方法。

除了上表之外，当 Heuristics 选项为 'basic'、'intermediate' 或 'advanced' 时，还会运行以下启发式方法。

ZI 取整 - 每当算法解出一个松弛 LP 时，此启发式方法就会运行。该启发式方法遍历每个带小数的整数变量并尝试将其移至邻近的整数，不影响在其他约束下的可行性。有关详细信息，请参阅 Hendel [7]。
1-opt - 每当算法找到新的整数可行解时，此启发式方法就会运行。该启发式方法在降低目标函数值的同时遍历每个整数变量并尝试将其移至邻近的整数，不影响在其他约束下的可行性。
2-opt - 每当算法找到新的整数可行解时，此启发式方法就会运行。2-opt 查找影响同一约束的所有整数变量对组，这意味着它们在 A 或 Aeq 约束矩阵的同一行中有非零项。对于每个对组，2-opt 选取一个整数可行解，并使用所有四个可能的移动（上-上、上-下、下-上和下-下）向上或向下移动变量对组的值，寻找具有更好的目标函数值的可行相邻解。该算法通过计算满足约束的对组中每个变量的最大移位大小（相同幅值）来测试每个整数变量对组，同时改进目标函数值。

在启发式方法阶段开始时，intlinprog 运行平凡启发式方法，除非 Heuristics 是 'none' 或者您在 x0 参数中提供初始整数可行点。平凡启发式方法检查以下各点的可行性：

全部为零
上界
下界（如果非零）
“锁定”点

“锁定”点仅针对所有变量的上界和下界有限的问题定义。每个变量的“锁定”点是其上界或下界，按如下方式进行选择。对于每个变量 j，计算线性约束矩阵 A(:,j) 中对应正条目的数量，并减去对应负条目的数量。如果结果为正值，则使用该变量的下界 lb(j)。否则，使用该变量的上界 ub(j)。“锁定”点尝试满足每个变量的最大数量的线性不等式约束，但不一定可行。

在每个启发式方法都获得一个可行解后，intlinprog 调用输出函数和绘图函数。请参阅intlinprog Output Function and Plot Function Syntax。

如果包含 x0 参数，intlinprog 将在 'rins' 中和引导潜水启发式方法中使用该值，直到找到更好的整数可行点。因此，当您提供 x0 时，您可以通过将 'Heuristics' 选项设置为 'rins-diving' 或使用 'rins' 的其他设置来获得满意的结果。

分支定界

分支定界方法构造一系列尝试收敛到 MILP 解的子问题。子问题给出关于解 f^Tx 的一系列上界和下界。第一个上界是任何可行解，第一个下界是松弛问题的解。关于上界的讨论，请参阅使用启发式方法求出可行解。

如线性规划中所述，线性规划松弛问题的任何解都比 MILP 解具有更低的目标函数值。此外，任何可行点 x_feas 都满足

f^Tx_feas ≥ f^Tx，

因为 f^Tx 是所有可行点中最小的。

在这种情况下，节点是具有以下特征的 LP：它具有与原始问题相同的目标函数、边界和线性约束，但没有整数约束，并对线性约束或边界作出特定改变。根节点是没有整数约束且线性约束或边界也不发生变化的原始问题，这意味着根节点是初始松弛 LP。

从起始边界开始，分支定界方法通过从根节点产生分支来构造新子问题。分支步骤采用启发式方法，遵循以下规则之一。每个规则都基于这样的理念：约束一个变量，使之小于或等于整数 J，或者大于或等于 J+1，从而拆分问题。当 x_LP 中对应于 intcon 所指定整数的项不是整数时，会出现这两个子问题。在此处，x_LP 是松弛问题的解。以 J 为变量的下限（向下舍入），以 J+1 为上限（向上舍入）。由此产生的两个问题的解大于或等于 f^Tx_LP，因为这两个问题具有更多约束。因此，此过程可能会提高下界。

分支定界方法的性能取决于基于何种规则选择要拆分的变量（分支规则）。您可以在 BranchRule 选项中设置算法所用的规则，这些规则如下：

'maxpscost' - 选择具有最大伪代价的小数变量。
伪代价
变量 i 的伪代价基于以下两个量：当 i 被选为分支变量时下界变化量的经验估计值，以及当前点 x 的分量 i 的小数部分。小数部分 p 分为两块，即下部和上部：
p_i^– = x(i) – ⌊x(i)⌋
p_i⁺ = 1 – p_i^–。
将 x_i^– 设置为受 x(i) ≤ ⌊x(i)⌋ 约束的线性规划的解，并让目标函数中的变化表示为
Δ_i^– = f^Tx_i^– – f^Tx。
类似地，Δ_i⁺ 是当问题受 x(i) ≥ ⌈x(i)⌉ 约束时目标函数的变化量。
变量 x_i 的每单位变化引起的目标函数增量为
$d_{i}^{-} = \frac{Δ_{i}^{-}}{p_{i}^{-}} or d_{i}^{+} = \frac{Δ_{i}^{+}}{p_{i}^{+}} .$
将分支定界算法过程中到该点为止的 s_i^– 和 s_i⁺ 的经验平均值设为 d_i^– 和 d_i⁺。在获得任何观测值之前，项的经验值将初始化为目标系数 f(i) 的绝对值。然后，'maxpscost' 规则对符合以下条件的节点 i 进行分支：该节点可基于正权重 w⁺ 和 w^– 最大化以下量
w^– * p_i^– * s_i^– + w⁺ * p_i⁺ * s_i⁺。
大体而言，此规则将选择有可能最大程度地增加下界的系数。
'strongpscost' - 类似于 'maxpscost'，但对于每个变量，求解器不会将伪代价初始化为 1，而是仅在伪代价具有更可靠的估计值后，才尝试对变量进行分支。为了获得更可靠的估计值，求解器执行以下操作（请参阅 Achterberg、Koch 和 Martin 的论著·[1]）。
- 对所有潜在分支变量（当前为小数但应为整数的变量）排序，依据是这些变量当前基于伪代价的得分。
- 从得分最高的变量开始，基于当前分支变量运行两个松弛线性规划（前提是该变量尚未用于分支计算）。求解器使用这两个解来更新当前分支变量的伪代价。求解器可以提前停止此过程，以节省选择分支的时间。
- 继续在列表中选择变量，直到当前基于伪代价的最高得分经连续 k 个变量保持不变，其中 k 是内部选择的值，通常在 5 到 10 之间。
- 对基于伪代价得分最高的变量进行分支。求解器可能已在先前的伪代价估算过程中基于此变量计算了松弛线性规划。
由于额外的线性规划解，'strongpscost' 分支的每次迭代都比默认的 'maxpscost' 花费更长时间。但是，分支定界迭代的数量通常会减少，因此 'strongpscost' 方法总体上可以节省时间。
'reliability' - 类似于 'strongpscost'，但 'reliability' 不是仅对未初始化的伪代价分支运行松弛线性规划，而是对每个变量运行多达 k2 次规划，其中 k2 是小整数，如 4 或 8。因此，与 'strongpscost' 相比，'reliability' 的分支速度更慢，但分支定界迭代可能更少。
'mostfractional' - 选择小数部分最接近 1/2 的变量。
'maxfun' - 选择目标向量 f 中最大的绝对值所对应的变量。

在算法分支后，有两个新节点需要探查。算法使用以下规则之一在所有可用节点中选择要探查的节点：

'minobj' - 选择目标函数值最低的节点。
'mininfeas' - 选择整数不可行性之和最小的节点。这意味着对于节点中的每个整数不可行分量 x(i)，选择 p_i^– 和 p_i⁺ 中较小的一个进行累加，其中
p_i^– = x(i) – ⌊x(i)⌋
p_i⁺ = 1 – p_i^–。
'simplebestproj' - 选择具有最佳投影的节点。
最佳投影
将目前为止找到的最佳整数可行点记作 x_B，根节点处的 LP 松弛解记作 x_R，要检查的节点记作 x。将节点 x 处整数不可行性的总和记作 in(x)（请参阅 'mininfeas'）。最佳投影规则的目标是最小化
$f^{T} x + \frac{f^{T} x_{B} - f^{T} x_{R}}{i n (x_{R})} i n (x) .$
如果到当前为止未找到整数可行点，请设置 f^Tx_B = 0。

intlinprog 将原始问题的信息（例如目标函数的最大公约数，即 GCD）纳入考虑，从而跳过一些子问题的分析。

分支定界过程继续按部就班地生成子问题以进行分析，并丢弃那些不会改进目标函数上界或下界的子问题，直到满足以下停止条件之一：

算法超过 MaxTime 选项指定的运行时间。
目标函数的上下界之差小于 AbsoluteGapTolerance 或 RelativeGapTolerance 容差。
探查的节点数量超过 MaxNodes 选项指定的数量。
整数可行点数超过 MaxFeasiblePoints 选项指定的数量。

有关分支定界过程的详细信息，请参阅 Nemhauser 和 Wolsey 的 [9] 和 Wolsey 的 [11]。

参考

[1] Achterberg, T., T. Koch and A. Martin. Branching rules revisited. Operations Research Letters 33, 2005, pp. 42–54. Available at https://opus4.kobv.de/opus4-zib/files/788/ZR-04-13.pdf.

[2] Andersen, E. D., and Andersen, K. D. Presolving in linear programming. Mathematical Programming 71, pp. 221–245, 1995.

[3] Atamtürk, A., G. L. Nemhauser, M. W. P. Savelsbergh. Conflict graphs in solving integer programming problems. European Journal of Operational Research 121, 2000, pp. 40–55.

[4] Berthold, T. Primal Heuristics for Mixed Integer Programs. Technischen Universität Berlin, September 2006. Available at https://www.zib.de/groetschel/students/Diplom-Berthold.pdf.

[5] Cornuéjols, G. Valid inequalities for mixed integer linear programs. Mathematical Programming B, Vol. 112, pp. 3–44, 2008.

[6] Danna, E., Rothberg, E., Le Pape, C. Exploring relaxation induced neighborhoods to improve MIP solutions. Mathematical Programming, Vol. 102, issue 1, pp. 71–90, 2005.

[7] Hendel, G. New Rounding and Propagation Heuristics for Mixed Integer Programming. Bachelor's thesis at Technische Universität Berlin, 2011. PDF available at https://opus4.kobv.de/opus4-zib/files/1332/bachelor_thesis_main.pdf.

[8] Mészáros C., and Suhl, U. H. Advanced preprocessing techniques for linear and quadratic programming. OR Spectrum, 25(4), pp. 575–595, 2003.

[9] Nemhauser, G. L. and Wolsey, L. A. Integer and Combinatorial Optimization. Wiley-Interscience, New York, 1999.

[10] Savelsbergh, M. W. P. Preprocessing and Probing Techniques for Mixed Integer Programming Problems. ORSA J. Computing, Vol. 6, No. 4, pp. 445–454, 1994.

[11] Wolsey, L. A. Integer Programming. Wiley-Interscience, New York, 1998.