solnpool.gms : Cplex Solution Pool for a Simple Facility Location Problem

Description

A simple version of a facility location problem is used to show how the
solution pool and the tools associated with it work. This example is taken
from the CPLEX 11 User's Manual (ILOG, CPLEX 11 User's Manual, 2007)

A company is considering opening as many as four warehouses in order to serve
nine different regions. The goal is to minimize the sum of fixed costs
associated with opening warehouses as well as the various transportation
costs incurred to ship goods from the warehouses to the regions.

Whether or not to open a warehouse is represented by binary variable ow.
Whether or not to ship goods from warehouse i to region j is represented
by binary variable oa.

Each region needs a specified amount of goods, and each warehouse can store
only a limited quantity of goods. In addition, each region must be served
by exactly one warehouse.


Small Model of Type : MIP


Category : GAMS Model library


Main file : solnpool.gms

$title CPLEX Solution Pool for a Simple Facility Location Problem (SOLNPOOL,SEQ=326)

$onText
A simple version of a facility location problem is used to show how the
solution pool and the tools associated with it work. This example is taken
from the CPLEX 11 User's Manual (ILOG, CPLEX 11 User's Manual, 2007)

A company is considering opening as many as four warehouses in order to serve
nine different regions. The goal is to minimize the sum of fixed costs
associated with opening warehouses as well as the various transportation
costs incurred to ship goods from the warehouses to the regions.

Whether or not to open a warehouse is represented by binary variable ow.
Whether or not to ship goods from warehouse i to region j is represented
by binary variable oa.

Each region needs a specified amount of goods, and each warehouse can store
only a limited quantity of goods. In addition, each region must be served
by exactly one warehouse.


The following GAMS program demonstrates a number of different approaches to
collecting solution pools. GAMS will store the individual solutions in GDX
containers/files which can then be further used by other programs or the same
GAMS run. CPLEX will name these GDX containers 'soln_loc_pNN.gdx', where NN
will be the serial number of the solution found. To manage the different
solutions, the names of the GDX containers created by CPLEX will be stored
in solnpool.gdx in the set 'index' using the set elements file*.

Eight examples are solved in this gams run.

Keywords: mixed integer linear programming, facility location problem, CPLEX
          solution pool
$offText

Set
   i 'warehouses' / w1*w4 /
   j 'regions'    / r1*r9 /;

Parameter
   f(i) 'fixed costs' / w1 130, w2 150, w3 170, w4 180 /
   c(i) 'capacity'    / w1  90, w2 110, w3 130, w4 150 /
   d(j) 'demand'      / r1 10, r2 10, r3 12, r4 15, r5 15
                        r6 15, r7 20, r8 20, r9 30        /;

Table t(j,i) 'transport costs'
        w1  w2  w3  w4
   r1   10  30  25  55
   r2   10  25  25  45
   r3   20  23  30  40
   r4   25  10  26  40
   r5   28  12  20  29
   r6   36  19  16  22
   r7   40  39  22  27
   r8   75  65  55  35
   r9   34  43  41  62;

Variable
   totcost 'total cost'
   fcost   'fixed cost'
   tcost   'transportation cost'
   ow(i)   'indicator for open warehouse'
   oa(i,j) 'indicator for open shipment arc warehouse to region';

Binary Variable ow, oa;

Equation
   deftotcost 'definition total cost'
   deffcost   'definition fixed cost'
   deftcost   'definition transportation cost'
   defwcap(i) 'limit utilization of warehouse by its capacity'
   onew(j)    'only one warehouse per region'
   defow(i,j) 'warehouse open if shipment from i to j';

deftotcost.. totcost =e= fcost + tcost;

deffcost..   fcost   =e= sum(i, f(i)*ow(i));

deftcost..   tcost   =e= sum((i,j), t(j,i)*oa(i,j));

defwcap(i).. sum(j, d(j)*oa(i,j)) =l= c(i);

onew(j)..    sum(i, oa(i,j)) =e= 1;

defow(i,j).. ow(i) =g= oa(i,j);

Model loc / all /;

* Define sets, parameters and files to hold solutions
Set
   soln           'possible solutions in the solution pool' / file1*file1000 /
   solnpool(soln) 'actual solutions';

Scalar cardsoln 'number of solutions in the pool';

Alias (soln,s1,s2), (*,u);

Parameter
   owX(soln,i)    'warehouse indicator by solution'
   oaX(soln,i,j)  'arc indicator by solution'
   xcostX(soln,*) 'cost structure by solution';

File fcpx / cplex.opt /;

option limRow = 0, limCol = 0, optCr = 0, mip = cplex;
loc.optFile   = 1;
loc.solPrint  = %solPrint.quiet%;
loc.savePoint = 1;

* The code to load different solution from gdx files will be used
* several times in this program and we therefore copy it into an include file.
$onEcho > readsoln.gms
execute_load 'solnpool.gdx', solnpool = index;
cardsoln = card(solnpool);

display cardsoln;
oaX(soln,i,j)  = 0;
owX(soln,i)    = 0;
xcostX(soln,u) = 0;

loop(solnpool(soln),
   put_utility 'gdxin' / solnpool.te(soln);
   execute_loadpoint;
   oaX(soln,i,j)             = round(oa.l(i,j));
   owX(soln,i)               = round(ow.l(i));
   xcostX(soln,'totcost')    = totcost.l;
   xcostX(soln,'tcost')      = tcost.l;
   xcostX(soln,'fcost')      = fcost.l;
   xcostX(soln,'fcost^0.96') = fcost.l**0.96;
);
* Restore the solution reported to GAMS
execute_loadpoint 'loc_p.gdx';
$offEcho

* 1. Collect the incumbents found during the regular optimize procedure
*    The Cplex option 'solnpool' triggers the collection of solutions in
*    the GDX container solnpool.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1';
solve loc min totcost using mip;
$include readsoln
display xcostX;

* 2. Use the populate procedure instead of regular optimize procedure (option
*    'solnpoolpop 2'). By default we will generate 20 solutions determined by
*    the default of option populatelim (only valid for one thread). This is a 
*    simple model which is quickly solved with heuristics and cuts, so we need 
*    to let Cplex retain sufficient exploration space to find alternative solutions. 
*    This is done with option 'solnpoolintensity 4'. With solutions where the optimal 
*    solution can not so quickly  be found, the default for this option should be suffict.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2';
solve loc min totcost using mip;
$include readsoln
display xcostX;

* 3. Now it gets more complicated. After the first call to populate we instruct
*    the GAMS/Cplex link to execute a GAMS program 'simple.gms' that decides if
*    the populate procedure should be called again. Moreover, we will instruct
*    Cplex to delete the solution 11 to 15 from the pool. After this two rounds
*    we will have 35 solutions in the pool. Note the use of %ncall% which expands
*    to an integer counting the number of previous execution of 'simple.gms'.
*    A nonzero return code of the gams run called by Cplex will signal to
*    Cplex that this is was the final call.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
            / 'solnpoolpoprepeat simple.gms'
            / 'solnpoolpopdel delsol.txt';
$onEchoV > simple.gms
$ifE %ncalls%>0 abort 'Terminate search'
File f / delsol.txt /;
put  f '11 12 13 14 15';
$offEcho

solve loc min totcost using mip;
$include readsoln
abort$(cardsoln <> 35) 'Expected to get 35 solutions';
display xcostX;

* 4. Next we call the populate procedure, but we want solutions that are
*    within 3% of the optimum. If we find more than 20, we are willing to
*    call the populate procedure once more. It turns out that there are 11
*    solutions of this quality, so no need to call populate again.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
            / 'solnpoolpoprepeat simple.gms'
            / 'solnpoolgap 0.03';
solve loc min totcost using mip;
$include readsoln
display xcostX;

* 5. Lets look at the diversity of the solution by counting the differences
*    between the shipment indicator variables. Lets limit the number of
*    solutions in the pool by 10 and require solution within 10% of the
*    optimum.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
            / 'solnpoolcapacity 10'   / 'solnpoolgap 0.1';
solve loc min totcost using mip;
$include readsoln

Scalar aggdiff 'aggregated differences' / 0 /;

loop((s1,s2)$(not sameas(s1,s2) and solnpool(s1) and solnpool(s2)),
   aggdiff = aggdiff + sum((i,j), oaX(s1,i,j) xor oaX(s2,i,j));
);
display aggdiff;

* 6. We repeat the experiment by now setting the solution pool replacement
*    strategy to 'diversity' and let the populate procedure find many more
*    solutions, we should see an increase in the aggregated difference
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
            / 'solnpoolcapacity 10'   / 'solnpoolgap 0.1'
            / 'populatelim 10000'     / 'solnpoolreplace 2';
solve loc min totcost using mip;
$include readsoln

Scalar aggdiffX 'aggregated differences' / 0 /;

loop((s1,s2)$(not sameas(s1,s2) and solnpool(s1) and solnpool(s2)),
   aggdiffX = aggdiffX + sum((i,j), oaX(s1,i,j) xor oaX(s2,i,j));
);
display aggdiffX;
abort$(aggdiffX < aggdiff) 'We expected *more* diversity';

* 7. We can fine tune diversity by using a diversity filter. Suppose that
*    facilities w1 and w2 are open. Let a solution keeping those two facilities
*    open be the reference. We use a diversity filter to stipulate that any
*    solution added to the solution pool must differ from the reference by
*    decisions to open or close at least two other facilities. The following
*    filter enforces this diversity by specifying a minimum diversity of 2.
*    Note that the reference solution becomes the incumbent and is reported as
*    the first solution in the pool.
put fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
       / 'divfltlo 2'            / 'writeflt ow2.flt';

Parameter owrefsol(i) / w1 1, w2 1, w3 0, w4 0 /;

loop(i, put / 'ow.divflt("' i.tl:0 '") ' owrefsol(i):0:0;);
putClose fcpx;

solve loc min totcost using mip;
$include readsoln
display owX;
loop(solnpool(soln)$(ord(soln) > 1),
   abort$(sum(i,abs(owX(soln,i) - owrefsol(i))) < 2) 'solution differs to little';);

* 8. We can also implement more complicated filters using the incumbent filter.
*    For example, we want to enforce that transportation cost is less than
*    fixedcost**0.96. The incumbent filter is implemented as part of the BCH
*    facility (http://www.gams.com/docs/bch.htm). Whenever a new incumbent is
*    found by Cplex, Cplex will run the GAMS program 'incbflt'. The incumbent
*    is provided in a GDX container called 'bchout_i'. The GAMS program can
*    inspect and decide if the incumbent should be accepted or rejected by
*    Cplex. A nonzero return code of the gams run called by Cplex will signal
*    to Cplex that the incumbent is accecpted.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2'
            / 'userincbcall incbflt.gms';
$onEcho > incbflt.gms
Variable tcost, fcost;
$gdxIn bchout_i
$load tcost fcost
abort$(tcost.l < fcost.l**0.96) 'Accept';
$offEcho

solve loc min totcost using mip;
$include readsoln
display xcostX;
loop(solnpool(soln),
   abort$(xcostX(soln,'tcost') >= xcostX(soln,'fcost')) 'tcost too big';);

* 9. We might sometimes want to order the solutions, for example, best
*    to worst. For this we can use the $libInclude rank utility. Note that rank
*    returns the sorted indices.
putClose fcpx 'solnpool solnpool.gdx' / 'threads 1' / 'solnpoolintensity 4' / 'solnpoolpop 2';
solve loc min totcost using mip;
$include readsoln

Parameter
   uval(soln) 'unsorted objective values'
   idx(soln)  'sorted index position';

uval(soln) = xcostX(soln,'totcost');

* sort solutions via rank
$libInclude rank uval solnpool idx
File fsorted / 'sorted.txt' /;
put  fsorted 'Sorted solutions:';

* Get the best five in a file
loop(solnpool(soln),
   if(idx(soln) <= 5,
      put #(idx(soln)+1) @1 soln.tl:0 ':' @10 uval(soln);
   );
);
putClose;

* Get all solutions sorted in a new parameter
Parameter sval(soln) 'sorted objective values';
sval(soln + (idx(soln) - ord(soln)))$solnpool(soln) = uval(soln);
display sval;