Hybridization of CP and MP

Topics covered in this chapter:

This chapter introduces the concept of LP/MIP relaxations for problems formulated as CP models with Xpress Kalis. By means of examples we show how to

• use automated linear programming relaxations,
• work with linear programming relaxations,
• configure CP search to be guided by information from LP/MIP relaxations.

Note: The functionality described in this chapter requires Xpress Optimizer to be installed and licensed in addition to Xpress Kalis.

Linear relaxations

Constraint Programming and Mathematical Programming (here used for Linear and Mixed Integer Programming, LP/MIP) are two different approaches to solving optimization problems that have a number of complementary features. Some problems, or parts of problems, are better solved by CP, others by LP/MIP techniques. Based on this observation it is possible to devise various approaches for hybrid CP-MP problem solving. The FICO Xpress Optimization whitepaper 'Hybrid MIP/CP solving with Xpress Optimizer and Xpress Kalis' describes several decomposition schemes where coordination between MIP and CP models is done on the Mosel level. Xpress Kalis now provides a built-in linear programming relaxation functionality where communication and cooperation take place directly on the solver level, with optional configuration and user interaction from the Mosel CP model.

In the context of Xpress Kalis we use the term linear relaxation for sets of linear (in)equality constraints that are solved by an external LP/MIP solver. In the case of an automatic relaxation this will be a reformulation of (a subset of) the constraints in the CP problem. A user-defined 'linear relaxation', however, may also include additional information, in which case this name choice might be somewhat misleading.

Xpress Kalis can build automatically several linear (or mixed integer linear) relaxations of a CP problem and use Xpress Optimizer to solve them. Xpress Kalis uses a double modeling approach (as in [Hei98]) by exchanging automatically and in a bidirectional way, information such as objective bounds, infeasibility, optimal relaxed solutions and reduced costs between the CP solver and the linear relaxations solver(s) during the search for a solution.

1. Linear relaxations can be generated fully automatically.
2. Alternatively, relaxations can be configured by selecting sets of constraints to form subproblems and by defining their solving frequency.
3. Relaxation problems can also be defined freely by the user; they may even include constraints that do not occur in the original CP problem.

These three ways of defining linear relaxations are presented in the following sections.

Automatic relaxation

When developing a CP model we recommend that you always try the effect of using the automatic relaxation. Indeed, in many of the examples presented in the previous chapters (such as examples a4sugar_ka.mos from Section occurrence: Sugar production or j5tax_ka.mos from Section equiv: Location of income tax offices) we can simplify the model formulation by removing the definition of the search strategy and enabling the automated relaxation by adding the following line before the start of the optimization:

 setparam("kalis_auto_relax",true)

whereby in many cases you will observe that the time to the first or optimal solution is shortened, or that optimality is proven faster than with CP-only search.

Integer knapsack problem with side constraint

The following integer knapsack problem with an additional 'all_different' constraint on its variables will serve us in this and the following sections to show how to work with linear relaxations:

maximize 5·x₁ + 8·x₂ + 4·x₃ + x₄
s.t. x_j ∈ {1,3,8,9,12} for j = 1,2,3,4
3·x₁ + 5·x₂ + 2·x₃ ≤ 50
2·x₁ + x₃ + 5·x₄ ≤ 75
all-different(x₁, x₂, x₃, x₄)

Implementation

The model implementation shown below is fairly straightforward, with the exception of the setting of the control parameter KALIS_AUTO_RELAX.

model "Knapsack with side constraints"
 uses "kalis"

 declarations
  VALUES = {1,3,8,9,12}
  R = 1..4
  x: array(R) of cpvar              ! Decision variables
  benefit: cpvar                    ! The objective to minimize
 end-declarations

! Enable output printing
 setparam("kalis_verbose_level", 1)

! Setting name of variables for pretty printing
 forall(i in R) setname(x(i),"x"+i)
 setname(benefit,"benefit")

! Set initial domains for variables
 forall(i in R) setdomain(x(i), VALUES)

! Knapsack constraints
 3*x(1) + 5*x(2) + 2*x(3) <= 50
 2*x(1) + x(3) + 5*x(4) <= 75

! Additional global constraint
 all_different(union(i in R) {x(i)})

! Objective function
 benefit = 5*x(1) + 8*x(2) + 4*x(3) + x(4)

! Initial propagation
 if not cp_propagate: exit(1)

! Display bounds on objective after constraint propagation
 writeln("Constraint propagation objective ", benefit)


! **** Linear relaxation ****

! Enable automatic linear relaxations
 setparam("kalis_auto_relax", true)

! Solve the problem
 if cp_maximize(benefit):
  cp_show_sol                      ! Output optimal solution to screen

end-model

Results

Enabling the automatic relaxation for this problem reduces the number of nodes spent by the search from 17 to 4. However, more time is spent in every node due to the overhead incurred by the LP solving. As solving times in this small problem are negligible in all formulations time measurements are not really meaningful; in larger problems the effect of using linear relaxations needs to be tested carefully.

Configuring automatic relaxations

This section and the remainder of the chapter are addressed at more expert readers who, besides experience with CP, also have some knowledge of LP/MIP modeling and solving techniques.

Configuration choices

The formulation of a linear relaxation (sub)problem with Xpress Kalis has several components, namely the problem definition (=statement of variables and constraints), the configuration of a linear relaxation solver, and its insertion into the search strategy of the original CP problem. Here is a summary of the different steps and configuration choices.

1. Choosing the constraint type
   

   The following constraints can be relaxed automatically with Xpress Kalis:

   • linear constraints
   • all-different
   • occurrence
   • distribute
   • min/max
   • absolute value
   • distance
   • element
   • cycle
   • logical (implies, or, and, equiv)
   To obtain an automatic relaxation of the whole problem, Xpress Kalis forms the intersection of the relaxations of the constraints of the CP model. This relaxation can be retrieved by a call to the function cp_get_linrelax. An optional argument to this function lets you specify the (set of) constraint types to relax. The relaxation for a specific constraint—as opposed to constraint type—is obtained with get_linrelax.
2. Choosing the relaxation type
   

   Xpress Kalis generally provides several distinct relaxations for a given constraint. Methods cp_get_linrelax and get_linrelax take an integer parameter to choose the type of the relaxation (0 for an LP oriented relaxation and 1 for a MIP oriented relaxation). Besides the automatic relaxations, the user can also freely define his own relaxations of these or other constraints of his CP models (see Section User-defined relaxations below).

3. Defining relaxation solvers
   

   To solve a linear relaxation, Xpress Kalis uses relaxation solvers (model objects of the type cplinrelaxsolver) that define for the relaxation:

   1. an objective variable,
   2. an optimization sense (either minimization or maximization),
   3. a configuration of when and what to do with the relaxation,
   4. a flag indicating whether to solve the relaxation as a pure LP problem or as a MIP.

   The following predefined configurations are implemented within kalis:

   KALIS_TOPNODE_RELAX_SOLVER

   This configuration solves the relaxation only at the top node and provides bound for the objective variable.

   KALIS_TREENODE_RELAX_SOLVER

   This configuration solves the relaxation at each node of the search tree and provides bounds for the objective variable and performs reduced costs propagation.

   KALIS_BILEVEL_RELAX_SOLVER

   This configuration solves the relaxation whenever all the variables (by default all the discrete variables such as cpvar) of a user defined set are instantiated. After the resolution of the relaxation, all the other variables are instantiated to the optimal value of the relaxation. The principal use of the bilevel configuration is to decompose simply and automatically the CP model into a main problem and a subproblem.

   In addition to the predefined schemes, the user can take full control of the relaxation solver by means of callbacks (see the Xpress Kalis Mosel Reference Manual for details).

   Note: a model can contain several linear relaxations solvers, e.g., for solving different subproblems or with different configuration options (such as one solver for the top node and a different relaxation solver for the tree nodes of the CP search).

4. Branching with relaxation information
   

   Certain predefined branching schemes and value and variable selection heuristics in kalis are based on the optimal solution of a relaxation. The information obtained from the relaxations can be used to explore quickly the most promising portions of the search space. For example, a MIP style branching scheme is defined by

   cp_set_branching(assign_and_forbid(KALIS_LARGEST_REDUCED_COST(mysolver),
                                      KALIS_NEAREST_RELAXED_VALUE(mysolver)))

   where the KALIS_LARGEST_REDUCED_COST variable selection heuristic selects the variable with the largest reduced cost in the optimal solution of the relaxation in the specified relaxation solver, and the KALIS_NEAREST_RELAXED_VALUE value selection heuristic selects the value in the domain of the variable that is the nearest (using L1-Norm) to the value of this variable in the optimal solution of the relaxation in the relaxation solver passed in the argument.

Implementation

The model shown below implements a similar relaxation as the automatic version in the previous section. This model has two parts: in the beginning we define a standard CP model for our Knapsack problem and trigger constraint propagation to display the initial bounds on the objective function. The second part shows how to define and configure a linear relaxation solver. At first, the model generates the automatic relaxation (and displays it on screen). The argument value 0 in cp_get_linrelax indicates that we want an LP-oriented relaxation (use 1 for a MIP-oriented relaxation). The relaxation solver defined with get_linrelax_solver is configured to maximize the objective benefit, solving the problem as a MIP, and this only once, at the first node of the CP search. The linear relaxation solver is added to the search process by a call to cp_add_linrelax_solver (there may be several solvers in a search process). The model definition is completed by a `MIP style' branching scheme that branches first on the variables with largest reduced cost, and tests first the values nearest to the optimal solution of the relaxation.

model "Knapsack with side constraints"
 uses "kalis", "mmsystem"

 declarations
  VALUES = {1,3,8,9,12}
  R = 1..4
  x: array(R) of cpvar              ! Decision variables
  benefit: cpvar                    ! The objective to minimize
 end-declarations

! Enable output printing
 setparam("kalis_verbose_level", 1)

! Setting name of variables for pretty printing
 forall(i in R) setname(x(i),"x"+i)
 setname(benefit,"benefit")

! Set initial domains for variables
 forall(i in R) setdomain(x(i), VALUES)

! Knapsack constraints
 3*x(1) + 5*x(2) + 2*x(3) <= 50
 2*x(1) + x(3) + 5*x(4) <= 75

! Additional global constraint
 all_different(union(i in R) {x(i)})

! Objective function
 benefit = 5*x(1) + 8*x(2) + 4*x(3) + x(4)

! Initial propagation
 if not cp_propagate: exit(1)

! Display bounds on objective after constraint propagation
 writeln("Constraint propagation objective ", benefit)


! **** Linear relaxation ****

 declarations
  myrelaxall: cplinrelax
 end-declarations

! Build an automatic 'LP' oriented linear relaxation
 myrelaxall:= cp_get_linrelax(0)

! Output the relaxation to the screen
 cp_show_relax(myrelaxall)

 mysolver:= get_linrelax_solver(myrelaxall, benefit, KALIS_MAXIMIZE,
                KALIS_SOLVE_AS_LP, KALIS_TREENODE_RELAX_SOLVER)

! Define the linear relaxation
 cp_add_linrelax_solver(mysolver)

! Define a 'MIP' style branching scheme using the solution of the relaxation
 cp_set_branching(assign_var(KALIS_LARGEST_REDUCED_COST(mysolver),
                             KALIS_NEAREST_RELAXED_VALUE(mysolver)))

! Solve the problem
 starttime:= gettime
 if cp_maximize(benefit) then
  write(gettime-starttime, "sec. ")
  cp_show_sol                      ! Output optimal solution to screen
 end-if

end-model

Results

The formulation of a CP model for this problem is straightforward. However, solving this type of problems with pure CP search methods tends to be very time consuming due to the poor quality of the bounds on the objective variable derived through constraint propagation. Formulating a MIP model for this problem is a perhaps somewhat more complicated task due to the need to linearize the 'all_different' relation. Here the automatic reformulation provided by Xpress Kalis is certainly helpful. The LP relaxation provides valuable bounds and insight on how to direct the search. Combining both makes the CP search go straight to the optimal solution and optimality is proven immediatly by the relaxation bounds. The configuration options allow you to tune the behavior of the combined solution algorithm.

Try to experiment with other settings:

• LP or MIP style relaxations (first parameter of cp_get_linrelax),
• relax only certain constraints, e.g., by specifying KALIS_LINEAR_CONSTRAINTS in cp_get_linrelax:

  cp_get_linrelax(0, {KALIS_LINEAR_CONSTRAINTS})

• solving frequency (tree nodes or top node): last argument of get_linrelax_solver,
• subproblem type (LP or MIP): 4^th argument of get_linrelax_solver.

Note: Instead of displaying the linear relaxation on screen you may also choose to export it to a file in LP format using procedure export_prob, specifying a linear relaxation solver and a matrix filename (the extension ".lp" will be appended to the name) as arguments.

export_prob(mysolver, "mymatfile")

User-defined relaxations

As an alternative or in addition to the built-in relaxation definitions, linear relaxations (or parts thereof) can be built up from scratch, similarly to the definition of linear constraints in a pure CP model. For the representation of linear relaxations, Xpress Kalis provides the types cpauxvar and cpauxlinexp (for variables and linear expressions that do not occur in the original CP model). The linearization of certain global constraints requires a mapping from integer-valued CP decision variables to sets of 0-1 indicator variables in the LP/MIP formulation, where each indicator variable is associated with a value in the domain of the CP variable. These indicator variables (of type cpauxvar) are obtained with function get_indicator. Besides these special variable types, linear relaxation constraints can also contain variables of type cpvar or cpfloatvar.

Implementation

The implementations below replace the configured definition of a relaxation by a relaxation built up 'manually' by stating the constraints one by one. Just like the previous model version, a relaxation solver is defined with get_linrelax_solver. This linear relaxation solver is then added to the search process by a call to cp_add_linrelax_solver and the model definition is completed by a `MIP style' branching scheme that branches first on the variables with largest reduced cost, and tests first the values nearest to the optimal solution of the relaxation.

This model extract selects one by one the constraints to be added to the relaxation:

 declarations
  myrelax: cplinrelax
 end-declarations

! Build an 'LP' oriented linear relaxation
 myrelax:= get_linrelax(all_different(union(i in R) {x(i)}), 0)
 myrelax += get_linrelax(3*x(1) + 5*x(2) + 2*x(3) <= 50, 0)
 myrelax += get_linrelax(2*x(1) + x(3) + 5*x(4) <= 75, 0)
 myrelax += get_linrelax(benefit - (5*x(1) + 8*x(2) + 4*x(3) + x(4)) = 0, 0)

 mysolver:= get_linrelax_solver(myrelax, benefit, KALIS_MAXIMIZE,
                KALIS_SOLVE_AS_LP, KALIS_TREENODE_RELAX_SOLVER)

! Output the relaxation to the screen (after creation of solver to see all!)
 cp_show_relax(myrelax)

! Define the linear relaxation
 cp_add_linrelax_solver(mysolver)

! Define a 'MIP' style branching scheme using the solution of the relaxation
 cp_set_branching(assign_var(KALIS_LARGEST_REDUCED_COST(mysolver),
                             KALIS_NEAREST_RELAXED_VALUE(mysolver)))

! Solve the problem
 starttime:= gettime
 if cp_maximize(benefit) then
  write(gettime-starttime, "sec. ")
  cp_show_sol                      ! Output optimal solution to screen
 end-if

When adding linear constraints to the relaxation we have the choice to add them directly, or by using the function get_linrelax as in the model extract above. When using this function we can specify through its second argument whether to disregard the integrality condition on any variables the linear constraint.

If we do not want to use any of the built-in relaxation functionality, we can also state the relaxations 'manually'. In the example below the linear relaxation is formulated entirely with the binary indicators for domain values, leaving out the variables x_i that are present in all the automated relaxation versions we have seen so far.

 declarations
  myrelax: cplinrelax
 end-declarations

! Build a user-defined linear relaxation

! Formulation of 'alldifferent':
 forall(val in VALUES)
  myrelax += sum(i in R | contains(x(i),val)) get_indicator(x(i), val) <= 1
 forall(i in R)
  myrelax += sum(val in VALUES | contains(x(i),val)) get_indicator(x(i), val) = 1

! Reformulation of linear constraints
 myrelax +=
  3*sum(val in VALUES | contains(x(1),val)) val*get_indicator(x(1), val) +
  5*sum(val in VALUES | contains(x(2),val)) val*get_indicator(x(2), val) +
  2*sum(val in VALUES | contains(x(3),val)) val*get_indicator(x(3), val) <= 50
 myrelax +=
   2*sum(val in VALUES | contains(x(1),val)) val*get_indicator(x(1), val) +
   sum(val in VALUES | contains(x(3),val)) val*get_indicator(x(3), val) +
   5*sum(val in VALUES | contains(x(4),val)) val*get_indicator(x(4), val) <= 75
 myrelax += benefit -
  (5*sum(val in VALUES | contains(x(1),val)) val*get_indicator(x(1), val) +
   8*sum(val in VALUES | contains(x(2),val)) val*get_indicator(x(2), val) +
   4*sum(val in VALUES | contains(x(3),val)) val*get_indicator(x(3), val) +
   sum(val in VALUES | contains(x(4),val)) val*get_indicator(x(4), val)) = 0

 mysolver:= get_linrelax_solver(myrelax, benefit, KALIS_MAXIMIZE,
                KALIS_SOLVE_AS_LP, KALIS_TREENODE_RELAX_SOLVER)

Results

When comparing the program output obtained from the two formulations above, it appears that the automatic formulation provides better guidance to the CP search (that is, less nodes are required to terminate the enumeration) than the 'manual' formulation that works exclusively with the binary indicator variables. There often are different possibilities for formulating linear relaxations and some experimentation is required to determine which is the best choice.

General recommendations for the formulation of linear relaxations are (a) to be careful to define a relaxation that makes sense and (b) to make sure that the relaxation is connected to the original problem (through variables shared by both formulations) so to obtain the desired benefits from information exchange between the two views of the problem.

Note: kalis also makes accessible additional configuration options for the underlying LP/MIP solver; the interested reader is referred to the documentation of set_linrelax_solver_attribute in the Xpress Kalis Mosel Reference Manual.