$ontext

   Transportation problem with stochastic demands.
   OSL-SE formulation

   Erwin Kalvelagen, January 2003


$offtext

sets
  i 'factories' /f1*f3/
  j 'distribution centers' /d1*d5/
  s 'individual scenarios' /lo,mid,hi/
;

parameter capacity(i) /f1 500, f2 450, f3 650/;

table demand(j,s) 'possible outcomes for demand'
      lo  mid   hi
  d1 150  160  170
  d2 100  120  135
  d3 250  270  300
  d4 300  325  350
  d5 600  700  800
;

table prob(j,s) 'probabilities for table demand'
      lo  mid   hi
  d1 .25  .50  .25
  d2 .25  .50  .25
  d3 .25  .50  .25
  d4 .30  .40  .30
  d5 .30  .40  .30
;

loop(j, abort$(abs(sum(s,prob(j,s))-1)>0.001) "probabilities don't add up");


*
* set up joint probabilities
* total scenarios =  3**5 = 243
*
set n 'nodes' /root,js1*js243/;
set js(n) /js1*js243/;
parameter jdemand(j,js) 'joint distribution: outcomes';
parameter jprob(js)   'joint distribution: probabilities';

alias (s,s1,s2,s3,s4,s5);
set current(js);
current('js1') = yes;
loop((s1,s2,s3,s4,s5),
   jdemand('d1',current) = demand('d1',s1);
   jdemand('d2',current) = demand('d2',s2);
   jdemand('d3',current) = demand('d3',s3);
   jdemand('d4',current) = demand('d4',s4);
   jdemand('d5',current) = demand('d5',s5);
   jprob(current) = prob('d1',s1)*prob('d2',s2)*prob('d3',s3)*prob('d4',s4)*prob('d5',s5);
   current(js) = current(js-1);
);
display jdemand,jprob;

scalar sumprob;
sumprob = sum(js, jprob(js));
display sumprob;
abort$(abs(sumprob-1)>0.00001) "joint probabilities don't add up";

table transcost(i,j) 'unit transportation cost'
       d1    d2    d3    d4    d5
  f1   2.49  5.21  3.76  4.85  2.07
  f2   1.46  2.54  1.83  1.86  4.76
  f3   3.26  3.08  2.60  3.76  4.45
;

scalar prodcost  'unit production cost' /14/;
scalar price     'sales price' /24/;
scalar wastecost 'cost of removal of overstocked products' /4/;



*-----------------------------------------------------------------------
* We form here the deterministic equivalent
*-----------------------------------------------------------------------

variables
   salesw(j,n)    'stochastic version of sales'
   wastew(j,n)    'stochastic version of waste'
   received(j,n)  'amount of product received in distribution center'
   ship(i,j,n)    'shipments'
   product(i,n)   'production'
   profit
;
positive variable ship,product,salesw,wastew;

equations
   objw
   sellingw(j,n)
   receive(j,n)
   production(i,n)
;

objw.. profit =e= sum((j,js),price*jprob(js)*salesw(j,js))
                 - sum((i,j), transcost(i,j)*ship(i,j,'root'))
                 - sum((j,js), wastecost*jprob(js)*wastew(j,js))
                 - sum(i,prodcost*product(i,'root'));

receive(j,'root')..     received(j,'root') =e= sum(i, ship(i,j,'root'));
sellingw(j,js)..        received(j,'root') =e= salesw(j,js)+wastew(j,js);
production(i,'root')..  product(i,'root') =e= sum(j, ship(i,j,'root'));

product.up(i,'root') = capacity(i);
salesw.up(j,js) = jdemand(j,js);


option lp=oslse;
model stoch /objw,production,receive,sellingw/;
solve stoch maximizing profit using lp;

display ship.l,product.l,salesw.l,wastew.l;