Scipy：使用稀疏矩阵进行线性规划

Question

我希望在Python中解决一个线性规划问题。变量的数量（我将其称为N）非常大（约50000），为了以scipy.optimize.linprog要求的方式来表达问题，我必须构建两个N x N矩阵（如下所示的A和B）。该LP可以写成：

minimize: c.x
subject to:
    A.x <= a
    B.x  = b
    x_i >= 0 for all i in {0, ..., n}

其中 . 表示点积，a、b 和 c 是长度为 N 的向量。

我的经验是构建这样的大矩阵（A 和 B 都有大约 50000x50000 = 25*10^8 个元素）会遇到一些问题：如果硬件不太强，NumPy 可能会拒绝构建如此大的矩阵（例如参见 Very large matrices using Python and NumPy），即使 NumPy 没有出现问题地创建矩阵，也存在巨大的性能问题。考虑到 NumPy 处理的数据量非常庞大，这是很自然的。

然而，尽管我的线性程序具有 N 个变量，但我处理的矩阵非常稀疏。其中一个矩阵只有第一行有条目，另一个矩阵只有前 M 行有条目，其中 M < N/2。当然，我希望利用这一事实。

据我所读（例如尝试使用稀疏矩阵解决Scipy优化问题并失败），scipy.optimize.linprog不能处理稀疏矩阵。因此，我有以下问题：

• SciPy是否真的没有提供使用稀疏矩阵解决线性规划问题的可能性？（如果不是，我该如何做？）
• 你知道任何替代库可以比使用非稀疏矩阵的SciPy更有效地解决问题吗？（在上面的thread中建议的library似乎对我的目的来说不够灵活-据我了解其文档）
• 是否可以期望使用纯Python（无C）实现的单纯形算法的新实现将比使用非稀疏矩阵的SciPy更有效？

- Samufi

Answer 1

我认为形成一个密集的矩阵（或两个矩阵）来解决大规模稀疏线性规划可能不是正确的做法。在解决大规模稀疏线性规划时，重要的是使用能够处理此类问题并以不显式创建任何这些零元素的方式生成模型的求解器。

编写一个稳定、快速、稀疏的Python Simplex LP求解器作为SciPy密集求解器的替代品并不是一项微不足道的任务。而且，用纯Python编写的求解器可能性能不佳。

对于您所指示的规模，虽然不是非常非常大（可能是大中等大小的模型），但您可能希望考虑像Cplex、Gurobi或Mosek这样的商业求解器。这些求解器非常快速和可靠（它们可以解决您提出的任何LP问题）。它们都有Python API。对于学术界来说，这些求解器是免费或非常便宜的。

如果你想使用开源求解器，可以考虑使用COIN CLP求解器。它还有Python接口。

如果你的模型更加复杂，你也可以考虑使用Python建模工具，比如Pulp或Pyomo（Gurobi在Python中也有很好的建模支持）。

- Erwin Kalvelagen

Answer 2

我无法相信没有人向你指出PuLP的方向！你将能够高效地创建你的问题，像这样：

import pulp
prob = pulp.LpProblem("test problem",pulp.LpMaximize)
x = pulp.LpVariable.dicts('x', range(5), lowBound=0.0)
prob += pulp.lpSum([(ix+1)*x[ix] for ix in range(5)]), "objective"
prob += pulp.lpSum(x)<=3, "capacity"

prob.solve()
for k, v in prob.variablesDict().iteritems():
    print k, v.value()

PuLP非常棒，自带一个非常不错的求解器（CBC），可以连接开源和商业求解器。我目前正在为一家林业公司使用它进行生产，并尝试使用Dippy来解决我们遇到的最难的（整数）问题。祝好运！

- Sergio Lucero