Constraint programming (CP) represents a different approach to optimization than that of classical optimization theory. CP is based on feasibility (finding feasible solutions to a problem), and focuses on the constraints and variables rather than the objective function. For many types of problems, CP can find optimal solutions much faster than a MIP solver.
The following sections show how to solve the same integer optimization problem that was previously solved with the MIP solver, but this time using the CP solver. The problem is shown below:
- Maximize x + 10y subject to the following constraints:
-
x + 7 y ≤ 17.5 x ≤ 3.5 x ≥ 0 y ≥ 0 x, y integers
The following sections present a Python program that solves this problem using the CP solver. Some sections of the program are quite similar to the program that uses the MIP solver, but there are some significant differences, particularly in defining the constraints and objective. If you're planning on using both solvers, it's worthwhile to compare these two examples.
Create the solver and variables
The following code, which creates the solver and variables for the problem, is almost the same as the corresponding MIP code.
from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import solver_parameters_pb2
def main():
# Instantiate a CP solver.
parameters = pywrapcp.Solver.DefaultSolverParameters()
solver = pywrapcp.Solver("simple_CP", parameters)
# x and y are integer non-negative variables.
x = solver.IntVar(0, 17, 'x')
y = solver.IntVar(0, 17, 'y')
The command solver.IntVar(0, 17, 'x') creates an integer variable whose values
are between 0 and 17. The intervals for the variables just need to be large enough to
include all the points in the feasible region. In this case (because the variables are
positive), you can set the upper bound to be the largest constant on the right hand side of
the constraint inequalities, which is 17.5, and round down to the nearest integer. (You
could choose a smaller bound for x, but this is not important to the solver.)
Define the constraints
The following code defines the constraints for the problem:
solver.Add(2*x + 14*y <= 35) solver.Add(2*x <= 7)
To increase computational speed, the CP solver works over the integers, so all of its inputs must be integers. To handle constraints that have non-integer constants, like the ones in this problem, you must first multiply the entire constraint by a sufficiently large positive integer so that the result has all integer constants. This does not change the feasible solutions to the problem.
In the current example, you can multiply the constraints
-
x + 7 y ≤ 17.5 x ≤ 3.5
-
2x + 14 y ≤ 35 2x ≤ 7
Note that unlike the MIP solver, it is not necessary to add the coefficients for the constraints individually.
Define the objective
The following code defines the objective for the problem:
obj_expr = solver.IntVar(0, 1000, "obj_expr") solver.Add(obj_expr == x + 10*y) objective = solver.Maximize(obj_expr, 1)
Decision builders
Before calling the solver, you need to create a decision builder for the problem. The
decision builder creates the search tree and
determines the order in which the solver searches solutions. The following code creates the
decision builder using the solver's Phase method. (The reason for the term "phase"
is that in more complicated problems,
the search can involve multiple phases, each of which employs different techniques for
finding solutions.)
decision_builder = solver.Phase([x, y],
solver.CHOOSE_FIRST_UNBOUND,
solver.ASSIGN_MIN_VALUE)
The Phase method has three input parameters:
- The first parameter contains the decision variables — the variables the solver uses to decide which node of the tree to visit next. In this case, the decision variables are simply x and y.
- The second parameter specifies how the solver chooses the next variable for
the search. Here the code uses
CHOOSE_FIRST_UNBOUND, which means the solver chooses the first unbound variable. Other options areCHOOSE_LOWEST_MINandCHOOSE_RANDOM. - The third parameter specifies how the solver chooses the next value of the
current search variable to check. Here, we've chosen
ASSIGN_MIN_VALUE, which means the solver chooses the smallest value that hasn't already been searched. We could also have chosenASSIGN_MAX_VALUEorASSIGN_RANDOM_VALUE.
In this example, the decision builder doesn't play a crucial role, because the search tree is so small. But in larger problems, choosing the right parameters for the decision builder can speed up the search significantly.
Create a solution collector
The CP solver normally returns many feasible solutions to an optimization problem. A useful tool for storing these results is the solution collector. The following code shows how to create one.
# Create a solution collector. collector = solver.LastSolutionCollector() # Add the decision variables. collector.Add(x) collector.Add(y) # Add the objective. collector.AddObjective(obj_expr)
In this example, the solution collector stores just the last solution returned by the solver. Since this is a maximization problem, each solution returned by the solver has a larger objective value than the previous one. As a result, the last solution returned by the solver is the optimal solution.
Another useful collector is the AllSolutionCollector, which stores all solutions
returned by the solver.
Invoke the solver
The following code invokes the solver and displays the results.
solver.Solve(decision_builder, [objective, collector])
if collector.SolutionCount() > 0:
best_solution = collector.SolutionCount() - 1
print("Objective value:", collector.ObjectiveValue(best_solution))
print()
print('x= ', collector.Value(best_solution, x))
print('y= ', collector.Value(best_solution, y))
The results of the search, shown below, are the same as the results returned by the MIP solver:
$ python my_projects/int_prog.py Optimal objective value = 23 x = 3 y = 2
If you want to see all solutions returned by the solver, rather than just the optimal one,
you can use the AllSolutionCollector as follows:
collector = solver.AllSolutionCollector()
collector.Add(x)
collector.Add(y)
collector.AddObjective(obj_expr)
solver.Solve(decision_builder, [objective, collector])
solution_count = collector.SolutionCount()
print("Number of solutions: ", solution_count)
print()
for i in range(solution_count):
print("x = ", collector.Value(i, x))
print("y = ", collector.Value(i, x))
print("Objective value: ", collector.Value(i, obj_expr))
print()
Here's the output of the modified code:
Number of solutions: 6 x = 0 y = 0 Objective value: 0 x = 0 y = 1 Objective value: 10 x = 0 y = 2 Objective value: 20 x = 1 y = 2 Objective value: 21 x = 2 y = 2 Objective value: 22 x = 3 y = 2 Objective value: 23
For maximization problems,
the solver only returns a solution when its objective value is greater than that of the
previous one. The variables x and y are returned in increasing lexicographic order because
of the setting ASSIGN_MIN_VALUE in the decision builder. If you change this to
ASSIGN_MAX_VALUE, the solver will begin at the optimal point x = 3,
y = 2, and not return any further solutions.
Finally, see what happens if you remove the objective from the call to the solver.
solver.Solve(decision_builder, collector)Now this is no longer a maximization problem, so the solver returns all feasible solutions:
Number of solutions: 12 x = 0 y = 0 Objective value: 0 x = 0 y = 1 Objective value: 10 x = 0 y = 2 Objective value: 20 x = 1 y = 0 Objective value: 1 x = 1 y = 1 Objective value: 11 x = 1 y = 2 Objective value: 21 x = 2 y = 0 Objective value: 2 x = 2 y = 1 Objective value: 12 x = 2 y = 2 Objective value: 22 x = 3 y = 0 Objective value: 3 x = 3 y = 1 Objective value: 13 x = 3 y = 2 Objective value: 23
The entire program
The entire program is shown below.
from __future__ import print_function
from ortools.constraint_solver import pywrapcp
from ortools.constraint_solver import solver_parameters_pb2
def main():
# Instantiate a CP solver.
parameters = pywrapcp.Solver.DefaultSolverParameters()
solver = pywrapcp.Solver("simple_CP", parameters)
# x and y are integer non-negative variables.
x = solver.IntVar(0, 17, 'x')
y = solver.IntVar(0, 17, 'y')
solver.Add(2*x + 14*y <= 35)
solver.Add(2*x <= 7)
obj_expr = solver.IntVar(0, 1000, "obj_expr")
solver.Add(obj_expr == x + 10*y)
objective = solver.Maximize(obj_expr, 1)
decision_builder = solver.Phase([x, y],
solver.CHOOSE_FIRST_UNBOUND,
solver.ASSIGN_MIN_VALUE)
# Create a solution collector.
collector = solver.LastSolutionCollector()
# Add the decision variables.
collector.Add(x)
collector.Add(y)
# Add the objective.
collector.AddObjective(obj_expr)
solver.Solve(decision_builder, [objective, collector])
if collector.SolutionCount() > 0:
best_solution = collector.SolutionCount() - 1
print("Objective value:", collector.ObjectiveValue(best_solution))
print()
print('x= ', collector.Value(best_solution, x))
print('y= ', collector.Value(best_solution, y))
if __name__ == '__main__':
main()