$ondollar title Hazardous Waste (multiILP) Example 8.3 of Rardin (1998)

$offsymxref offsymlist offuelxref offuellist offupper
option limrow = 0, limcol = 0;

option mip = zoom;

sets
i  "waste sources" /1*7/,
j  "disposal sites" /1*3/,
k  "routes" /1*2/,
cs "weighting cases" /1*7/,
o  "objective functions" /1*2/;

table d(i,j,k) "travel distance (miles) from i to j on k"
   1.1    1.2    2.1    2.2    3.1    3.2
1  200    280    850   1090    900   1100
2  400    530    730    860    450    600
3  600    735    550    600    210    240
4  900   1060    450    570    180    360
5  600    640    390    440    360    510
6  900   1240    100    120    640    800
7 1230   1410    400    460   1305   1500;

table p(i,j,k) "population (000) from i to j on k"
   1.1    1.2    2.1    2.2    3.1    3.2
1   50     15    300     80    400    190
2  105     60    380    210    350    160
3  300    130    520    220    270    140
4  620    410    700    430    800    280
5  205    180    440    370    680    330
6  390    125     80     30    800    410
7  465    310    180    105   1245    790;

parameter s(i) "required K tons outflow at i"
/1 1.2, 2 0.5, 3 0.3, 4 0.7, 5 0.6, 6 0.1, 7 0.2/;

table gam(cs,o) "weight on tonmiles in case cs"
      1     2
1     1     0
2     0     1
3    10     1
4    10     5
5    10    10
6     5    10
7     1    10;
scalar gamma1, gamma2;

scalar
numsite  "number of sites to be selected" /2/;

free variables
tonmi(k)   "ton miles on routes k",
tottm      "total ton miles of shipping",
tonpopl    "total ton population of shipping"
wgtobj     "weighted objective value for case cs";

positive variables x(i,j,k) "amount shipped i to j along k";

binary variables y(j) "whether to open site j";

equations
dist           "minimize total ton mileage",
popl           "minimize total ton population",
composite      "weighted objective",
tmdef(k)       "define ton miles on routes k"
supply(i)      "ship all at i",
sites          "number of sites",
switch(i,j,k)  "site switching for i, j and k";

dist..
tottm =e= sum((i,j,k), d(i,j,k)*x(i,j,k));

popl..
tonpopl  =e= sum((i,j,k), p(i,j,k)*x(i,j,k));

composite..
wgtobj =e= gamma1 * tottm + gamma2 * tonpopl;

tmdef(k)..
tonmi(k) =e= sum((i,j), d(i,j,k)*x(i,j,k));

supply(i)..
sum((j,k), x(i,j,k)) =e= s(i);

sites..
sum(j, y(j)) =e= numsite;

switch(i,j,k)..
x(i,j,k) =l= s(i)*y(j);

model hazwaste /all/;

loop(cs,
   gamma1 = gam(cs,'1');
   gamma2 = gam(cs,'2');
   display gamma1, gamma2;
   solve hazwaste using mip minimizing wgtobj;
);