Initializing help system before first use

# Combining CP and MIP

Topics covered in this section:

- Example: Machine assignment and sequencing
- Model formulation
  - MIP model
  - CP model
- Implementation
- Results
- Parallel solving of CP subproblems

Application problems often combine different subproblems that are solved better with one or another solver, making the complete problem difficult or unmanageable for a single solver. A typical example is production planning and scheduling applications. The long-term (planning) aspects are usually more easily handled by LP solvers whereas the short-term (scheduling) subproblems are better suited for CP solvers. Solving these two parts completely independent of each other may lead to infeasible scheduling subproblems or plans that do not correspond to the reality of production. A possible solution to this dilemma is to iteratively solve LP planning problems and CP scheduling problems, until a feasible schedule for the planned quantities is obtained.

Another method of combined MIP-CP problem solving that provides a tighter integration of the two techniques consists of solving CP subproblems for generating cuts at the nodes of a MIP Branch-and-Bound search. This technique has already been applied successfully to several large-scale planning and scheduling applications by PSA and BASF (see for instance the description of hybrid MIP-CP algorithms implemented with Mosel in [BP03] and [Sad04]). This type of combination is more technical than sequential CP and LP/MIP solving since it requires the developer of the algorithm to interact with the MIP search at every node.

Cut generation algorithms can be implemented with the help of the Xpress Optimizer callbacks (see the `Mosel Language Reference Manual' for the definition of callbacks with Mosel and the `Xpress Optimizer Reference Manual' for an explanation of the Optimizer callback functions).

The original description of the example in this section was published in [JG99]. A prototype implementation was developed by N. Pisaruk in the context of the EU-project LISCOS.

## Example: Machine assignment and sequencing

We need to produce 12 products on a set of three machines. Each machine may produce all of the products but processing times and costs vary (Table Machine-dependent production costs and durations). Furthermore, for every product we are given its release and due dates (Table Release dates and due dates of products). We wish to determine a production plan for all products that minimizes the total production cost.

Table 2: Machine-dependent production costs and durations
| Production costs | Durations
--- | --- | ---
Prod. \ Mach. |  | 1 | 2 | 3 |  | 1 | 2 | 3
1 |  | 12 | 6 | 7 |  | 10 | 14 | 13
2 |  | 13 | 6 | 10 |  | 7 | 9 | 8
3 |  | 10 | 4 | 6 |  | 11 | 17 | 15
4 |  | 8 | 4 | 5 |  | 6 | 9 | 12
5 |  | 12 | 6 | 7 |  | 4 | 6 | 10
6 |  | 10 | 5 | 6 |  | 2 | 3 | 4
7 |  | 7 | 4 | 5 |  | 10 | 15 | 16
8 |  | 9 | 5 | 5 |  | 8 | 11 | 12
9 |  | 10 | 5 | 7 |  | 10 | 14 | 13
10 |  | 8 | 4 | 5 |  | 8 | 11 | 14
11 |  | 15 | 8 | 9 |  | 9 | 12 | 16
12 |  | 13 | 7 | 7 |  | 3 | 5 | 6

Table 3: Release dates and due dates of products
Product | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12
--- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | --- | ---
Release | 2 | 4 | 5 | 7 | 9 | 0 | 3 | 6 | 11 | 2 | 3 | 4
Due date | 32 | 33 | 36 | 37 | 39 | 34 | 30 | 26 | 36 | 38 | 31 | 22

## Model formulation

We are going to represent this problem by two subproblems: the machine assignment problem and the sequencing of operations on machines. The former is implemented by a MIP model; the latter is formulated as a CP (single machine) problem.

#### MIP model

Let COSTpm denote the production cost and DURpm the processing time of product p (p ∈ PRODS, the set of products) on machine m (m ∈ MACH, the set of machines).

To formulate the machine assignment problem we introduce binary variables usepm that take the value 1 if product p is produced on machine m and zero otherwise. The objective function is then given as

minimize

∑
p ∈ PRODS

∑
m ∈ MACH

COSTpm · usepm

The assignment constraints expressing that each order needs exactly one machine for processing it are defined as follows:

∀ p ∈ PRODS:

∑
m ∈ MACH

usepm = 1

In addition to these constraints that already fully state the problem we may define some additional constraints to strengthen the LP relaxation. All production takes place between the earliest release date and the latest due date. If we denote the length of this interval by MAX_LOAD, we may formulate the following valid inequalities expressing that the total processing time of products assigned to a machine cannot exceed MAX_LOAD:

∀ m ∈ MACH:

∑
p ∈ PRODS

DURpm · usepm ≤ MAX_LOAD

#### CP model

Once the set of operations assigned to machine m, ProdMachm (ProdMachm ⊆ PRODS), is known, we obtain the following sequencing problem for this machine:

∀ p ∈ ProdMachm: startp ∈ {RELp, ..., DUEp-DURpm }
∀ p, q ∈ ProdMachm, p < q: startp + DURpm ≤ startq ∨ startq + DURqm ≤ startp

## Implementation

We are using the following algorithm for modeling and solving this problem:

Define the MIP machine assignment problem. Define the operations of the CP model. Start the MIP Branch-and-Bound search. At every node of the MIP search: while function generate_cuts adds at least one cut re-solve the LP-relaxation Function generate_cuts for all machines m call generate_cut_machine(m) if at least one cut has been generated Display statistics Function generate_cut_machine(m) Collect all operations assigned to machine m if more than one operation assigned to m Solve the CP sequencing problem for m if sequencing succeeds Save the solution otherwise Add an infeasibility cut for machine m to the MIP

The implementation of this model is split into two Mosel models: the first, sched_main.mos, contains the main MIP problem and the definition of the cut generation algorithm. The second model, sched_sub.mos, implements the CP single machine sequencing model.

The first part of the main model sets up the data arrays, compiles and loads the CP submodel, calls subroutines for the model definition and problem solving, and finally produces some summary result output. We have defined the filename of the data file as a parameter to be able to change the name of the data file at the execution of the model without having to change the model source. Correspondingly, all data, including the sizes of index sets, are read in from file. At first, we read in only the values of NP and NM. Subsequently, when declaring the sets and arrays that make use of these values, NP and NM are known and the arrays are created as fixed arrays. Otherwise, if their indexing sets are not known, these arrays would automatically be declared as dynamic arrays and for all but arrays of basic types (real, integer, etc.) we have to create their entries explicitly.

    model "Schedule (MIP + CP) main problem"
     uses "mmsystem", "mmxprs", "mmjobs"

     parameters
      DATAFILE = "Data/sched_3_12.dat"
      VERBOSE = 1
     end-parameters

     forward procedure define_MIP_model
     forward procedure setup_cutmanager
     forward function generate_cuts: boolean
     forward procedure print_solution

     declarations
      NP: integer                             ! Number of operations (products)
      NM: integer                             ! Number of machines
     end-declarations

     initializations from DATAFILE
      NP NM
     end-initializations

     declarations
      PRODS = 1..NP                           ! Set of products
      MACH = 1..NM                            ! Set of machines

      REL: array(PRODS) of integer            ! Release dates of orders
      DUE: array(PRODS) of integer            ! Due dates of orders
      MAX_LOAD: integer                       ! max_p DUE(p) - min_p REL(p)
      COST: array(PRODS,MACH) of integer      ! Processing cost of products
      DUR: array(PRODS,MACH) of integer       ! Processing times of products
      starttime: real                         ! Measure program execution time
      ctcut: integer                          ! Counter for cuts
      solstart: array(PRODS) of integer
                                              ! **** MIP model:
      use: array(PRODS,MACH) of mpvar         ! 1 if p uses machine m, otherwise 0
      Cost: linctr                            ! Objective function

      totsolve,totCP: real                    ! Time measurement
      ctrun: integer                          ! Counter of CP runs
      CPmodel: Model                          ! Reference to the CP sequencing model
      ev: Event                               ! Event
      EVENT_SOLVED=2                          ! Event codes sent by submodels
      EVENT_FAILED=3
     end-declarations

     ! Read data from file
     initializations from DATAFILE
      REL DUE COST DUR
     end-initializations

    ! **** Problem definition ****
     define_MIP_model                         ! Definition of the MIP model
     res:=compile("sched_sub.mos")            ! Compile the CP model
     load(CPmodel, "sched_sub.bim")           ! Load the CP model

    ! **** Solution algorithm ****
     starttime:= gettime
     setup_cutmanager                         ! Settings for the MIP search

     totsolve:= 0.0
     initializations to "raw:"
      totsolve as "shmem:solvetime"
      REL as "shmem:REL" DUE as "shmem:DUE"
     end-initializations

     minimize(Cost)                           ! Solve the problem

     writeln("Number of cuts generated: ", ctcut)
     writeln("(", gettime-starttime, "sec) Best solution value: ", getobjval)
     initializations from "raw:"
      totsolve as "shmem:solvetime"
     end-initializations
     writeln("Total CP solve time: ", totsolve)
     writeln("Total CP time: ", totCP)
     writeln("CP runs: ", ctrun)

The MIP model corresponds closely to the mathematical model that we have seen in the previous section.

     procedure define_MIP_model

     ! Objective: total processing cost
      Cost:= sum(p in PRODS, m in MACH) COST(p,m) * use(p,m)

     ! Each order needs exactly one machine for processing
      forall(p in PRODS) sum(m in MACH) use(p,m) = 1

     ! Valid inequalities for strengthening the LP relaxation
      MAX_LOAD:= max(p in PRODS) DUE(p) - min(p in PRODS) REL(p)
      forall(m in MACH) sum(p in PRODS) DUR(p,m) * use(p,m) <= MAX_LOAD

      forall(p in PRODS, m in MACH) use(p,m) is_binary

     end-procedure

The cut generation callback function generate_cuts is called at least once per MIP node. For every machine, it checks whether the assigned operations can be scheduled or whether an infeasibility cut needs to be added. If any cuts have been added, the LP relaxation needs to be re-solved and the cut generation function will be called again, until no more cuts are added.
The function generate_cut_machine first collects all tasks that have been assigned to the given machine m into the set ProdMach by calling the procedure products_on_machine. If there are still unassigned tasks the returned set is empty, otherwise, if the set has more than one element it tries to solve the sequencing subproblem (function solve_CP_problem). If this problem cannot be solved, then the function adds a cut to the MIP problem that makes the current assignment of operations to this machine infeasible.

     procedure products_on_machine(m: integer, ProdMach: set of integer)

      forall(p in PRODS) do
       val:=getsol(use(p,m))
       if (val > 0 and val < 1) then
       ProdMach:={}
        break
       elif val>0.5 then
        ProdMach+={p}
       end-if
      end-do

     end-procedure

    !-----------------------------------------------------------------
    ! Generate a cut for machine m if the sequencing subproblem is infeasible
     function generate_cut_machine(m: integer): boolean
      declarations
       ProdMach: set of integer
      end-declarations

     ! Collect the operations assigned to machine m
      products_on_machine(m, ProdMach)

     ! Solve the sequencing problem (CP model): if solved, save the solution,
     ! otherwise add an infeasibility cut to the MIP problem
      size:= getsize(ProdMach)
      returned:= false
      if (size>1) then
       if not solve_CP_problem(m, ProdMach, 1) then
        Cut:= sum(p in ProdMach) use(p,m) - (size-1)
        if VERBOSE > 2 then
         writeln(m,": ", ProdMach, " <= ", size-1)
        end-if
        addcut(1, CT_LEQ, Cut)
        returned:= true
       end-if
      end-if

     end-function

    !-----------------------------------------------------------------
    ! Cut generation callback function
     function generate_cuts: boolean
      returned:=false; ctcutold:=ctcut

      forall(m in MACH) do
       if generate_cut_machine(m) then
        ctcut+=1
       end-if
      end-do
      if ctcut-ctcutold>0 and VERBOSE>1 then
       writeln("Node ", getparam("XPRS_NODES"), ": ", ctcut-ctcutold,
               " cut(s) added")
      end-if

     end-function

The solving of the CP model is started from the function solve_CP_problem that writes out the necessary data to shared memory and starts the execution of the submodel contained in the file sched_sub.mos.

     function solve_CP_problem(m: integer, ProdMach: set of integer,
                               mode: integer): boolean
      declarations
       DURm: dynamic array(range) of integer
       sol: dynamic array(range) of integer
       solvetime: real
      end-declarations

     ! Data for CP model
      forall(p in ProdMach) DURm(p):= DUR(p,m)
      initializations to "raw:"
       ProdMach as "shmem:ProdMach"
       DURm as "shmem:DURm"
      end-initializations

     ! Solve the problem and retrieve the solution if it is feasible
      startsolve:= gettime
      returned:= false
      if(getstatus(CPmodel)=RT_RUNNING) then
        fflush
        writeln("CPmodel is running")
        fflush
        exit(1)
      end-if

      ctrun+=1
      run(CPmodel, "NP=" + NP + ",VERBOSE=" + VERBOSE + ",MODE=" + mode)
      wait                                  ! Wait for a message from the submodel
      ev:= getnextevent                     ! Retrieve the event
      if getclass(ev)=EVENT_SOLVED then
       returned:= true
       if mode = 2 then
        initializations from "raw:"
         sol as "shmem:solstart"
        end-initializations
        forall(p in ProdMach) solstart(p):=sol(p)
       end-if
      elif getclass(ev)<>EVENT_FAILED then
       writeln("Problem with Kalis")
       exit(2)
      end-if
      wait
      dropnextevent                         ! Ignore "submodel finished" event
      totCP+= (gettime-startsolve)
     end-function

We complete the MIP model with settings for the cut manager and the definition of the integer solution callback. The Mosel comparison tolerance is set to a slightly larger value than the tolerance applied by Xpress Optimizer. It is important to switch the LP presolve off since we interfere with the matrix during the execution of the algorithm (alternatively, it is possible to fine-tune presolve to use only non-destructive algorithms). Sufficiently large space for cuts and cut coefficients should be reserved in the matrix. We also enable output printing by the Optimizer and choose among different MIP log frequencies (depending on model parameter VERBOSE.

     procedure setup_cutmanager
      setparam("XPRS_CUTSTRATEGY", 0)           ! Disable automatic cuts
      feastol:= getparam("XPRS_FEASTOL")        ! Get Optimizer zero tolerance
      setparam("zerotol", feastol * 10)         ! Set comparison tolerance of Mosel
      setparam("XPRS_PRESOLVE", 0)              ! Disable presolve
      setparam("XPRS_MIPPRESOLVE", 0)           ! Disable MIP presolve
      command("KEEPARTIFICIALS=0")              ! No global red. cost fixing
      setparam("XPRS_SBBEST", 0)                ! Turn strong branching off
      setparam("XPRS_HEUREMPHASIS", 0)          ! Disable MIP heuristics
      setparam("XPRS_EXTRAROWS", 10000)         ! Reserve space for cuts
      setparam("XPRS_EXTRAELEMS", NP*30000)     ! ... and cut coefficients
      setcallback(XPRS_CB_OPTNODE, ->generate_cuts) ! Define the optnode callback
      setcallback(XPRS_CB_INTSOL, ->print_solution) ! Define the integer sol. cb.
      setparam("XPRS_COLORDER", 2)
      case VERBOSE of
      1: do
          setparam("XPRS_VERBOSE", true)
          setparam("XPRS_MIPLOG", -200)
         end-do
      2: do
          setparam("XPRS_VERBOSE", true)
          setparam("XPRS_MIPLOG", -100)
         end-do
      3: do                                     ! Detailed MIP output
          setparam("XPRS_VERBOSE", true)
          setparam("XPRS_MIPLOG", 3)
         end-do
      end-case

     end-procedure

The definition of the integer solution callback is, in parts, similar to the function generate_cut_machine. To obtain a detailed solution output we need to re-solve all CP subproblems, this time with run MODE two, meaning that the CP model writes its solution information to shared memory.

     procedure print_solution
      declarations
       ProdMach: set of integer
      end-declarations

      writeln("(",gettime-starttime, "sec) Solution ",
              getparam("XPRS_MIPSOLS"), ": Cost: ", getsol(Cost))

      if VERBOSE > 1 then
       forall(p in PRODS) do
        forall(m in MACH) write(getsol(use(p,m))," ")
        writeln
       end-do
      end-if

      if VERBOSE > 0 then
       forall(m in MACH) do
        ProdMach:= {}

      ! Collect the operations assigned to machine m
        products_on_machine(m, ProdMach)

        Size:= getsize(ProdMach)
        if Size > 1 then
       ! (Re)solve the CP sequencing problem
         if not solve_CP_problem(m, ProdMach, 2) then
          writeln("Something wrong here: ", m, ProdMach)
         end-if
        elif Size=1 then
         elem:=min(p in ProdMach) p
         solstart(elem):=REL(elem)
        end-if
       end-do

     ! Print out the result
       forall(p in PRODS) do
        msol:=sum(m in MACH) m*getsol(use(p,m))
        writeln(p, " -> ", msol,": ", strfmt(solstart(p),2), " - ",
                strfmt(DUR(p,round(msol))+solstart(p),2), "  [",
                REL(p), ", ", DUE(p), "]")
       end-do
       writeln("Time: ", gettime - starttime, "sec")
       writeln
       fflush
      end-if
     end-procedure

The following code listing shows the complete CP submodel. At every execution, the set of tasks assigned to one machine and the corresponding durations are read from shared memory. The disjunctions between pairs of tasks are posted explicitly to be able to stop the addition of constraints if an infeasibility is detected during the definition of the problem. The search stops at the first feasible solution. If a solution was found, it is passed back to the main model if the model parameter MODE has the value two. In every case, after termination of the CP search the submodel sends a solution status event back to the main model.

    model "Schedule (MIP + CP) CP subproblem"
     uses "kalis", "mmjobs" , "mmsystem"

     parameters
      VERBOSE = 1
      NP = 12                                 ! Number of products
      MODE = 1                                ! 1 - decide feasibility
                                              ! 2 - return complete solution
     end-parameters

     startsolve:= gettime

     declarations
      PRODS = 1..NP                           ! Set of products
      ProdMach: set of integer
     end-declarations

     initializations from "raw:"
      ProdMach as "shmem:ProdMach"
     end-initializations

     finalize(ProdMach)

     declarations
      REL: array(PRODS) of integer            ! Release dates of orders
      DUE: array(PRODS) of integer            ! Due dates of orders
      DURm: array(ProdMach) of integer        ! Processing times on machine m
      solstart: array(ProdMach) of integer    ! Solution values for start times

      start: array(ProdMach) of cpvar         ! Start times of tasks
      Disj: array(range) of cpctr             ! Disjunctive constraints
      Strategy: array(range) of cpbranching   ! Enumeration strategy
      EVENT_SOLVED=2                          ! Event codes sent by submodels
      EVENT_FAILED=3
      solvetime: real
     end-declarations

     initializations from "raw:"
      DURm as "shmem:DURm" REL as "shmem:REL" DUE as "shmem:DUE"
     end-initializations

    ! Bounds on start times
     forall(p in ProdMach) setdomain(start(p), REL(p), DUE(p)-DURm(p))

    ! Disjunctive constraint
     ct:= 1
     forall(p,q in ProdMach| p<q) do
      Disj(ct):= start(p) + DURm(p) <= start(q) or start(q) + DURm(q) <= start(p)
      ct+= 1
     end-do

    ! Post disjunctions to the solver
     nDisj:= ct; j:=1; res:= true
     while (res and j<nDisj) do
      res:= cp_post(Disj(j))
      j+=1
     end-do

    ! Solve the problem
     if res then
      Strategy(1):= settle_disjunction(Disj)
      Strategy(2):= assign_and_forbid(KALIS_SMALLEST_DOMAIN, KALIS_MIN_TO_MAX,
                                     start)
      cp_set_branching(Strategy)
      res:= cp_find_next_sol
     end-if

    ! Pass solution to main problem
     if res then
      forall(p in ProdMach) solstart(p):= getsol(start(p))
      if MODE=2 then
       initializations to "raw:"
        solstart as "shmem:solstart"
       end-initializations
      end-if
      send(EVENT_SOLVED,0)
     else
      send(EVENT_FAILED,0)
     end-if

    ! Update total running time measurement
     initializations from "raw:"
      solvetime as "shmem:solvetime"
     end-initializations
     solvetime+= gettime-startsolve
     initializations to "raw:"
      solvetime as "shmem:solvetime"
     end-initializations

    end-model

## Results

The best solution produced for the data set sched_3_12 is the following :

    Cost: 92
    1 -> 3:  2 - 15  [2, 32]
    2 -> 3: 15 - 23  [4, 33]
    3 -> 2: 15 - 32  [5, 36]
    4 -> 1: 24 - 30  [7, 37]
    5 -> 2: 32 - 38  [9, 39]
    6 -> 2:  0 -  3  [0, 34]
    7 -> 1:  3 - 13  [3, 30]
    8 -> 1: 16 - 24  [6, 26]
    9 -> 3: 23 - 36  [11, 36]
    10 -> 1: 30 - 38  [2, 38]
    11 -> 2:  3 - 15  [3, 31]
    12 -> 1: 13 - 16  [4, 22]

A total of 1604 cuts are added to the MIP problem by 2691 CP model runs and the Branch-and-Bound search explores 12295 nodes. Optimality is proven within a few seconds on a Pentium IV PC.

It is possible to implement this problem entirely either with Xpress Optimizer or with Xpress Kalis. However, already for this three machines – 12 jobs instance the problem is extremely hard for either technique on its own. With CP it is difficult to prove optimality and with MIP the formulation of the disjunctions makes the definition of a large number of binary variables necessary (roughly in the order of number_of_machines · number_of_products2) which makes the problem impracticable to deal with.

## Parallel solving of CP subproblems

Instead of solving the CP single-machine subproblems at every MIP node sequentially, we can modify our Mosel models to solve the subproblems in parallel—especially when working on a multiprocessor machine this may speed up the cut generation process and hence shorten the total run time. We modify the algorithm of Section Implementation as follows:

Define the MIP machine assignment problem. Define the operations of the CP model. Start the MIP Branch-and-Bound search. At every node of the MIP search: while function generate_cuts adds at least one cut re-solve the LP-relaxation Function generate_cuts Collect all machines that are fully assigned into set ToSolve for all machines m ∈ ToSolve call start_CP_model(m) Wait for the solution status messages from all submodels if submodel m is infeasible Add an infeasibility cut for machine m to the MIP if at least one cut has been generated Display statistics Procedure start_CP_model(m) Collect all operations assigned to machine m Write data for this machine to memory Start the submodel execution

The modified version of the function generate_cuts looks as follows. For the full example code the reader is referred to the set of User Guide examples provided with the Xpress Kalis distribution (files sched_mainp.mos and sched_subp.mos).

     procedure products_on_machine(m: integer)

      NumOp(m):=0
      forall(p in PRODS) do
       val:=getsol(use(p,m))
       if (! not isintegral(use(p,m)) !) (val > 0 and val < 1) then
        NumOp(m):=0
        break
       elif val>0.5 then
        NumOp(m)+=1
        OpMach(m,NumOp(m)):= p
       end-if
      end-do

     end-procedure

    !-----------------------------------------------------------------
    ! Add an infeasibility cut for machine m to the MIP problem
     procedure add_cut_machine(m: integer)

      Cut:= sum(p in 1..NumOp(m)) use(OpMach(m,p),m) - (NumOp(m)-1)
      if VERBOSE > 1 then
       write(m,": ")
       forall(p in 1..NumOp(m)) write(OpMach(m,p), " ")
       writeln(" <= ", NumOp(m)-1)
      end-if
      addcut(1, CT_LEQ, Cut)

     end-procedure

The implementation of the CP submodels remains largely unchanged, with the exception of the labels employed for passing data via shared memory: we append the machine index to every data item to be able to distinguish between the data used by the different subproblems running in parallel.

For the data set sched_3_12.dat we have observed only a few percent decrease of the total running time on a dual processor machine using the parallel implementation: in many nodes only a single CP subproblem is solved and if there are several subproblems to be solved their execution may be of quite different length. For instances with a larger number of machines the parallelization is likely to show more effect.

Parent Topic Hybrid MIP/CP Solving

© 2001-2026 Fair Isaac Corporation. All rights reserved. This documentation is the property of Fair Isaac Corporation (“FICO”). Receipt or possession of this documentation does not convey rights to disclose, reproduce, make derivative works, use, or allow others to use it except solely for internal evaluation purposes to determine whether to purchase a license to the software described in this documentation, or as otherwise set forth in a written software license agreement between you and FICO (or a FICO affiliate). Use of this documentation and the software described in it must conform strictly to the foregoing permitted uses, and no other use is permitted.