Mixed Integer Programming
This chapter extends the model developed in Chapter 10 to a Mixed Integer Programming (MIP) problem. It describes how to
- define different types of discrete variables,
- get the MIP solution status and understand the MIP optimization log produced by Xpress-Optimizer.
Chapter 6 shows how to formulate and solve the same example with Mosel and in Chapter 16 the problem is input and solved directly with Xpress-Optimizer.
Extended problem description
The investor is unwilling to have small share holdings. He looks at the following two possibilities to formulate this constraint:
- Limiting the number of different shares taken into the portfolio.
- If a share is bought, at least a minimum amount 10% of the budget is spent on the share.
We are going to deal with these two constraints in two separate models.
MIP model 1: limiting the number of different shares
To be able to count the number of different values we are investing in, we introduce a second set of variables buys in the LP model developed in Chapter 2. These variables are indicator variables or binary variables. A variable buys takes the value 1 if the share s is taken into the portfolio and 0 otherwise.
We introduce the following constraint to limit the total number of assets to a maximum of MAXNUM. It expresses the constraint that at most MAXNUM of the variables buys may take the value 1 at the same time.
s 
SHARES buys 
MAXNUM
We now still need to link the new binary variables buys with the variables fracs, the quantity of every share selected into the portfolio. The relation that we wish to express is `if a share is selected into the portfolio, then it is counted in the total number of values' or `if fracs > 0 then buys = 1'. The following inequality formulates this implication:
s 
SHARES: fracs buys
If, for some s, fracs is non-zero, then buys must be greater than 0 and hence 1. Conversely, if buys is at 0, then fracs is also 0, meaning that no fraction of share s is taken into the portfolio. Notice that these constraints do not prevent the possibility that buys is at 1 and fracs at 0. However, this does not matter in our case, since any solution in which this is the case is also valid with both variables, buys and fracs, at 0.
Implementation with BCL
We extend the LP model developed in Chapter 10 with the new variables and constraints. The fact that the new variables are binary variables (i.e. they only take the values 0 and 1) is expressed through the type XPRB_BV at their creation.
Another common type of discrete variable is an integer variable, that is, a variable that can only take on integer values between given lower and upper bounds. These variables are defined in BCL with the type XPRB_UI. In the following section (MIP model 2) we shall see yet another example of discrete variables, namely semi-continuous variables.
#include <iostream>
#include "xprb_cpp.h"
using namespace std;
using namespace ::dashoptimization;
#define MAXNUM 4 // Max. number of shares to be selected
#define NSHARES 10 // Number of shares
#define NRISK 5 // Number of high-risk shares
#define NNA 4 // Number of North-American shares
double RET[] = {5,17,26,12,8,9,7,6,31,21}; // Estimated return in investment
int RISK[] = {1,2,3,8,9}; // High-risk values among shares
int NA[] = {0,1,2,3}; // Shares issued in N.-America
int main(int argc, char **argv)
{
int s;
XPRBprob p("FolioMIP1"); // Initialize a new problem in BCL
XPRBlinExp Risk,Na,Return,Cap,Num;
XPRBvar frac[NSHARES]; // Fraction of capital used per share
XPRBvar buy[NSHARES]; // 1 if asset is in portfolio, 0 otherwise
// Create the decision variables (including upper bounds for `frac')
for(s=0;s<NSHARES;s++)
{
frac[s] = p.newVar("frac", XPRB_PL, 0, 0.3);
buy[s] = p.newVar("buy", XPRB_BV);
}
// Objective: total return
for(s=0;s<NSHARES;s++) Return += RET[s]*frac[s];
p.setObj(Return); // Set the objective function
// Limit the percentage of high-risk values
for(s=0;s<NRISK;s++) Risk += frac[RISK[s]];
p.newCtr(Risk <= 1.0/3);
// Minimum amount of North-American values
for(s=0;s<NNA;s++) Na += frac[NA[s]];
p.newCtr(Na >= 0.5);
// Spend all the capital
for(s=0;s<NSHARES;s++) Cap += frac[s];
p.newCtr(Cap == 1);
// Limit the total number of assets
for(s=0;s<NSHARES;s++) Num += buy[s];
p.newCtr(Num <= MAXNUM);
// Linking the variables
for(s=0;s<NSHARES;s++) p.newCtr(frac[s] <= buy[s]);
// Solve the problem
p.maxim("g");
// Solution printing
cout << "Total return: " << p.getObjVal() << endl;
for(s=0;s<NSHARES;s++)
cout << s << ": " << frac[s].getSol()*100 << "% (" << buy[s].getSol()
<< ")" << endl;
return 0;
}
Besides the additional variables and constraints, the choice of optimization algorithm needs to be adapted to the problem type: we now wish to solve a MIP problem via Branch-and-Bound, indicated by the flag "g" (for `global search') for method maxim.
Just as with the LP problem in the previous chapter, it is usually helpful to check the solution status before accessing the MIP solution—only if the MIP status is `unfinished (solution found)' or `optimal' will BCL print out a meaningful solution:
char *MIPSTATUS[] = {"not loaded", "not optimized", "LP optimized",
"unfinished (no solution)",
"unfinished (solution found)", "infeasible", "optimal"};
cout << "Problem status: " << MIPSTATUS[p.getMIPStat()] << endl;
Analyzing the solution
As the result of the execution of our program we obtain the following output:
Reading Problem FolioMIP1
Problem Statistics
14 ( 514 spare) rows
20 ( 0 spare) structural columns
49 ( 5056 spare) non-zero elements
Global Statistics
10 entities 0 sets 0 set members
Presolved problem has: 14 rows 20 cols 49 non-zeros
LP relaxation tightened
Its Obj Value S Ninf Nneg Sum Inf Time
0 42.600000 D 12 0 6.166667 0
14 14.066667 D 0 0 .000000 0
Optimal solution found
*** Heuristic solution found: 8.900000 Time: 0 ***
Generating cuts
Its Type BestSoln BestBound Sols Add Del Gap GInf Time
1 K 8.900000 13.300000 1 4 0 33.08% 2 0
2 K 8.900000 13.100000 1 2 1 32.06% 0
Cuts in the matrix : 5
Cut elements in the matrix : 34
Cuts in the cutpool : 6
Cut elements in the cutpool: 38
*** Heuristic solution found: 8.900000 Time: 0 ***
Branch Parent Solution Entity Value/Bound Active GInf Time
*** Solution found **** Time: 0
1 0 13.100000 0 0 0
*** Relative MIP gap less than MIPRELSTOP ***
*** Search completed *** Time: 0 Nodes: 1
Number of integer feasible solutions found is 2
Best integer solution found is 13.100000
Best bound is 13.100000
Uncrunching matrix
Problem status: optimal
Total return: 13.1
0: 20% (1)
1: 0% (0)
2: 30% (1)
3: 0% (0)
4: 20% (1)
5: 30% (1)
6: 0% (0)
7: 0% (0)
8: 0% (0)
9: 0% (0)
At the beginning we see the log of the execution of Xpress-Optimizer: the problem statistics (we now have 14 constraints and 20 variables, out of which 10 are MIP variables, refered to as `entities'), the log of the execution of the LP algorithm, the log of the built-in MIP heuristics (a solution with the value 8.9 has been found), the log of the automated cut generation (a total of 6 cuts of type `K' = knapsack have been generated), and the log of the Branch-and-Bound search. Since this problem is very small, it is solved by the addition of cuts (additional constraints that cut off parts of the LP solution space, but no MIP solution) that tighten the LP formulation in such a way that the solution to the LP relaxation becomes integer feasible. The Branch-and-Bound search therefore stops at the first node.
The output printed by our program tells us that the problem has been solved to optimality (i.e. the MIP search has been completed and at least one integer feasible solution has been found). The maximum return is now lower than in the original LP problem due to the additional constraint. As required, only four different shares are selected to form the portfolio.
MIP model 2: imposing a minimum investment in each share
To formulate the second MIP model, we start again with the LP model from Chapters 2 and 10. The new constraint we wish to formulate is `if a share is bought, at least a minimum amount 10% of the budget is spent on the share.' Instead of simply constraining every variable fracs to take a value between 0 and 0.3, it now must either lie in the interval between 0.1 and 0.3 or take the value 0. This type of variable is known as semi-continuous variable. In the new model, we replace the bounds on the variables fracs by the following constraint:
s SHARES: fracs = 0 or 0.1 fracs 0.3
Implementation with BCL
The following program implements the MIP model 2. The semi-continuous variables are defined by the type XPRB_SC. By default, BCL assumes a continuous limit of 1, se we need to set this value to 0.1 with the method setLim.
A similar type is available for integer variables that take either the value 0 or an integer value between a given limit and their upper bound (so-called semi-continuous integers): XPRB_SI. A third composite type is a partial integer which takes integer values from its lower bound to a given limit value and is continuous beyond this value (marked by XPRB_PI).
#include <iostream>
#include "xprb_cpp.h"
using namespace std;
using namespace ::dashoptimization;
#define NSHARES 10 // Number of shares
#define NRISK 5 // Number of high-risk shares
#define NNA 4 // Number of North-American shares
double RET[] = {5,17,26,12,8,9,7,6,31,21}; // Estimated return in investment
int RISK[] = {1,2,3,8,9}; // High-risk values among shares
int NA[] = {0,1,2,3}; // Shares issued in N.-America
int main(int argc, char **argv)
{
int s;
XPRBprob p("FolioSC"); // Initialize a new problem in BCL
XPRBlinExp Risk,Na,Return,Cap;
XPRBvar frac[NSHARES]; // Fraction of capital used per share
// Create the decision variables
for(s=0;s<NSHARES;s++)
{
frac[s] = p.newVar("frac", XPRB_SC, 0, 0.3);
frac[s].setLim(0.1);
}
// Objective: total return
for(s=0;s<NSHARES;s++) Return += RET[s]*frac[s];
p.setObj(Return); // Set the objective function
// Limit the percentage of high-risk values
for(s=0;s<NRISK;s++) Risk += frac[RISK[s]];
p.newCtr(Risk <= 1.0/3);
// Minimum amount of North-American values
for(s=0;s<NNA;s++) Na += frac[NA[s]];
p.newCtr(Na >= 0.5);
// Spend all the capital
for(s=0;s<NSHARES;s++) Cap += frac[s];
p.newCtr(Cap == 1);
// Solve the problem
p.maxim("g");
// Solution printing
cout << "Total return: " << p.getObjVal() << endl;
for(s=0;s<NSHARES;s++)
cout << s << ": " << frac[s].getSol()*100 << "%" << endl;
return 0;
}
When executing this program we obtain the following output (leaving out the part printed by the Optimizer):
Total return: 14.0333 0: 30% 1: 0% 2: 20% 3: 0% 4: 10% 5: 26.6667% 6: 0% 7: 0% 8: 13.3333% 9: 0%
Now five securities are chosen for the portfolio, each forming at least 10% and at most 30% of the total investment. Due to the additional constraint, the optimal MIP solution value is again lower than the initial LP solution value.
If you have any comments or suggestions about these pages, please send mail to docs@dashoptimization.com.