运筹系列77：开源线性规划软件clp使用简介

Posted 2022-11-04 IE06

从本期开始阅读开源求解器coin_or的源码。整体项目介绍见https://blog.csdn.net/kittyzc/article/details/100182189
第一篇从clp开始。

1. 上手

1.1 快速使用

首先是简单的调用测试，在mac上首先安装clp的库：brew install coin-or-tools/coinor/clp，然后新建项目进行调用，各项配置如下：

使用这种方法没有办法debug细节。

1.2 命令行方式

>>> Clp
>>> Clp:import p0033.mps
>>> Clp: dualS
>>> Clp: solu sol.txt

或者干脆写成一行：

>>> clp p0033.mps -primalS
Coin LP version 1.17.3, build Dec  4 2020
command line - clp p0033.mps -primalS 
At line 15 NAME          P0033
At line 16 ROWS
At line 34 COLUMNS
At line 109 RHS
At line 118 BOUNDS
At line 152 ENDATA
Problem P0033 has 16 rows, 33 columns and 98 elements
Model was imported from ./p0033.mps in 0.000231 seconds
Presolve 15 (-1) rows, 32 (-1) columns and 97 (-1) elements
0  Obj 0 Primal inf 18.49249 (10) Dual inf 5.5987499e+11 (32)
26  Obj 2520.5717
Optimal - objective value 2520.5717
After Postsolve, objective 2520.5717, infeasibilities - dual 0 (0), primal 0 (0)
Optimal objective 2520.571739 - 26 iterations time 0.002, Presolve 0.00

关键词清单如下：

allS(lack) barr(ier) basisI(n) basisO(ut) directory dualS(implex)
 either(Simplex) end exit export gsolu(tion) guess import max(imize) min(imize)
 para(metrics) primalS(implex) printM(ask) quit restoreS(olution) saveS(olution)
 solu(tion) solv(e) stat(istics) stop

2. 基本函数

最基本的结构如下：

2.1 获取结果

调用求解的函数为primal()，求解后的信息可以通过下列函数获取：

2.2 通用求解配置

2.3 单纯形法求解配置

2.4 presolve

调用方法如下，首先加ClpPresolve.hpp

2.5 状态函数

调用和设置函数如下：

// Get row status...
Status status=model.getRowStatus(sequenceNumber)
// ... or get column status.
Status status=model.getColumnStatus(sequenceNumber)
// Set row status to basic (for example)...
model.setRowStatus(sequenceNumber,ClpSimplex::basic)
// ... or column status to basic.
model.setColumnStatus(sequenceNumber,ClpSimplex::basic)

3. 进阶函数

3.1 pivot选项

在prime simplex中，任何非基变量都可以用作pivot；在dual simplex中，任何基变量都可以用于pivot。
对于dual simplex，基础类叫做ClpDualRowPivot，实现的两个instances是ClpDualRowDantzig 和 ClpDualRowSteepest。
调用方法如下：

ClpDualRowSteepest steep(1); // 0 uninitialized, 1 compute weights
model.setDualRowPivotAlgorithm(steep);

3.2 参数矩阵

如果全是正负1，可以用ClpPlusMinusOneMatrix。建立一个CoinPackedMatrix的instance：matrix，调用方法如下：

ClpPlusMinusOneMatrix plusMinus(matrix);
assert (plusMinus.getIndices()); // would be zero if not +- one
model.loadProblem(plusMinus, lowerColumn,upperColumn, objective, lower,upper);

对于网络流问题，每一列只有一个参数，可以用ClpNetworkMatrix，调用方法如下：

ClpNetworkMatrix network(numberColumns,head,tail);
model.loadProblem(network,  lowerColumn, upperColumn, objective,  lower, upper);

此时clp会使用网络流的方法去求解。

3.3 打印日志

如果什么都不想打印，可以设置model.setLogLevel(int value)中的value为0。
也可以自定义文件日志：

FILE * fp; // assumed open
CoinMessageHandler handler(fp);
model.passInMessageHandler(&handler);

4. 示例文件说明

按照说明文档进行分类说明：

4.1 基本示例文件

• MINIMUMCPP：读取mps文件并求解
• DEFAULTSCPP：参数设置与读取示例，包括求解参数、结果参数；如何设置warm start；

/* Status of problem:
        0 - optimal
        1 - primal infeasible
        2 - dual infeasible
        3 - stopped on iterations etc
        4 - stopped due to errors
*/
std::cout << "Model status is " << model.status() 

• DRIVERCPP：设置SolveType、开启presolve

ClpSolve solvectl;
solvectl.setSolveType(ClpSolve::useDual);
# solvectl.setSolveType(ClpSolve::useBarrier);
solvectl.setPresolveType(ClpSolve::presolveOn);
model.initialSolve(solvectl);

• NETWORKCPP：如何读取非mps文件

ClpSimplex  model;
CoinPackedMatrix matrix;
int * lengths = NULL;
matrix.assignMatrix(true, numberRows, numberColumns, 2 * numberColumns, element, row, start, lengths);
ClpNetworkMatrix network(matrix);
model.loadProblem(network, lowerColumn, upperColumn, objective,
                       lower, upper);

• TESTBARRIERCPP：使用barrier方法的例子

ClpInterior  model;
int status = model.readMps( "../p0033.mps", true);
ClpCholeskyDense * cholesky = new ClpCholeskyDense();
model.setCholesky(cholesky);
model.primalDual();
// crossover的后处理算法，保证用内点法求解到的是一个顶点，而且是最近的那个顶点。
ClpSimplex model2(model);
model2.createStatus();
model2.primal(1);

4.2 其他示例文件

• DUALCUTSCPP：cut的基本使用方法
• DECOMPOSECPP：Dantzig-Wolfe decomposition的示例。
• SPRINTCPP：求解long，thin的模型。