Solvers

In this user guide we have explained the basics of working with Mosel, focussing on the formulation of Linear and Mixed Integer Programming problems and related solution techniques. However, the Mosel language is not limited to certain types of Mathematical Programming problems.

An extension to standard LP and MIP is the possibility to work with quadratic objective functions (commonly referred to as Quadratic Programming and Mixed Integer Quadratic Programming). This functionality is provided by the module mmquad and documented in the Mosel Language Reference Manual.

The module mmnl extends the Mosel language with functionality for handling general non-linear constraints. This module (equally documented in the Mosel Language Reference Manual) does not contain any solver on its own. In combination with mmxprs you can use it to formulate and solve QCQP (Quadratically Constrained Quadratic Programming) problems through the QCQP solvers within Xpress Optimizer.

All other solvers of FICO Xpress Optimization (Xpress NonLinear for solving non-linear problems, and Xpress Kalis for Constraint Programming, CP) are provided with separate manuals and their own sets of examples. Please see the Xpress website for an overview of the available products.

With Mosel it is possible to combine several solvers to formulate hybrid solution approaches for solving difficult application problems. The whitepaper Hybrid MIP/CP solving with Xpress Optimizer and Xpress Kalis, available for download from the Xpress website, gives several examples of hybrid solving with LP/MIP and Constraint Programming.

Below follow a few examples for some of the solvers mentioned above.

QCQP solving with Xpress Optimizer

In Section Recursion we have solved a quadratically constrained problem by a recursion algorithm. One might be tempted to try solving this problem with a QCQP solver. Unfortunately, it does not meet the properties required by the definition of QCQP problems: problems solvable by QCQP must not contain equality constraints with quadratic terms (reformulation as two inequalities will not work either) since such constraints do not satisfy the convexity condition. We shall see below how to solve this problem with the general non-linear solver Xpress NonLinear (using SLP).

Let us therefore take a look at a different problem: the task is to distribute a set of points represented by tuples of x-/y-coordinates on a plane minimizing the total squared distance between all pairs of points. For each point i we are given a target location (CX_i,CY_i) and the (square of the) maximum allowable distance to this location, the (squared) radius R_i around this location.

In mathematical terms, we have two decision variables x_i and y_i for the coordinates of every point i. The objective to minimize the total squared distance between all points is expressed by the following sum.

For every point i we have the following quadratic inequality.

The resulting Mosel model (file airport_qp.mos) looks thus.

model "airport"
 uses "mmxprs", "mmnl"

 declarations
  RN: range                      ! Set of airports
  R: array(RN) of real           ! Square of max. distance to given location
  CX,CY: array(RN) of real       ! Target location for each point
  x,y: array(RN) of mpvar        ! x-/y- coordinates
  LimDist: array(RN) of nlctr
 end-declarations

 initialisations from "airport.dat"
  CY  CX  R
 end-initialisations

! Set bounds on variables
 forall(i in RN) do
  -10<=x(i); x(i)<=10
  -10<=y(i); y(i)<=10
 end-do

! Objective: minimise the total squared distance between all points
 TotDist:= sum(i,j in RN | i<j) ((x(i)-x(j))^2+(y(i)-y(j))^2)

! Constraints: all points within given distance of their target location
 forall(i in RN)
  LimDist(i):= (x(i)-CX(i))^2+(y(i)-CY(i))^2 <= R(i)

 setparam("XPRS_verbose", true);
 minimise(TotDist);

 writeln("Solution: ", getobjval);
 forall(i in RN) writeln(i, ": ", getsol(x(i)), ", ", getsol(y(i)))

end-model

A QCQP matrix can be exported to a text file (option "" for MPS or "l" LP format) through the writeprob function of Xpress Optimizer. That is, you need to add the following lines to your model after the problem definition:

 setparam("XPRS_loadnames", true)
 loadprob(TotDist)
 writeprob("airport.mat","l")

A graphical representation of the result with mmsvg, obtained with the following lines of Mosel code, is shown in Figure Result graphic in SVG format.

! Set the size of the displayed graph
 svgsetgraphviewbox(-10,-10,10,10)
 svgsetgraphscale(20)

! Draw the target locations
 svgaddgroup("T", "Target area", SVG_SILVER)
 svgsetstyle(SVG_FILL,SVG_CURRENT)
 forall(i in RN) svgaddcircle(CX(i), CY(i), sqrt(R(i)))

! Draw the solution points
 svgaddgroup("S", "Solution", SVG_BLUE)
 forall(i in RN) svgaddpoint(x(i).sol, y(i).sol)

! Output to file
 svgsave("airport.svg")

! Update the display
 svgrefresh
 svgwaitclose

Figure 21.4: Result graphic in SVG format

Xpress NonLinear

The following example solves the non-linear financial planning problem from Section Recursion with Xpress NonLinear (using the SLP solver). The definition of constraints is straightforward (nonlinear constraints have the type nlctr with mmxnlp). Other functionality contributed by the module mmxnlp that appears in this model are the setinitval subroutine—used for setting start values for the decision variables— and the overloaded minimize subroutine for loading and solving nonlinear problems. Here, we simply need to become feasible and therefore use 0 as argument to minimize. Xpress NonLinear automatically selects a solver among the installed solvers of the Xpress suite (Simplex, Barrier, SLP, or Knitro) depending on the detected problem type. We know that this problem can be solved by recursion and therefore preselect the SLP solver by setting the parameter XNLP_SOLVER.

model "Recursion (NLP)"
 uses "mmxnlp"                       ! Use Xpress NonLinear

 declarations
  NT=6                               ! Time horizon
  QUARTERS=1..NT                     ! Range of time periods
  M,P,V: array(QUARTERS) of real     ! Payments

  interest: array(QUARTERS) of mpvar ! Interest
  net: array(QUARTERS) of mpvar      ! Net
  balance: array(QUARTERS) of mpvar  ! Balance
  rate: mpvar                        ! Interest rate
 end-declarations

 M:: [-1000, 0, 0, 0, 0, 0]
 P:: [206.6, 206.6, 206.6, 206.6, 206.6, 0]
 V:: [-2.95, 0, 0, 0, 0, 0]

 setinitval(rate, 0)                 ! Set initial values for variables
 forall(t in QUARTERS) setinitval(balance(t), 1)

                                     ! net = payments - interest
 forall(t in QUARTERS) net(t) = (M(t)+P(t)+V(t)) - interest(t)	

                                     ! Money balance across periods
 forall(t in QUARTERS) balance(t) = if(t>1, balance(t-1), 0) - net(t)	

                                     ! Interest rate
 forall(t in 2..NT) -(365/92)*interest(t) + balance(t-1) * rate = 0

 interest(1) = 0                     ! Initial interest is zero
 forall(t in QUARTERS) net(t) is_free
 forall(t in 1..NT-1) balance(t) is_free
 balance(NT) = 0                     ! Final balance is zero

 ! setparam("XNLP_VERBOSE",true)     ! Uncomment to see detailed output
 setparam("XNLP_SOLVER", 0)          ! Use the SLP solver

 minimize(0)                         ! Solve the problem (get feasible)

                                     ! Print the solution
 writeln("\nThe interest rate is ", getsol(rate))
 write(strfmt("t",5), strfmt(" ",4))
 forall(t in QUARTERS) write(strfmt(t,5), strfmt(" ",3))
 write("\nBalances ")
 forall(t in QUARTERS) write(strfmt(getsol(balance(t)),8,2))
 write("\nInterest ")
 forall(t in QUARTERS) write(strfmt(getsol(interest(t)),8,2))
 writeln

end-model

The results displayed by this model are exactly the same as those we have obtained from the recursion algorithm in Section Recursion.

Xpress Kalis

Constraint Programming is an entirely different method of representing and solving problems compared to the Mathematical Programming approaches we have employed so far in this manual. Consequently, the Mosel module kalis defines its own set of types that behave differently from their Mathematicl Programming counterparts, including cpvar, cpfloatvar and cpctr for decision variables and constraints, but also cpbranching for search strategies (a standard constituent of CP models) and aggregate modeling objects such as cptask and cpresource.

Below we show the CP implementation of a binpacking example from the book `Applications of optimization with Xpress-MP' (Section 9.4 `Backing up files'). The problem is to save sixteen files of different sizes onto empty disks of the same fixed capacity minimizing the total number of disks that are used.

Our model formulation remains close to the Mixed Integer Programming formulation and yet shows some specifities of Constraint Programming. Two sets of decision variables are used, save_f indicating the choice of disk for file f and use_fd the amount of diskspace used by a file on a particular disk. Whereas the variables save_f simply take integer values in a specified interval, each of the use_fd variables may take only two values, 0 or SIZE_f. For every file / disk combination we establish the logical relation

a constraint that cannot be stated in this form in a Mathematical Programming model. A second, so-called global constraint relation is used to state the maximum relation for calculating the number of disks used:

And finally, the constraint limiting the capacity available per disk is a linear inequality that could occur in the same form in a Mathematical Programming model.

This is the complete Mosel CP model for the binpacking problem.

model "D-4 Bin packing (CP)"
 uses "kalis"

 declarations
  ND: integer                        ! Number of floppy disks
  FILES = 1..16                      ! Set of files
  DISKS: range                       ! Set of disks

  CAP:  integer                      ! Floppy disk size
  SIZE: array(FILES) of integer      ! Size of files to be saved
 end-declarations

 initializations from 'd4backup.dat'
  CAP SIZE
 end-initializations

! Provide a sufficiently large number of disks
 ND:= ceil((sum(f in FILES) SIZE(f))/CAP)
 DISKS:= 1..ND
 finalize(DISKS)

 setparam("kalis_default_lb", 0)

 declarations
  save: array(FILES) of cpvar       ! Disk a file is saved on
  use: array(FILES,DISKS) of cpvar  ! Space used by file on disk
  diskuse: cpvar                    ! Number of disks used
 end-declarations

! Set variable domains
 forall(f in FILES) setdomain(save(f), DISKS)
 forall(f in FILES, d in DISKS) setdomain(use(f,d), {0, SIZE(f)})

! Correspondence between disk choice and space used
 forall(f in FILES, d in DISKS) equiv(save(f)=d, use(f,d)=SIZE(f))

! Limit the number of disks used
 diskuse = maximum(save)

! Capacity limit of disks
 forall(d in DISKS) sum(f in FILES) use(f,d) <= CAP

! Minimize the total number of disks used
 if not cp_minimize(diskuse) then
  writeln("Problem infeasible")
 end-if

! Solution printing
 writeln("Number of disks used: ", getsol(diskuse))
 forall(d in 1..getsol(diskuse)) do
  write(d, ":")
  forall(f in FILES) write( if(getsol(save(f))=d , " "+SIZE(f), ""))
  writeln("  space used: ", getsol(sum(f in FILES) use(f,d)))
 end-do

end-model

This implementation is one of several possible formulations of this problem as a CP model. An alternative and generally more efficient model being the formulation as a cumulative scheduling problem, where the disks are represented by a single resource of discrete capacity, the files to save correspond to tasks of duration 1 with a resource requirement defined by the file size. The objective in this case is to minimize the duration of the schedule (= number of disks used). The interested reader is refered to the Xpress Kalis User Guide for a detailed discussion of this problem.