Python - 规划问题求解器概览

(图: 装箱问题图示)

一. 前言

在【学界】运筹学数学规划|离散优化求解器大搜罗 - 知乎 (zhihu.com)这篇文章中, 作者罗列了当前市面上的大部分主流(付费)的求解器, 缺少对开源(免费)的(python)库的介绍, 这里做一个补充.

几大扛把子(付费)就不需要多提了, 这里主要看看开源(免费)的几个库, 目前找到的几个start/fork数量较多, 文档较为完善的库.

就上述的状态来看, Google出品的ortools和微软的z3-solver得到了较多的关注, 但是z3针对的问题和其他的库有一点差异.

多数情况下, 除非特别大的模型(如数以万计的变量)或者是一些冷僻模型, 解决问题的核心多卡在数学建模(基础) - 模型转为代码(主体框架)这个过程型(即如何将数学模型转为代码表示的形式, 这看上去没什么难度, 但是极为容易犯错的一环)上, 不用太过于纠结使用哪个库. 至于使用哪个库作为主力, 主要考虑库api设计的使用难度, 提供文档的友好性(核心部分), 库维护状态. 同时需要看该库的使用用户数量, 问题提交数等, 这些间接反映库质量和维护的状态.

需要注意的是上述(开源)库, 均未有中文文档提供(google的ortools提供的也是机器翻译, 未经人工修正的中文文档).

补充: 求解器：助力智能决策的利器 - 知乎, 这篇文章更详细的罗列来了各类求解器.

1 - 商用求解器

2 - 开源求解器

3 - 软件集成

二. 问题

下面全部的代码测试, 均使用此题目.

2.1 线性规划

$max\;z = 2x_1 + 3x_2 - 5x_3\\ s.t = \begin{cases} x_1 + x_2 + x_3 = 7\\ 2x_1 - 5x_2 + x_3 \geq 10\\ x_1 + 3x_2 + x_3 \leq 12\\ x_1, x_2, x_3 \geq 0 \end{cases}\\$

示例, 使用的是Python关于数学的处理 | Lian (kyouichirou.github.io)该文的案例, 在此文中使用scipy来求解.

from scipy import optimize

# 矩阵形式
# * -1, 转换为min
z = [-2, -3, 5]

a = [
    [1, 1, 1],
    # * -1, 转换为小于等于
    [-2, 5, -1],
    [1, 3, 1]
]

b = [
    7, -10, 12
]

# 约束
x1b = (0, None)
x2b = (0, None)
x3b = (0, None)

optimize.linprog(
    z,
    A_ub=a, b_ub=b,
    bounds=[x1b, x2b, x3b],
    method='revised simplex'
)

     con: array([], dtype=float64)
     fun: -14.571428571428571
 message: 'Optimization terminated successfully.'
     nit: 2
   slack: array([0.        , 0.        , 3.85714286])
  status: 0
 success: True
       x: array([6.42857143, 0.57142857, 0.        ])

可以对比以下几个库使用和scipy之间的差异.

scipy.optimize模块提供的求解支持.

三. 安装与使用

conda的下载源速度较慢, 建议将下载源改成国内的镜像站, 清华大学这个开源镜像站做的很好.

在.condarc写入(位于用户文件夹之下, windows)

channels:
  - defaults
show_channel_urls: true
default_channels:
  - https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/main
  - https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/r
  - https://mirrors.tuna.tsinghua.edu.cn/anaconda/pkgs/msys2
custom_channels:
  conda-forge: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  msys2: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  bioconda: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  menpo: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  pytorch: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  pytorch-lts: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  simpleitk: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud
  deepmodeling: https://mirrors.tuna.tsinghua.edu.cn/anaconda/cloud/

# 删除目录缓存
conda clean -i

建议在虚拟环境下测试.

3.1 ortools

个人使用或更为青睐的, 放在前面.

Google OR-Tools Reference | Google for Developers

OR-Tools is an open source software suite for optimization, tuned for tackling the world's toughest problems in vehicle routing, flows, integer and linear programming, and constraint programming.

After modeling your problem in the programming language of your choice, you can use any of a half dozen solvers to solve it: commercial solvers such as Gurobi or CPLEX, or open-source solvers such as SCIP, GLPK, or Google's GLOP and award-winning CP-SAT.

这个库没有提供conda的安装方式, 比较奇怪的是, 在Ortools Python :: Anaconda.org找到的库名称发生变更, 但是安装却出现找不到库的问题.

@Lian ➜ ~\Desktop ( org_env 3.10.13) 1.569s conda install -c conda-forge ortools-python
Collecting package metadata (current_repodata.json): done
Solving environment: failed with initial frozen solve. Retrying with flexible solve.
Collecting package metadata (repodata.json): done
Solving environment: failed with initial frozen solve. Retrying with flexible solve.

PackagesNotFoundError: The following packages are not available from current channels:

  - ortools-python

Current channels:

  - https://conda.anaconda.org/conda-forge/win-64
  - https://conda.anaconda.org/conda-forge/noarch
  - https://repo.anaconda.com/pkgs/main/win-64
  - https://repo.anaconda.com/pkgs/main/noarch
  - https://repo.anaconda.com/pkgs/r/win-64
  - https://repo.anaconda.com/pkgs/r/noarch
  - https://repo.anaconda.com/pkgs/msys2/win-64
  - https://repo.anaconda.com/pkgs/msys2/noarch

To search for alternate channels that may provide the conda package you're
looking for, navigate to

    https://anaconda.org

and use the search bar at the top of the page.

改回原站点源还是出现找不到包错误.

3.1.1 安装

ortools - PyPI提供好了编译的文件, 直接下载下来安装即可.

conda create -n org_env python=3.10

@Lian ➜ ~\Desktop ( base 3.9.12)  conda activate org_env

pip install ortools-9.7.2996-cp310-cp310-win_amd64.whl

@Lian ➜ D:\common_download ( org_env 3.10.13) 14.885s python
Python 3.10.13 | packaged by Anaconda, Inc. | (main, Sep 11 2023, 13:24:38) [MSC v.1916 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import ortools as cr
>>> cr.__version__
'9.7.2996'

3.1.2 示例

ortools的使用比较灵活

from ortools.linear_solver import pywraplp

model = pywraplp.Solver.CreateSolver("CLP")
x1 = model.NumVar(0, model.infinity(), name='x1')
x2 = model.NumVar(0, model.infinity(), name='x2')
x3 = model.NumVar(0, model.infinity(), name='x3')

model.Add(x1 + x2 + x3 == 7.0)
model.Add(2 * x1 - 5 * x2 + x3 >= 10)
model.Add(x1 + 3 * x2 + x3 <= 12)

model.Maximize(2 * x1 + 3*x2 - 5*x3)

status = model.Solve()

model.Objective().Value()

# 14.57142857142857

二者是等价的, 但是下面的看起来很复杂, 实际上只是代码多了点而已.

from ortools.linear_solver import pywraplp

model = pywraplp.Solver('LinearExample', pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)
x1 = model.NumVar(0, model.infinity(), name='x1')
x2 = model.NumVar(0, model.infinity(), name='x2')
x3 = model.NumVar(0, model.infinity(), name='x3')

c1 = model.Constraint(7, 7)
c1.SetCoefficient(x1, 1)
c1.SetCoefficient(x2, 1)
c1.SetCoefficient(x3, 1)

c2 = model.Constraint(10, model.infinity())
c2.SetCoefficient(x1, 2)
c2.SetCoefficient(x2, -5)
c2.SetCoefficient(x3, 1)

c3 = model.Constraint(-model.infinity(), 12)
c3.SetCoefficient(x1, 1)
c3.SetCoefficient(x2, 3)
c3.SetCoefficient(x3, 1)

objective = model.Objective()
objective.SetCoefficient(x1, 2)
objective.SetCoefficient(x2, 3)
objective.SetCoefficient(x3, -5)
objective.SetMaximization()

model.Solve()

print(objective.Value())

# 14.57142857142857

3.1.3 小结

Google出品, 为这个库提供了"强大/可靠"的背书, 从使用上来看, 这个库的api设计也比较灵活, 更为重要的是, 这个同时整合了其他的求解器.

如MIP(Mixed Integer Programming)(混合整数规划)问题, 使用的是scip的求解器

提供了绝大部分运筹学常见问题的求解支持.

总体而言, 这是一个很不错的库.

但是, 这库存在一个较大的问题: 文档友好性.

如: 搜索的结果, 找不到, 找不准, 幸亏示例还算丰富.

例如想快速定位某个函数的参数解析, 颇为不便.

3.2 pyscipopt

PySCIPOpt: Overview

文档应该是这几个库中最为完善和严谨的了和这个库的开发背景为研究机构很大关系啊!

This project provides an interface from Python to the SCIP Optimization Suite. Please review SCIP's license restrictions before installing PySCIPOpt.

SCIP is currently one of the fastest non-commercial solvers for mixed integer programming (MIP) and mixed integer nonlinear programming (MINLP). It is also a framework for constraint integer programming and branch-cut-and-price. It allows for total control of the solution process and the access of detailed information down to the guts of the solver.

需要搭配scip使用, 但是需要注意并不一定需要安装scip的exe文件, 使用conda可以直接安装, 需要注意对应的版本.

PySCIPOpt/INSTALL.md at master - scipopt/PySCIPOpt - GitHub提供的安装方式, pip的安装似乎有问题, 尽管按照文档提供的方法, 也设置了环境变量等操作了, 但是一直报错.

# 提示编译环境出问题
# 安装编译器
# 不想安装大块头的 vc studio c++
conda install libpython m2w64-toolchain -c msys2
# 尽管添加还是不行

文档的建议也是使用conda来安装, 还是让conda来解决问题好一点.

3.2.1 安装

创建一个虚拟环境

conda create -n scip_env python=3.10

conda下并不需要手动安装scip文件

conda install --channel conda-forge pyscipopt

默认状态下, 安装的版本并不是最新, 需要手动指定版本号:

conda install --channel conda-forge pyscipopt==4.3.0

@Lian ➜ ~\Desktop ( scip_env 3.10.13) 1.24s python
Python 3.10.13 | packaged by Anaconda, Inc. | (main, Sep 11 2023, 13:24:38) [MSC v.1916 64 bit (AMD64)] on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import pyscipopt as sp
>>> sp.__version__
'4.3.0'

3.2.2 示例

from pyscipopt import Model

# 模型初始化
model = Model('tets')

# 创建变量
x1 = model.addVar(vtype='C', name='x1')
x2 = model.addVar(vtype='C', name='x2')
x3 = model.addVar(vtype='C', name='x3')

# 添加约束
model.addCons(x1 + x2 + x3 == 7)
model.addCons(2 * x1 - 5 * x2 + x3 >= 10)
model.addCons(x1 + 3 * x2 + x3 <= 12)
model.addCons(x1 >= 0)
model.addCons(x2 >= 0)
model.addCons(x3 >= 0)

# 设置目标函数
model.setObjective(2 * x1 + 3*x2 - 5*x3, 'maximize')
# 求解
model.optimize()
# 读取最佳解
sol = model.getBestSol()
# 最优目标值
target=model.getObjVal(model)

# 14.571428571428571

可以看到, 这种api的设计和代码逻辑和上述的数学模型执行逻辑非常相近, 同时代码的风格和sympy非常之相像, 近似于数学符号的使用, 个人认为这种风格的代码是比较友好, 在数学模型转换为代码时这个过程更为自然.

3.3 cvxpy

CVXPY is a Python-embedded modeling language for convex optimization problems. It automatically transforms the problem into standard form, calls a solver, and unpacks the results.

Welcome to CVXPY 1.4 - CVXPY 1.4 documentation

文档是python用户熟悉的味道, 这一点非常友好.

3.3.1 安装

conda install -c conda-forge cvxpy

安装没什么特别的地方, 直接conda安装即可.

3.3.2 示例

import cvxpy as cp
import numpy as np

# 创建决策变量
x = cp.Variable(3)

# 其他的和scipy的基本一样
z = np.array([-2, -3, 5])

a = np.array([
    [1, 1, 1],
    # * -1, 转换为小于等于
    [-2, 5, -1],
    [1, 3, 1]
])

b = np.array([
    7, -10, 12
])

# 约束的使用相对麻烦
c = np.array([
    [1, 0, 0],
    [0, 1, 0],
    [0, 0, 1]
])

d = np.array([0, 0, 0])

prob = cp.Problem(cp.Minimize(z.T @ x), [a @ x <= b,  c @ x >= d])
print(prob.solve(solver = cp.ECOS_BB))

# -14.571428571127495

>>> import cvxpy as cp
(CVXPY) Oct 20 07:31:43 PM: Encountered unexpected exception importing solver OSQP:
ImportError('DLL load failed while importing qdldl: The specified module could not be found.')

感觉用起来没有scipy的感觉好, 但是很明显cvxpy的api设计更为严谨, 前后一致(都是使用矩阵的形式, 尽管会在使用上造成一定的麻烦, 但是却可以做到前后的统一).

3.4 pulp

PuLP is an LP modeler written in Python. PuLP can generate MPS or LP files and call GLPK, COIN-OR CLP/CBC, CPLEX, GUROBI, MOSEK, XPRESS, CHOCO, MIPCL, SCIP to solve linear problems.

Optimization with PuLP - PuLP 2.7.0 documentation (coin-or.github.io)

对这个库不是很感兴趣, 暂未测试

3.5 z3-solver

Python: package z3

Z3 is an SMT solver and supports the SMTLIB format.

Z3 is an efficient Satisfiability Modulo Theories (SMT) solver from Microsoft Research. Z3 is a solver for symbolic logic, a foundation for many software engineering tools. SMT solvers rely on a tight integration of specialized engines of proof. Each engine owns a piece of the global puzzle and implements specialized algorithms. For example, Z3’s engine for arithmetic integrates Simplex, cuts and polynomial reasoning, while an engine for strings are regular expressions integrate methods for symbolic derivatives of regular languages. A theme shared among many of the algorithms is how they exploit a duality between finding satisfying solutions and finding refutation proofs. The solver also integrates engines for global and local inferences and global propagation. Z3 is used in a wide range of software engineering applications, ranging from program verification, compiler validation, testing, fuzzing using dynamic symbolic execution, model-based software development, network verification, and optimization.

微软这个和其他有有点差异, 同时文档有点狂野, 这里不做测试.

>>> x = Int('x')
>>> y = Int('y')
>>> s = Solver()
>>> s.add(x > 0)
>>> s.add(x < 2)
>>> s.add(y == x + 1)
>>> s.check()
sat
>>> m = s.model()
>>> m[x]
1
>>> m[y]
2

看了文档的示例, 还算容易理解.

乎此之外, z3支持JavaScript[Z3 JavaScript | Online Z3 Guide (microsoft.github.io)](https://microsoft.github.io/z3guide/programming/Z3 JavaScript Examples), 这是比较少见的.

四. 总结

技之道, 顺时应物; 术之道, 万变不离其中.