Further modeling topics
Data input and index sets
BCL requires the user to read data into their own structures or data arrays by using standard C functions for accessing data files. The functions XPRBreadarrline and XPRBreadline read data from data files in the diskdata format (see the documentation of the module mmetc in the Xpress-Mosel Language Reference Manual for details). The first function reads (dense) data tables with all entries of the same type, the second reads tables with items of different types (such as text strings and numbers). In particular, XPRBreadline is well suited to read sparse data tables that are indexed by so-called index sets. Roughly speaking, an index set is a set of items such as text strings that index data tables and other objects in the model in a clearer way than numerical values (for details refer to the Xpress-Mosel Reference Manual).
- Data input from file
FILE *datafile;char name[50];double dval, dvals[5];XPRBreadlinecb(XPRB_FGETS,datafile,200,"T,d",name,&dval);XPRBreadarrlinecb(XPRB_FGETS,datafile,200,"d;",dvals,5);- Create a new index set
- XPRBidx set1;set1 = XPRBnewidxset(prob,"Set1",100);
- Add index to a set
- XPRBaddidxel(set1,"Prob1");
- Accessing index sets
- int size, ind;ind = XPRBgetidxel(set1,"Prod1");name = XPRBgetidxelname(set1,14);name = XPRBgetidxsetname(set1);size = XPRBgetidxsetsize(set1);
Figure 3.1: Data input from file and accessing index sets: creation of sets, addition of elements, retrieving elements, and the index set size.
A new index set is created by calling function XPRBnewidxset. Set elements are added with function XPRBaddidxel. An element of a set can be retrieved either by its name (XPRBgetidxel) or by its order number within the set (using the function XPRBgetidxelname). A data item may be part of several index sets. Function XPRBgetidxsetsize returns the current size (i.e. the number of set elements) of an index set.
The definition of index sets may be nested, that is while reading a data file the user may fill up several index sets at a time. The size of index sets grows automatically as required. The user sets some initial size at the creation of the set, but if less elements are added the size returned by XPRBgetidxsetsize will be smaller than this value and if more elements are added the size is increased accordingly.
Example
Taking the program example from the previous chapter, we now assume that we want to give names to the jobs, such as ABC14, DE45F, GH9IJ99, KLMN789. We further assume that these names, together with the durations, are given in a separate data file, durations.dat:
ABC14, 3 DE45F, 4 GH9IJ99, 2 KLMN789, 2
If data is read with function XPRBreadline, it is possible to use comments (preceded by !) and line continuation signs (&) in the data file. (Note that single strings and numbers may not be written over several lines.) The input function also skips blanks and empty lines. If separator signs other than blanks are used, the value 0 may be omitted in the data file (for instance, a data line 0, 0, 0 could as well be written as , , or, using blanks as separators, 0 0 0). The following is functionally equivalent to the contents of durations.dat:
ABC14, 3 ! product1, duration1
DE45F, & ! this line is continued
4 ! in the next line
GH9IJ99, 2 ! blanks are skipped
! as well as empty lines
KLMN789, 2
Separating the input data from the definition allows the same model to be rerun with different data sets without having to recompile the program code. To accommodate data in this form the model program must be written or edited appropriately. In the following program, a function for data input is added to the code seen in the previous chapter. The space for the decision variable arrays is allocated once the array sizes are known. Notice that we use the job names as the names of the decision variables.
#include <stdio.h>
#include <stdlib.h>
#include "xprb.h"
#define MAXNJ 4 /* Maximum number of jobs */
#define NT 10 /* Time limit */
int NJ=0; /* Number of jobs read in */
double DUR[MAXNJ]; /* Durations of jobs */
XPRBidxset Jobs; /* Job names */
XPRBvar *start; /* Start times of jobs */
XPRBvar **delta; /* Binaries for start times */
XPRBvar z; /* Max. completion time */
XPRBprob prob; /* BCL problem */
void read_data(void)
{
char name[100];
FILE *datafile;
Jobs = XPRBnewidxset(prob,"jobs",MAXNJ);
/* Create a new index set */
datafile=fopen("durations.dat","r");
/* Open data file for read access */
while(NJ<MAXNJ) &&
XPRBreadlinecb(XPRB_FGETS, datafile, 99, "T,d", name, &DUR[NJ]))
{ /* Read in all (non-empty) lines up to the end of the file */
XPRBaddidxel(Jobs,name); /* Add job to the index set */
NJ++;
}
fclose(datafile); /* Close the input file */
printf("Number of jobs read: %d\n", XPRBgetidxsetsize(Jobs));
}
void jobs_model(void)
{
XPRBctr ctr;
int j,t;
/* Create start time variables with bounds */
start = (XPRBvar *)malloc(NJ * sizeof(XPRBvar));
if(start==NULL)
{ printf("Not enough memory for 'start' variables.\n");
exit(0); }
for(j=0;j<NJ;j++)
start[j] = XPRBnewvar(prob,XPRB_PL,"start",0,NT-DUR[j]+1);
z = XPRBnewvar(prob,XPRB_PL,"z",0,NT); /* Makespan var. */
/* Declare binaries for each job */
delta = (XPRBvar **)malloc(NJ * sizeof(XPRBvar*));
if(delta==NULL)
{ printf("Not enough memory for 'delta' variables.\n");
exit(0); }
for(j=0;j<NJ;j++)
{
delta[j] = (XPRBvar *)malloc(NT* sizeof(XPRBvar));
if(delta[j]==NULL
{ printf("Not enough memory for 'delta_j' variables.\n");
exit(0); }
delta[j][t] = XPRBnewvar(XPRB_BV,
XPRBnewname("delta%s_%d",XPRBgetidxelname(Jobs,j),t+1),
0,1);
}
for(j=0;j<NJ;j++) /* Calculate max. completion time */
XPRBnewprec(prob,"Makespan",start[j],DUR[j],z);
/* Precedence relation betw. jobs */
XPRBnewprec(prob,"Prec",start[0],DUR[0],start[2]);
for(j=0;j<NJ;j++) /* Linking start times & binaries */
{
ctr = XPRBnewctr(prob,"Link",XPRB_E);
for(t=0;t<(NT-DUR[j]+1);t++)
XPRBaddterm(ctr,delta[j][t],t+1);
XPRBaddterm(ctr,start[j],-1);
}
for(j=0;j<NJ;j++) /* Unique start time for each job */
{
ctr = XPRBnewctr(prob,"One",XPRB_E);
for(t=0;t<(NT-DUR[j]+1);t++) XPRBaddterm(ctr,delta[j][t],1);
XPRBaddterm(ctr,NULL,1);
}
ctr = XPRBnewctr(prob,"OBJ",XPRB_N);
XPRBaddterm(ctr,z,1);
XPRBsetobj(prob,ctr); /* Set objective function */
jobs_solve(); /* Solve the problem */
free(start);
for(j=0;j<NJ;j++) free(delta[j]);
free(delta);
}
int main(int argc, char **argv)
{
prob=XPRBnewprob("Jobs"); /* Initialization */
read_data(); /* Read data from file */
jobs_model(); /* Define & solve the problem */
return 0;
}
Special Ordered Sets
Basic functions
Special Ordered Sets of type n (n=1,2) are sets of variables of which at most n may be non-zero at an integer feasible solution. Associated with each set member is a real number (weight), establishing an ordering among the members of the set. In SOS of type 2, any positive variables must be adjacent in the sense of this ordering.
- XPRBsos set1, set2;XPRBarrvar s;
- Immediate (ref. constraint)
- XPRBctr c;set1=XPRBnewsosrc(prob,"sA",XPRB_S2,s,c);
- Immediate (coefficients)
- double C[] = {1,2,3,4};set2=XPRBnewsosw(prob,"sB",XPRB_S1,s,C);
- Consecutive definition
- set2=XPRBnewsos(prob,"sB",XPRB_S1);XPRBaddsosarrel(set2,s,C);
- Delete set definition
- XPRBdelsos(set2);
- Accessing sets
- XPRBaddsosel(set2,s[2],4,5);XPRBdelsosel(set1,s[0]);XPRBgetsosname(set1);XPRBgetsostype(set2);
Figure 3.2: Defining and accessing SOS: immediate (single function) by indicating a reference constraint; or consecutive definition by adding coefficients for all
members.
In BCL, Special Ordered Sets may be defined in different ways as illustrated in Figure Defining and accessing SOS. As with arrays and constraints, they may be created either with a call to a single function (see Section Array-based SOS definition), or by adding coefficients consecutively.
In the basic, incremental definition, function XPRBnewsos marks the beginning of the definition of a set. Single members are added by function XPRBaddsosel and arrays by function XPRBaddsosarrel, each time indicating the corresponding coefficients. Single elements, or an entire set definition, can be deleted with functions XPRBdelsosel and XPRBdelsos respectively. BCL also has functions to retrieve the name of a SOS and its type (XPRBgetsosname and XPRBgetsostype). It is also possible to set branching directives for a SOS (function XPRBsetsosdir), including priorities, choice of the preferred branching direction and definition of pseudo costs.
Array-based SOS definition
BCL provides two functions for creating Special Ordered Sets with a single function call: XPRBnewsosrc and XPRBnewsosw. With both functions, a new SOS is created by indicating the type (1 or 2), an array of variables and the corresponding weight coefficients for establishing an ordering among the set elements. With XPRBnewsosrc, these coefficients are taken from the variables' coefficients in the indicated reference constraint. When using function XPRBnewsosw, the user directly provides an array of weight coefficients.
Example
In the previous examples, instead of defining the delta variables as binaries, the problem can also be formulated using SOS of type 1. In this case, the delta variables are defined to be continuous as the SOS1 property and their unit sum ensure that one and only one takes the value one.
XPRBprob prob; /* BCL problem */
XPRBvar delta[NJ][NT]; /* Variables for start times */
XPRBsos set[NJ];
void jobs_model(void)
{
...
for(j=0;j<NJ;j++) /* Declare a variable for each job */
for(t=0;t<NT-DUR[j]+1;t++) /* and for each start time */
delta[j][t] = XPRBnewvar(prob, XPRB_PL,
XPRBnewname("delta%d%d",j+1,t+1), 0, 1);
for(j=0;j<NJ;j++)
{ /* Create a new SOS1 */
set[j] = XPRBnewsos(prob, "sosj", XPRB_S1);
for(t=0;t<NT-D[j]+1;t++) /* Add variables to the SOS */
XPRBaddsosel(set[j], delta[j][t], t+1);
}
}
In order to simplify the definition of the SOS one can use the model formulation with variable arrays presented in the previous chapter. The constraints Link are employed as the reference constraints to determine the weight coefficient for each variable (the constraints need to be stored in an array, Link).
XPRBprob prob; /* BCL problem */
XPRBarrvar delta[NJ]; /* Sets of var.s for start times */
XPRBsos set[NJ];
void jobs_model(void)
{
XPRBctr Link[NJ]; /* "Link" constraints */
...
for(j=0;j<NJ;j++) /* Declare a set of var.s for each job */
delta[j] = XPRBnewarrvar(prob, (NT-(int)DUR[j]+1), XPRB_PL,
XPRBnewname("delta%d",j+1), 0, 1);
for(j=0;j<NJ;j++) /* Linking start times & binaries */
{
Link[j] = XPRBnewsumc(prob,"Link",delta[j],1,XPRB_E,0);
XPRBaddterm(Link[j],start[j],-1);
}
/* Create a SOS1 for each job using constraints "Link" as
reference constraints */
for(j=0;j<NJ;j++)
set[j] = XPRBnewsosrc(prob, "sosj", XPRB_S1, delta[j], Link[j]);
}
Instead of setting directives on the binary variables, we may now define branching directives for the SOS1.
for(j=0;j<NJ;j++) XPRBsetsosdir(set[j],XPRB_DN,0);
/* First branch downwards on sets */
Output and printing
BCL provides printing functions for variables, constraints, Special Ordered Sets, and index sets (XPRBprintvar, XPRBprintarrvar, XPRBprintctr, XPRBprintsos, XPRBprintidxset) as well as the entire model definition (XPRBprintprob). Any program output may be printed with XPRBprintf in a similar way to the C function printf. The output of all functions mentioned above is intercepted by the callback XPRBdefcbmsg if this function has previously been defined by the user.
It is also possible to output the problem to a file in extended LP format or as a matrix in extended MPS format (function XPRBexportprob). Note that unlike standard LP format, the extended LP format supports Special Ordered Sets and non-standard variable types (semi-continuous, semi-integer, or partial integers). Like the standard LP format it requires the sense of the objective function to be defined.
- File output
- XPRBexportprob(prob,XPRB_MPS,"expl2");
- Print model objects
- XPRBvar y;XPRBprintvar(y);
- XPRBarrvar av;XPRBprintarrvar(av);
- XPRBctr c;XPRBprintctr(c);
- XPRBsos s;XPRBprintsos(s);
- XPRBidxset is;XPRBprintidxset(is);
- Print a given problem
- XPRBprintprob(prob);
- Print program output
- XPRBprintf("Print this text");
- Compose a name string
- int i = 3;XPRBnewname("abc%d",i);
Figure 3.3: File output and printing.
Example
We may now augment the last few lines of the model definition (cmodel or cmodel_array) of our example with some output functions. Note that these output functions may be added at any time to print the current problem definition in BCL. The function XPRBprintprob prints the complete BCL problem definition to the standard output. The function XPRBexportprob writes the problem definition in LP format or as a matrix in extended MPS format to the indicated file.
XPRBprintprob(prob); /* Print out the problem definition */
XPRBexportprob(prob,XPRB_MPS,"expl1");
/* Output matrix to MPS file */
Instead of printing the entire problem with function XPRBprintprob, it is also possible to display single variables or constraints as soon as they have been defined. The following modified extract of the model definition may serve as an example.
#include <stdio.h>
#include "xprb.h"
#define NJ 4 /* Number of jobs */
#define NT 10 /* Time limit */
double DUR[] = {3,4,2,2}; /* Durations of jobs */
XPRBvar start[NJ]; /* Start times of jobs */
XPRBprob prob; /* BCL problem */
...
void cmodel(void)
{
XPRBctr ctr;
int j,t;
prob=XPRBnewprob("Jobs"); /* Initialization */
for(j=0;j<NJ;j++) /* Create start time variables */
{
start[j] = XPRBnewvar(prob,XPRB_PL,"start",0,NT);
XPRBprintvar(start[j]);
XPRBprintf(", ");
}
...
/* Precedence relation betw. jobs */
ctr = XPRBnewprec(prob,"Prec",start[0],DUR[0],start[2]);
XPRBprintctr(ctr);
...
}
Quadratic Programming with BCL
As an extension to LP and MIP, BCL also provides support for formulating and solving Quadratic Programming (QP) problems, that is, problems with linear constraints with a quadratic objective function of the form
where x is the vector of decision variables, c is the cost vector, and Q is the quadratic cost coefficient matrix. The matrix Q must be symmetric. It should also be positive semi-definite if the problem is to be minimized, and negative semi-definite if it is to be maximized, because the Xpress-Optimizer solves convex QP problems. If the problem is not convex, the solution algorithms may not converge at all, or may only converge to a locally optimal solution.
- Add quadratic term
- XPRBvar x1;XPRBaddqterm(prob,x1,x1,3);
- Set quadratic term
- XPRBvar x2;XPRBsetqterm(prob,x1,x2,-7.2);
- Delete all quadratic terms
- XPRBdelqobj(prob);
Figure 3.4: Defining and accessing objective function terms in BCL.
In BCL, the quadratic part of the objective function is defined termwise, much like constraint definition. The coefficient of a quadratic term is either set to a given value (XPRBsetqterm) or its value is augmented by the given value (XPRBaddqterm). Since the quadratic part of the objective function is treated separately from the linear part, it has to be deleted with a call to XPRBdelqobj if the objective is deleted or changed to a different constraint. Note that the definition of the quadratic objective terms should always be preceded by the definition of the corresponding variables.
Unless BCL is used in Student Mode, functions XPRBprintprob, XPRBexportprob, and XPRBprintctr (if applied to the objective function) will print or output the complete problem definition to a file, including the quadratic part of the objective.
Example
The following code fragment defines the quadratic objective function
It is assumed that the variables xi have been created earlier in the program.
XPRBctr cobj;
XPRBvar x[2];
...
cobj = XPRBnewctr(prob,"OBJ",XPRB_N); /* Create the objective
constraint */
XPRBaddterm(prob,cobj, x[0],-6); /* Add the linear term */
XPRBsetobj(prob,cobj); /* Choose the obj. function */
XPRBaddqterm(prob,x[0],x[0],2); /* Add quadratic terms */
XPRBaddqterm(prob,x[0],x[1],-2);
XPRBaddqterm(prob,x[1],x[1],1);
User error handling
In this section we use a small, infeasible problem to demonstrate how the error handling and all printed messages produced by BCL can be intercepted by the user's program. This is done by defining the corresponding BCL callback functions and changing the error handling flag. If error handling by BCL is disabled, then the definition of the error callback replaces the necessity to check for the return values of the BCL functions called by a program.
User error handling may be required if a BCL program is embedded in some larger application or if the program is run under Windows from an application with windows. In all other cases it will usually be sufficient to use the error handling provided by BCL.
#include <stdio.h>
#include <setjmp.h>
#include <string.h>
#include "xprb.h"
jmp_buf model_failed; /* Marker for the longjump */
void modinf(XPRBprob prob)
{
XPRBvar x[3];
XPRBctr ctr[2], cobj;
int i;
for(i=0;i<2;i++) /* Create two integer variables */
x[i]=XPRBnewvar(prob, XPRB_UI, XPRBnewname("x_%d",i),0,100);
/* Create the constraints:
C1: 2x0 + 3x1 >= 41
C2: x0 + 2x1 = 13 */
ctr[0]=XPRBnewctr(prob,"C1",XPRB_G);
XPRBaddterm(ctr[0],x[0],2);
XPRBaddterm(ctr[0],x[1],3);
XPRBaddterm(ctr[0],NULL,41);
ctr[1]=XPRBnewctr(prob,"C2",XPRB_E);
XPRBaddterm(ctr[1],x[0],1);
XPRBaddterm(ctr[1],x[1],2);
XPRBaddterm(ctr[1],NULL,13);
/* Uncomment the following line to cause an error in the model
that triggers the user error handling: */
/* x[2]=XPRBnewvar(prob, XPRB_UI, "x_2", 10, 1); */
/* Objective: minimize x0+x1 */
cobj = XPRBnewctr(prob,"OBJ",XPRB_N);
for(i=0;i<2;i++) XPRBaddterm(cobj, x[i], 1);
XPRBsetobj(prob,cobj); /* Select objective function */
XPRBsetsense(prob,XPRB_MINIM); /* Obj. sense: minimization */
XPRBprintprob(prob); /* Print current problem */
XPRBsolve(prob,""); /* Solve the LP */
XPRBprintf(prob,
"problem status: %d LP status: %d MIP status: %d\n",
XPRBgetprobstat(prob), XPRBgetlpstat(prob),
XPRBgetmipstat(prob));
/* This problem is infeasible, that means the following command
will fail. It prints a warning if the message level is at
least 2 */
XPRBprintf(prob,"Objective: %g\n",XPRBgetobjval(prob));
for(i=0;i<2;i++) /* Print solution values */
XPRBprintf(prob,"%s:%g, ", XPRBgetvarname(x[i]),
XPRBgetsol(x[i]));
XPRBprintf(prob,"\n");
}
/**** User error handling function ****/
void XPRB_CC usererror(XPRBprob prob, void *vp, int num,
int type, const char *t)
{
printf("BCL error %d: %s\n", num, t);
if(type==XPRB_ERR) longjmp(model_failed,1);
}
/**** User printing function ****/
void XPRB_CC userprint(XPRBprob prob, void *vp, const char *msg)
{
static int rtsbefore=1;
/* Print 'BCL output' whenever a new output line starts,
otherwise continue to print the current line. */
if(rtsbefore)
printf("BCL output: %s", msg);
else
printf("%s",msg);
rtsbefore=(msg[strlen(msg)-1]=='\n');
}
int main(int argc, char **argv)
{
XPRBprob prob;
XPRBseterrctrl(0); /* Switch to error handling by the
user's program */
XPRBsetmsglevel(NULL,2); /* Set the printing flag to
printing errors and warnings */
XPRBdefcbmsg(NULL,userprint, NULL);
/* Define the printing callback func. */
if((prob=XPRBnewprob("ExplInf"))==NULL)
{ /* Initialize a new problem in BCL */
fprintf(stderr,"I cannot create the problem\n");
return 1;
}
else
if(setjmp(model_failed)) /* Set a marker at this point */
{
fprintf(stderr,"I cannot build the problem\n");
XPRBdelprob(prob); /* Delete the part of the problem
that has been created */
XPRBdefcberr(prob,NULL, NULL);
/* Reset the error callback */
return 1;
}
else
{
XPRBdefcberr(prob,usererror, NULL);
/* Define the error handling callback */
modinf(prob); /* Formulate and solve the problem */
XPRBdefcberr(prob,NULL, NULL);
/* Reset the error callback */
return 0;
}
}
Since this example defines the printing level and the printing callback function before creating the problem (that is, before BCL is initialized), we pass NULL as first argument.
If you have any comments or suggestions about these pages, please send mail to docs@dashoptimization.com.