CP Approach to Integer Optimization

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 faster than a MIP solver.

The following sections show how to solve an integer optimization problem using the CP-SAT solver and the original CP solver. The problem is shown below:

Maximize 2x + 2y + 3z subject to the following constraints:
2x + 7 y + 3z50
3x - 5y + 7z45
5x + 2y - 6z37
x, y, z0
x, y, z integers

Solution using the CP-SAT solver

The following sections present a Python program that solves this problem using the CP-SAT 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.

Declare the model

The following code declares the model for the problem.


from ortools.sat.python import cp_model


def main():
  model = cp_model.CpModel()

Create the variables

The following code creates the variables for the problem.

  x = model.NewIntVar(0, cp_model.INT32_MAX, 'x')
  y = model.NewIntVar(0, cp_model.INT32_MAX, 'y')
  z = model.NewIntVar(0, cp_model.INT32_MAX, 'z')

cp_model.INT32_MAX is the largest integer that can be represented on a 32-bit computer. So the command solver.NewIntVar(0, cp_model.INT32_MAX, 'x') creates an integer variable whose lower bound is 0, and which essentially has no upper bound.

Define the constraints

The following code defines the three linear constraints for the problem:

  model.Add(2*x + 7*y + 3*z <= 50)
  model.Add(3*x - 5*y + 7*z <= 45)
  model.Add(5*x + 2*y - 6*z <= 37)

Define the objective function

The following code defines the objective function for the problem and declares it to be a maximization problem:

  model.Maximize(2*x + 2*y + 3*z)

Call the solver

The following code calls the solver and displays the results.

  solver = cp_model.CpSolver()
  status = solver.Solve(model)

  if status == cp_model.OPTIMAL:
    print('Maximum of objective function: %i' % solver.ObjectiveValue())
    print()
    print('x value: ', solver.Value(x))
    print('y value: ', solver.Value(y))
    print('z value: ', solver.Value(z))

The output is shown below:

Maximum of objective function: 35

x value:  7
y value:  3
z value:  5

The entire program

The entire program is shown below.

from __future__ import print_function

from ortools.sat.python import cp_model


def main():
  model = cp_model.CpModel()
  # Set variable upper bound to max. of the upper bounds on the right hand
  # side of the of the constraints below. This does not affect solutions.
  x = model.NewIntVar(0, cp_model.INT32_MAX, 'x')
  y = model.NewIntVar(0, cp_model.INT32_MAX, 'y')
  z = model.NewIntVar(0, cp_model.INT32_MAX, 'z')

  model.Add(2*x + 7*y + 3*z <= 50)
  model.Add(3*x - 5*y + 7*z <= 45)
  model.Add(5*x + 2*y - 6*z <= 37)

  model.Maximize(2*x + 2*y + 3*z)

  solver = cp_model.CpSolver()
  status = solver.Solve(model)

  if status == cp_model.OPTIMAL:
    print('Maximum of objective function: %i' % solver.ObjectiveValue())
    print()
    print('x value: ', solver.Value(x))
    print('y value: ', solver.Value(y))
    print('z value: ', solver.Value(z))


if __name__ == '__main__':
  main()

Solution using the original CP solver

The following sections present a Python program that solves the integer optimization problem using the original CP solver. (However, we recommend the faster CP-SAT solver.)

Declare the solver

the following code declares the solver.

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)

Create the variables

The following code creates the variables.

  var_upper_bound = max(50, 45, 37)
  x = solver.IntVar(0, var_upper_bound, 'x')
  y = solver.IntVar(0, var_upper_bound, 'y')
  z = solver.IntVar(0, var_upper_bound, 'z')

The command solver.NewIntVar(0, var_upper_bound, 'x') creates an integer variable whose values are between 0 and var_upper_bound. Although the problem doesn't impose upper bounds on the individual variables, we still need to supply an upper bound to the solver. So we just choose a value that is large enough that it doesn't change the solutions. Setting the value to the maximum of the constraint upper bounds (shown below) is sufficient.

Define the constraints

The following code defines the constraints for the problem.

  solver.Add(2*x + 7*y + 3*z <= 50)
  solver.Add(3*x - 5*y + 7*z <= 45)
  solver.Add(5*x + 2*y - 6*z <= 37)

Define the objective function

The following code defines the objective for the problem and declares it to be a maximization problem:

  objective = solver.Maximize(2*x + 2*y + 3*z, 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.

  decision_builder = solver.Phase([x, y, z],
                                  solver.CHOOSE_FIRST_UNBOUND,
                                  solver.ASSIGN_MAX_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 x, y and z.
  • 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 are CHOOSE_LOWEST_MIN and CHOOSE_RANDOM.
  • The third parameter specifies how the solver chooses the next value of the current search variable to check. Here, we've chosen ASSIGN_MAX_VALUE, which means the solver chooses the largest value that hasn't already been searched. We could also have chosen ASSIGN_MIN_VALUE or ASSIGN_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.

  collector = solver.LastSolutionCollector()
  # Add the decision variables.
  collector.Add(x)
  collector.Add(y)
  collector.Add(z)
  # Add the objective.
  collector.AddObjective(2*x + 2*y +3*z)

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.

Call the solver

The following code calls the solver and displays the results.

  solver.Solve(decision_builder, [objective, collector])
  if collector.SolutionCount() > 0:
    best_solution = collector.SolutionCount() - 1
    print("Maximum of objective function:", collector.ObjectiveValue(best_solution))
    print()
    print('x= ', collector.Value(best_solution, x))
    print('y= ', collector.Value(best_solution, y))
    print('z= ', collector.Value(best_solution, z))
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])
    Here's the output of the program:
Maximum of objective function: 35

x value:  7
y value:  3
z value:  5

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)
  var_upper_bound = max(50, 45, 37)
  x = solver.IntVar(0, var_upper_bound, 'x')
  y = solver.IntVar(0, var_upper_bound, 'y')
  z = solver.IntVar(0, var_upper_bound, 'z')
  solver.Add(2*x + 7*y + 3*z <= 50)
  solver.Add(3*x - 5*y + 7*z <= 45)
  solver.Add(5*x + 2*y - 6*z <= 37)
  objective = solver.Maximize(2*x + 2*y + 3*z, 1)
  decision_builder = solver.Phase([x, y, z],
                                  solver.CHOOSE_FIRST_UNBOUND,
                                  solver.ASSIGN_MAX_VALUE)
    # Create a solution collector.
  collector = solver.LastSolutionCollector()
  # Add the decision variables.
  collector.Add(x)
  collector.Add(y)
  collector.Add(z)
  # Add the objective.
  collector.AddObjective(2*x + 2*y +3*z)
  solver.Solve(decision_builder, [objective, collector])
  if collector.SolutionCount() > 0:
    best_solution = collector.SolutionCount() - 1
    print("Maximum of objective function:", collector.ObjectiveValue(best_solution))
    print()
    print('x= ', collector.Value(best_solution, x))
    print('y= ', collector.Value(best_solution, y))
    print('z= ', collector.Value(best_solution, z))

if __name__ == '__main__':
  main()

Send feedback about...