注：中文非直接翻译英文，而是理解加工后的笔记，记录英文仅为学其专业表述。
目录
概述
Constraint optimization, or constraint programming(CP)，约束优化(规划)，用于在一个非常大的候选集合中找到可行解，其中的问题可以用任意约束来建模。
CP基于可行性(寻找可行解)而不是优化(寻找最优解)，它更关注约束和变量，而不是目标函数。
SP-SAT Solver
CP-SAT求解器技术上优于传统CP求解器。
The CP-SAT solver is technologically superior to the original CP solver and should be preferred in almost all situations.
为了增加运算速度，CP求解器处理的都是整数。
如果有非整数项约束，可以先将其乘一个整数，使其变成整数项。
If you begin with a problem that has constraints with non-integer terms, you need to first multiply those constraints by a sufficiently large integer so that all terms are integers
来看一个简单例子：寻找可行解。
- 变量x，y，z，每个只能取值0，1，2。
- eg:
x = model.NewIntVar(0, num_vals - 1, 'x')
- eg:
- 约束条件：x ≠ y
- eg:
model.Add(x != y)
- eg:
核心步骤：
- 声明模型
- 创建变量
- 创建约束条件
- 调用求解器
- 展示结果
from ortools.sat.python import cp_model
def SimpleSatProgram():
"""Minimal CP-SAT example to showcase calling the solver."""
# Creates the model.
model = cp_model.CpModel()
# Creates the variables.
num_vals = 3
x = model.NewIntVar(0, num_vals - 1, 'x')
y = model.NewIntVar(0, num_vals - 1, 'y')
z = model.NewIntVar(0, num_vals - 1, 'z')
# Creates the constraints.
model.Add(x != y)
# Creates a solver and solves the model.
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.FEASIBLE:
print('x = %i' % solver.Value(x))
print('y = %i' % solver.Value(y))
print('z = %i' % solver.Value(z))
SimpleSatProgram()
运行得
x = 1
y = 0
z = 0
status = solver.Solve(model)返回值状态含义:
- OPTIMAL 找到了最优解。
- FEASIBLE 找到了一个可行解，不过不一定是最优解。
- INFEASIBLE 无可行解。
- MODEL_INVALID 给定的CpModelProto没有通过验证步骤。可以通过调用
ValidateCpModel(model_proto)获得详细的错误。 - UNKNOWN 由于达到了搜索限制，模型的状态未知。
加目标函数，寻找最优解
假设目标函数是求x + 2y + 3z的最大值，在求解器前加
model.Maximize(x + 2 * y + 3 * z)
并替换展示条件
if status == cp_model.OPTIMAL:
运行得
x = 1
y = 2
z = 2
加回调函数，展示所有可行解
去掉目标函数，加上打印回调函数
from ortools.sat.python import cp_model
class VarArraySolutionPrinter(cp_model.CpSolverSolutionCallback):
"""Print intermediate solutions."""
def __init__(self, variables):
cp_model.CpSolverSolutionCallback.__init__(self)
self.__variables = variables
self.__solution_count = 0
def on_solution_callback(self):
self.__solution_count += 1
for v in self.__variables:
print('%s=%i' % (v, self.Value(v)), end=' ')
print()
def solution_count(self):
return self.__solution_count
def SearchForAllSolutionsSampleSat():
"""Showcases calling the solver to search for all solutions."""
# Creates the model.
model = cp_model.CpModel()
# Creates the variables.
num_vals = 3
x = model.NewIntVar(0, num_vals - 1, 'x')
y = model.NewIntVar(0, num_vals - 1, 'y')
z = model.NewIntVar(0, num_vals - 1, 'z')
# Create the constraints.
model.Add(x != y)
# Create a solver and solve.
solver = cp_model.CpSolver()
solution_printer = VarArraySolutionPrinter([x, y, z])
status = solver.SearchForAllSolutions(model, solution_printer)
print('Status = %s' % solver.StatusName(status))
print('Number of solutions found: %i' % solution_printer.solution_count())
SearchForAllSolutionsSampleSat()
运行得
x=0 y=1 z=0
x=1 y=2 z=0
x=1 y=2 z=1
x=1 y=2 z=2
x=1 y=0 z=2
x=1 y=0 z=1
x=2 y=0 z=1
x=2 y=1 z=1
x=2 y=1 z=2
x=2 y=0 z=2
x=0 y=1 z=2
x=0 y=1 z=1
x=0 y=2 z=1
x=0 y=2 z=2
x=1 y=0 z=0
x=2 y=0 z=0
x=2 y=1 z=0
x=0 y=2 z=0
Status = OPTIMAL
Number of solutions found: 18
展示 intermediate solutions
与展示所有可行解的异同主要在:
- 加
model.Maximize(x + 2 * y + 3 * z) - 换
status = solver.SolveWithSolutionCallback(model, solution_printer)
输出结果为:
x=0 y=1 z=0
x=0 y=2 z=0
x=0 y=2 z=1
x=0 y=2 z=2
x=1 y=2 z=2
Status = OPTIMAL
Number of solutions found: 5
与官网结果不一致。https://developers.google.cn/optimization/cp/cp_solver?hl=es-419#run_intermediate_sol
应该与python和or-tools版本有关。
Original CP Solver
一些特殊的小问题场景，传统CP求解器性能优于CP-SAT。
区别于CP-SAT的包，传统CP求解器用的是from ortools.constraint_solver import pywrapcp
打印方式也有区别，传统CP求解器打印方式是用while solver.NextSolution():来遍历打印。
构造器下例子中只用了一个phase，复杂问题可以用多个phase，以便于在不同阶段用不同的搜索策略。Phase()方法中的三个参数:
- vars 包含了该问题的变量，下例是 [x, y, z]
- IntVarStrategy 选择下一个unbound变量的规则。默认值 CHOOSE_FIRST_UNBOUND，表示每个步骤中，求解器选择变量出现顺序中的第一个unbound变量。
- IntValueStrategy 选择下一个值的规则。默认值 ASSIGN_MIN_VALUE，表示选择还没有试的变量中的最小值。另一个选项 ASSIGN_MAX_VALUE 反之。
import sys
from ortools.constraint_solver import pywrapcp
def main():
solver = pywrapcp.Solver("simple_example")
# Create the variables.
num_vals = 3
x = solver.IntVar(0, num_vals - 1, "x")
y = solver.IntVar(0, num_vals - 1, "y")
z = solver.IntVar(0, num_vals - 1, "z")
# Create the constraints.
solver.Add(x != y)
# Call the solver.
db = solver.Phase([x, y, z], solver.CHOOSE_FIRST_UNBOUND, solver.ASSIGN_MIN_VALUE)
solver.Solve(db)
print_solution(solver, x, y, z)
def print_solution(solver, x, y, z):
count = 0
while solver.NextSolution():
count += 1
print("x =", x.Value(), "y =", y.Value(), "z =", z.Value())
print("
Number of solutions found:", count)
if __name__ == "__main__":
main()
打印结果:
x = 0 y = 1 z = 0
x = 0 y = 1 z = 1
x = 0 y = 1 z = 2
x = 0 y = 2 z = 0
x = 0 y = 2 z = 1
x = 0 y = 2 z = 2
x = 1 y = 0 z = 0
x = 1 y = 0 z = 1
x = 1 y = 0 z = 2
x = 1 y = 2 z = 0
x = 1 y = 2 z = 1
x = 1 y = 2 z = 2
x = 2 y = 0 z = 0
x = 2 y = 0 z = 1
x = 2 y = 0 z = 2
x = 2 y = 1 z = 0
x = 2 y = 1 z = 1
x = 2 y = 1 z = 2
Number of solutions found: 18
Cryptarithmetic Puzzles
Cryptarithmetic Puzzles 是一种数学问题，每个字母代表唯一数字，目标是找到这些数字，使给定方程成立。
由
CP
+ IS
+ FUN
--------
= TRUE
找到
23
+ 74
+ 968
--------
= 1065
这是其中一个答案，还可以有其他解。
用 CP-SAT 或者 传统CP 求解器均可。
建模
约束构建如下:
- 等式 CP + IS + FUN = TRUE
- 每个字母代表一个数字
- C、I、F、T不能为0，因为开头不写0
这个问题官方解法代码有误，回头研究。
The N-queens Problem
n皇后问题是个很经典的问题。如果用深度优先搜索的解法见 link
In chess, a queen can attack horizontally, vertically, and diagonally. The N-queens problem asks:
How can N queens be placed on an NxN chessboard so that no two of them attack each other?
在国际象棋中，皇后可以水平、垂直和对角攻击。N-queens问题即: 如何把N个皇后放在一个NxN棋盘上，使它们不会互相攻击?
这并不是最优化问题，我们是要找到所有可行解。
CP求解器是尝试所有的可能来找到可行解。
传播和回溯
约束优化有两个关键元素:
- Propagation，传播可以加速搜索过程，因为减少了求解器继续做无谓的尝试。每次求解器给变量赋值时，约束条件都对为赋值的变量增加限制。
- Backtracking，回溯发生在继续向下搜索也肯定无解的情况下，故需要回到上一步尝试未试过的值。
构造约束:
- 约束禁止皇后区位于同一行、列或对角线上
queens[i] = jmeans there is a queen in column i and row j.- 设置在不同行列，
model.AddAllDifferent(queens) - 设置两个皇后不能在同一对角线，更为复杂些。
- 如果对角线是从左到右下坡的，则斜率为-1（等同于有直线 y = -x + b），故满足
queens(i) + i具有相同值的所有 queens(i) 在一条对角线上。 - 如果对角线是从左到右上坡的，则斜率为 1（等同于有直线 y = x + b），故满足
queens(i) - i具有相同值的所有 queens(i) 在一条对角线上。
- 如果对角线是从左到右下坡的，则斜率为-1（等同于有直线 y = -x + b），故满足
import sys
from ortools.sat.python import cp_model
def main(board_size):
model = cp_model.CpModel()
# Creates the variables.
# The array index is the column, and the value is the row.
queens = [model.NewIntVar(0, board_size - 1, 'x%i' % i)
for i in range(board_size)]
# Creates the constraints.
# The following sets the constraint that all queens are in different rows.
model.AddAllDifferent(queens)
# Note: all queens must be in different columns because the indices of queens are all different.
# The following sets the constraint that no two queens can be on the same diagonal.
for i in range(board_size):
# Note: is not used in the inner loop.
diag1 = []
diag2 = []
for j in range(board_size):
# Create variable array for queens(j) + j.
q1 = model.NewIntVar(0, 2 * board_size, 'diag1_%i' % i)
diag1.append(q1)
model.Add(q1 == queens[j] + j)
# Create variable array for queens(j) - j.
q2 = model.NewIntVar(-board_size, board_size, 'diag2_%i' % i)
diag2.append(q2)
model.Add(q2 == queens[j] - j)
model.AddAllDifferent(diag1)
model.AddAllDifferent(diag2)
### Solve model.
solver = cp_model.CpSolver()
solution_printer = SolutionPrinter(queens)
status = solver.SearchForAllSolutions(model, solution_printer)
print()
print('Solutions found : %i' % solution_printer.SolutionCount())
class SolutionPrinter(cp_model.CpSolverSolutionCallback):
"""Print intermediate solutions."""
def __init__(self, variables):
cp_model.CpSolverSolutionCallback.__init__(self)
self.__variables = variables
self.__solution_count = 0
def OnSolutionCallback(self):
self.__solution_count += 1
for v in self.__variables:
print('%s = %i' % (v, self.Value(v)), end = ' ')
print()
def SolutionCount(self):
return self.__solution_count
class DiagramPrinter(cp_model.CpSolverSolutionCallback):
def __init__(self, variables):
cp_model.CpSolverSolutionCallback.__init__(self)
self.__variables = variables
self.__solution_count = 0
def OnSolutionCallback(self):
self.__solution_count += 1
for v in self.__variables:
queen_col = int(self.Value(v))
board_size = len(self.__variables)
# Print row with queen.
for j in range(board_size):
if j == queen_col:
# There is a queen in column j, row i.
print("Q", end=" ")
else:
print("_", end=" ")
print()
print()
def SolutionCount(self):
return self.__solution_count
if __name__ == '__main__':
# By default, solve the 8x8 problem.
board_size = 8
if len(sys.argv) > 1:
board_size = int(sys.argv[1])
main(board_size)
欲图形化打印，调用SolutionPrinter处换成DiagramPrinter。
Setting solver limits
时间限制
eg:
solver.parameters.max_time_in_seconds = 10.0
找到指定数量的解后停止搜索
eg:
from ortools.sat.python import cp_model
class VarArraySolutionPrinterWithLimit(cp_model.CpSolverSolutionCallback):
"""Print intermediate solutions."""
def __init__(self, variables, limit):
cp_model.CpSolverSolutionCallback.__init__(self)
self.__variables = variables
self.__solution_count = 0
self.__solution_limit = limit
def on_solution_callback(self):
self.__solution_count += 1
for v in self.__variables:
print('%s=%i' % (v, self.Value(v)), end=' ')
print()
if self.__solution_count >= self.__solution_limit:
print('Stop search after %i solutions' % self.__solution_limit)
self.StopSearch()
def solution_count(self):
return self.__solution_count
def StopAfterNSolutionsSampleSat():
"""Showcases calling the solver to search for small number of solutions."""
# Creates the model.
model = cp_model.CpModel()
# Creates the variables.
num_vals = 3
x = model.NewIntVar(0, num_vals - 1, 'x')
y = model.NewIntVar(0, num_vals - 1, 'y')
z = model.NewIntVar(0, num_vals - 1, 'z')
# Create a solver and solve.
solver = cp_model.CpSolver()
solution_printer = VarArraySolutionPrinterWithLimit([x, y, z], 5)
status = solver.SearchForAllSolutions(model, solution_printer)
print('Status = %s' % solver.StatusName(status))
print('Number of solutions found: %i' % solution_printer.solution_count())
assert solution_printer.solution_count() == 5
StopAfterNSolutionsSampleSat()