/********************************************************
Xpress-BCL C++ Example Problems
===============================
file xbpurch.cxx
````````````````
Purchasing problem using SOS-2.
(c) 2008 Fair Isaac Corporation
author: S.Heipcke, Jan. 2000, rev. Mar. 2011
********************************************************/
/* WARNING. This model is not for the novice, but it does contain many *
* useful ideas. *
* The formulation is tricky because the discounts are all-quantity, *
* so the graph of cost against quantity purchased is discontinuous. *
* To maintain sanity in the special ordered set formulation, we must *
* coarsen the discontinuity by stopping purchases just at the break *
* point. */
#include <iostream>
#include <cstdio>
#include "xprb_cpp.h"
using namespace std;
using namespace ::dashoptimization;
#define PARAMSFILE XPRBDATAPATH "/purchase/params.dat"
#define MAXPERCFILE XPRBDATAPATH "/purchase/maxperc.dat"
#define REQFILE XPRBDATAPATH "/purchase/required.dat"
#define PRICEFILE XPRBDATAPATH "/purchase/pricebk.dat"
/****TABLES****/
int NS; /* Number of suppliers */
int NB; /* Number of price breaks */
int NB2; /* Useful parameter */
double **UC; /* Unit cost */
double **BR; /* Breakpoints (quantities at which unit cost changes) */
double **X; /* Coarsened break points */
double **C; /* Total cost at break points */
double *DELTA; /* Coarsening factors */
double *MAXPERC; /* Maximum percentage from each supplier */
double REQ; /* Total quantity required */
double delta = 0.10; /* Base coarsening factor */
XPRBprob p("Purchase"); /* Initialize a new problem in BCL */
/***********************************************************************/
void modPurchase()
{
XPRBexpr lobj, *lr, lc;
int s,b;
XPRBvar *x,**lam;
/****VARIABLES****/
x = new XPRBvar[NS]; /* Quantity to purchase from supplier s */
for(s=0;s<NS;s++) x[s] = p.newVar("x");
lam = new XPRBvar*[NS];
for(s=0;s<NS;s++)
{
lam[s] = new XPRBvar[NB2];
for(b=0;b<NB2;b++)
lam[s][b] = p.newVar(XPRBnewname("lam_%d",s+1));
} /* Weights at breakpoint b for supplier s */
/****OBJECTIVE****/
for(s=0;s<NS;s++) /* Minimize the sum of costs*weights */
for(b=0;b<NB2;b++) lobj += C[s][b]*lam[s][b];
p.setObj(p.newCtr("OBJ",lobj)); /* Set objective function */
/****CONSTRAINTS****/
/* Define x and also order the lam variables by breakpoint quantities */
lr = new XPRBexpr[NS];
for(s=0;s<NS;s++)
{
for(b=0;b<NB2;b++) lr[s] += X[s][b]*lam[s][b];
p.newCtr("DefX", lr[s] == x[s]);
}
/* The convexity constraint (lam sum to 1) */
for(s=0;s<NS;s++)
{
lc=0;
for(b=0;b<NB2;b++) lc += lam[s][b];
p.newCtr("Conv", lc == 1);
}
/* The minimum quantity that must be bought */
lc=0;
for(s=0;s<NS;s++) lc += x[s];
p.newCtr("MustBuy", lc >= REQ);
/****BOUNDS****/
/* No more than the maximum percentage from each supplier */
for(s=0;s<NS;s++) x[s].setUB(MAXPERC[s]*REQ/100.0);
/****SETS****/
/* Define the lam as SOS2 as we can linearly interpolate between the
* breakpoints. The weights are the DefX coefficients augmented by 1
because BCL does not accept the weight 0. */
for(s=0;s<NS;s++)
{
for(b=0;b<NB2;b++) lr[s] += lam[s][b];
p.newSos("SetLam", XPRB_S2, lr[s]);
}
/****SOLVING + OUTPUT****/
p.exportProb(XPRB_MPS,"purc"); /* Output to MPS file */
p.setSense(XPRB_MINIM); /* Choose the sense of the optimization */
p.mipOptimize(""); /* Solve the MIP-problem */
cout << "Objective: " << p.getObjVal() << endl; /* Get objective value */
for(s=0;s<NS;s++) /* Print out the solution values */
cout << x[s].getName() << ":" << x[s].getSol() << " ";
cout << endl;
for(s=0;s<NS;s++)
for(b=0;b<NB2;b++)
cout << lam[s][b].getName() << ":" << lam[s][b].getSol() << " ";
cout << endl;
delete [] BR;
delete [] UC;
delete [] lr;
for(s=0;s<NS;s++) delete [] lam[s];
delete [] lam;
delete [] x;
}
/***********************************************************************/
/**** Read data from files ****/
void readData()
{
FILE *datafile;
int s,b;
double val1,val2;
/* Read the parameter file */
datafile=fopen(PARAMSFILE,"r");
XPRBreadlinecb(XPRB_FGETS, datafile, 200, "i", &NS);
XPRBreadlinecb(XPRB_FGETS, datafile, 200, "i", &NB);
fclose(datafile);
NB2 = 2*NB;
/* Now define the tables that we will use, using the sizes we
have just read. */
UC= new double*[NS];
for(s=0;s<NS;s++) UC[s]=new double[NB];
BR = new double*[NS];
for(s=0;s<NS;s++) BR[s]=new double[NB];
X=new double*[NS];
for(s=0;s<NS;s++) X[s]=new double[NB2];
C=new double*[NS];
for(s=0;s<NS;s++) C[s]=new double[NB2];
DELTA=new double[NB2];
MAXPERC=new double[NS];
/* Define coarsening factors. */
DELTA[0]=0;
DELTA[1] = -1*delta;
for(b=2;b<NB2;b++) DELTA[b] = -1*DELTA[b-1];
/* Read the price and breakpoint data file */
for(s=0;s<NS;s++)
for(b=0;b<NB;b++)
{
UC[s][b]=0;
BR[s][b]=0;
}
datafile=fopen(PRICEFILE,"r");
while(XPRBreadlinecb(XPRB_FGETS, datafile, 200, "i,i,g,g", &s, &b,
&val1, &val2)==4)
{
UC[s-1][b-1]=val1;
BR[s-1][b-1]=val2;
}
fclose(datafile);
/* Read the percentages data file */
datafile=fopen(MAXPERCFILE,"r");
XPRBreadarrlinecb(XPRB_FGETS, datafile, 200, "g,", MAXPERC, NS);
fclose(datafile);
/* Read the required amount data file */
datafile=fopen(REQFILE,"r");
XPRBreadlinecb(XPRB_FGETS, datafile, 200, "g", &REQ);
fclose(datafile);
/* We now pick out the data from the tables into which they have
been read. */
for(s=0;s<NS;s++)
{
C[s][0] = 0; /* First total cost and (new) breakpoint are (0,0) */
X[s][0] = 0;
for(b=1;b<NB2;b++)
{ /* Rest calculated... */
C[s][b] = UC[s][b/2] * BR[s][(b-1)/2]; /* ...unit price*quantity */
X[s][b] = BR[s][(b-1)/2] + DELTA[b]; /* ...coarsened grids */
}
}
}
/***********************************************************************/
int main(int argc, char **argv)
{
readData(); /* Data input from file */
modPurchase(); /* Formulate and solve the problem */
return 0;
}
|