如何用Python解决具有大量决策变量的线性问题
由于涉及大量变量，某些组合问题引起的线性规划（LP）问题变得难以处理（Gilmore＆Gomory，1961）。在这种问题中，延迟列生成是一种潜在的解决方案，因为它不会从一开始就在问题中包含所有可能的决策变量。一个经典的应用是切割库存问题，在其中必须决定如何将给定宽度的卷切成更小的块，以满足确定的切割尺寸的需求。
本文将介绍列生成的一些主要理论方面，并用一维切割库存问题进行说明。在解决方案中，将使用Python库scipy。有兴趣的读者也可以在此git存储库中找到完整的代码。
如果您还不熟悉线性规划，则在阅读有关该主题的先前介绍后，您可能会更好地理解本文。
线性规划：理论与应用
线性优化的主要概念和Python实现
towardsdatascience.com
如果您已准备好深入探讨概念，并且需要一个实际的示例，请加入我们。
列生成：概述
当使用延迟列生成来解决LP时，每次迭代都将考虑仅具有子集列的问题。它被称为受限主问题（RMP）。在每次迭代中，解决RMP并获得其相应的最优对偶决策变量向量π。然后，算法使用此结果来解决子问题或定价问题，其旨在找到具有负减少成本的新决策变量。让我们回顾一下减少成本的定义。
对于任何非基本变量，变量的减少成本是非基本变量的目标函数系数必须改善的金额，以使该变量成为LP某个最优解的基本变量。Winston＆Goldberg（2004年）
因此，通过包括具有负减少成本的变量，我们可能会期望边际贡献提高目标值。回顾一下对偶变量作为影子成本的经济解释。新变量可能会导致节省，其节省与其对应的双重物成本的贡献成比例，同时通过其成本系数增加总成本。
通过解决子问题，我们确定了一组具有给定双重成本的最低减少成本的列S。如果未识别出具有负减少成本的列，则停止，因为π是原始问题的最优对偶解，与RMP的最优原始解一起，我们有一个最优原始/对偶对。否则，我们将列追加到RMP中并进行迭代（Klabjan，2005）。此过程由下图表示。
请注意，子问题是问题特定的，并且在某些情况下，将其制定为混合整数规划问题可能相当耗时。然后考虑可以返回每次迭代中多个负减少成本列的启发式和/或动态规划方法可能有用。在车辆路径问题的情况下，子问题通常是受约束的最短路径问题。有兴趣深入了解的读者可以参考Irnich＆Desaulniers（2005）以获得解决方案技术的洞察力。
现在请记住，所提出的延迟列生成方法解决了具有实值决策变量的线性规划问题。为了解决大型整数或混合整数规划问题，可以在解决松弛后使用一些启发式方法，或者使用Branch＆Price来施加整数约束。在后一种方法中，应考虑不仅在解决初始LP时已生成的那些列，还要在树的整个LP解决过程中生成新列。 Barnhart等人对Branch＆Price进行了更深入的讨论。 （1998）。作者研究了常见问题，案例研究和关于分支规则的有趣见解。
切割库存问题
假设我们有一组需求I，每个需求需要d个宽度为w的零件。此外，假设我们有一卷宽度为W的卷材，从中进行切割。已知切割模式集合将被表示为P。每个切割模式p的每个零件消耗一个单位的卷材（c = 1），并产生aᵢₚ个wᵢ的宽度单元。我们的目标是确定应使用每个模式p的切割次数x，以满足需求并最小化消耗的单位数。我们可以如下公式化这个问题。
与需求约束相关的对偶决策变量π然后用于定价问题（子问题）中寻找负减少成本的新模式。在切割库存问题的情况下，我们必须找到一种新的模式，将不同宽度的零件组合在一起，适合总宽度W，并通过帮助满足材料需求，将导致比新成本更多的节省。这是一个背包问题，可以如下公式化。
其中yᵢ对应于新切割模式中生产的宽度为wᵢ的零件数量。
由于我们知道每个新模式都有一个单位成本c = 1，我们可以通过以下方式计算新模式的减少成本：
在非负值的情况下，应停止列生成算法。
为了获得切割库存问题的整数解，一个简单的启发式方法是将线性松弛中获得的分数值向上取整。或者，也可以在强制实现整数约束的线性松弛中使用产生的模式集合来解决线性问题。我们将在本文中使用这两种策略。对于线性松弛接近完整整数模型的实例，这些策略可能非常成功。在需求量相对较小的其他实例中可能会出现一些差异。
如果要获得精确的整数解，分支定价法可能是一个不错的选择。在此方法中，在某些初始变量的分支后，可能会包括当前节点中具有降低成本的新列。对此感兴趣的人可以参考Carvalho（1998）和Vance（1998）获取更多详细信息。
现在让我们动手实践一些吧！
解决方案
让我们开始解决切割库存问题的Python实现，其中LP松弛被解决为最优解，并且已经产生的模式的整数模型被解决。我们将使用numpy进行线性代数运算，pandas处理数据框，scipy进行优化算法，以及matplotlib将切割模式可视化到图表中。
import numpy as npimport pandas as pdfrom scipy.optimize import linprogimport matplotlib.pyplot as plt让我们导入一个数据集，这是JuMP文档中提出的问题的修改版本。
dataset = pd.read_csv("data.txt", sep=" ")让我们实例化问题参数
# Total widthW = 100.0# Width and amount associated with each demandw = dataset.w.valuesd = dataset.d.values# LP parametersA = np.eye(dataset.shape[0]) * (W // w)c = np.ones_like(w)请注意，为了初始化A矩阵，我引入了产生每个宽度的最大可行滚动数的切割简单模式。假设有一些宽度为24的卷的需求。鉴于总宽度W为100，它将导致具有系数4的初始切割模式。
现在，让我们定义一个函数来解决子问题，给定总宽度 W，每个需求 w 关联的宽度向量和当前的双重变量 duals 向量。
def solve_knapsack(W, w, duals): return linprog( -duals, A_ub=np.atleast_2d(w), b_ub=np.atleast_1d(W), bounds=(0, np.inf), integrality=1, )在主循环中，对于受限制的主问题的每个解，我们将使用来自 scipy 的 linprog 函数。A 矩阵以及需求向量都以它们的负数作为传递，因为按照惯例，scipy 只考虑“小于或等于”不等式。
linprog(c, A_ub=-A, b_ub=-d, bounds=(0, None))解决方案应该有与需求相关的双重变量的负值存储在 ineqlin.marginals 属性中。
通过探索对偶概念，我们还可以通过解决以下 LP 来获得双重变量。
linprog(-d, A_ub=A.T, b_ub=c, bounds=(0, None))让我们在主循环中将它们放在一起。
# Initial solutionsol = linprog(c, A_ub=-A, b_ub=-d, bounds=(0, None))sol_dual = linprog(-d, A_ub=A.T, b_ub=c, bounds=(0, None))diff = np.abs(sol_dual.x + sol.ineqlin.marginals).sum()print(f"Compare duality difference: {diff}")# Iteratefor _ in range(1000): duals = -sol.ineqlin.marginals price_sol = solve_knapsack(W, w, duals) y = price_sol.x if 1 + price_sol.fun < -1e-4: print(f"Iteration: {_}; Reduced cost: {(1 + price_sol.fun):.3f}") A = np.hstack((A, y.reshape((-1, 1)))) c = np.append(c, 1) sol = linprog(c, A_ub=-A, b_ub=-d, bounds=(0, None)) else: break最后，让我们尝试将线性松弛解的解向上取整，并仅考虑在 LP 中产生的列来获得整数解。
sol_round = linprog(c, A_ub=-A, b_ub=-d, bounds=(0, np.inf), integrality=0)print(f"Rounding solution {np.ceil(sol_round.x).sum()}")sol = linprog(c, A_ub=-A, b_ub=-d, bounds=(0, np.inf), integrality=1)print(f"Integer solution: {sol.x.sum()}")应该返回：
- Rounding solution 339.0
- Integer solution: 334.0
在这种情况下，仅通过对 LP 强制实数约束而不是仅将松弛向上取整就可以将结果提高近 2%。对于那些想尝试分支定价的人，这是这个实例的精确解。
最后，让我们尝试可视化新的切割模式：
fig, ax = plt.subplots(figsize=[7, 3], dpi=100)hmap = ax.imshow(A > 1e-6, cmap="Blues")plt.axis('off')fig.tight_layout()plt.show()
我们使用列生成解决了切割库存问题。
进一步阅读
容量车辆路径问题（CVRP）最初由 Dantzig & Ramser（1959）提出，由于其组合性质而特别具有挑战性。作者在他们的原始论文中显示，即使对于小问题，可能的路径数量也非常大。例如，具有 15 个需求点的对称问题具有超过 6×10⁹ 个可能的路径。我发现看到随着时间的推移如何在这个和相关问题中探索列生成特别有趣。
虽然针对具有严格时间窗口限制的车辆路径问题，列生成自Desrochers等人的工作（1992年）以来已经得到了很好的建立，但分支定价算法在较少约束的情况下往往会失败。因此，纯列生成并不被认为是CVRP的一种有前途的方法。但是，Fukasawa等人（2006）将列生成与分支剪枝算法相结合，证明了文献中几个实例的最优性。其他作者进一步改进了CVRP的分支剪枝定价，我相信Pecin等人（2017年）的工作对感兴趣的读者尤其有吸引力。
结论
本文介绍了延迟列生成作为解决具有大量决策变量的线性规划的策略，而不需要明确考虑所有变量。经典的一维切割问题用来说明该方法，并在Python中实现了一种解决方案替代方案。完整代码可在此git存储库中获取。
参考文献
Barnhart，C.，Johnson，E. L.，Nemhauser，G. L.，Savelsbergh，M. W.和Vance，P. H.，1998。分支定价：用于解决巨大整数规划的列生成。运营研究，46（3），316-329。
de Carvalho，J. V.，1998。使用列生成和分支定界精确解决切割问题。国际运营研究交易，5（1），35-44。
Dantzig，G. B。和Ramser，J. H.，1959。卡车调度问题。管理科学，6（1），80-91。
Desrochers，M.，Desrosiers，J.和Solomon，M.，1992。一种新的优化算法，用于具有时间窗口的车辆路径问题。运营研究，40（2），342-354。
Fukasawa，R.，Longo，H.，Lysgaard，J.，Aragão，M. P. D.，Reis，M.，Uchoa，E。和Werneck，R. F.，2006。用于容量车辆路径问题的鲁棒分支定价。数学规划，106，491-511。
Gilmore，P. C.和Gomory，R. E.，1961。用于切割问题的线性规划方法。运营研究，9（6），849-859。
Irnich，S。和Desaulniers，G.，2005。具有资源约束的最短路径问题（第33-65页）。Springer US。
Klabjan，D.，2005。航空业中的大规模模型。列生成，163-195。
Pecin，D.，Pessoa，A.，Poggi，M。和Uchoa，E.，2017。用于容量车辆路径的改进的分支剪枝定价。数学规划计算，9，61-100。
Vance，P. H.，1998。用于一维切割问题的分支定价算法。计算优化和应用，9，211-228。
Winston，W. L。＆Goldberg，J. B.，2004。运筹学：应用和算法。第4版。加利福尼亚州贝尔蒙特：汤姆森布鲁克/科尔贝尔蒙特。