! Benders decomposition applied to            (BendersX.lng)
  multi-product\location inventory stocking.
 The problem:
    Stage 0: Decide how much to buy\make\stock, Y(i),
             at\of various origins or stock points, i.
    Stage 1: a)Random demands, DEM(s,j), occurs under scenario s, 
                 at various demand points\products,j.
             b)We then choose most profitable amount to 
               ship, X(s,i,j) under scenario s from supply point i
               to  demand point j.
               Any unsatisfied demand is lost, i.e., no sale;
! Benders decomposition is a solution method that may work
   well when the problem has the following 2 level structure:
      Master level consists of "complicating\linking" variables. 
         Some of these variables may be integer variables.
      Subproblem level has the features:
        a)if the master variables are fixed, the subproblem(s)
          are very easy to solve, e.g. separate into independent problems.
        b) subproblem is a continuous convex problem,
           e.g., an LP, so that dual prices are available for
           giving optimistic linear bounds on how the subproblem
           objective value will change as master variables are changed.
    The general iteration is:
      While( solution is not good enough:
        1) Solve master alone except for cuts from bounding subproblem.
        2) Fix master variables and solve subproblem, given the values of master variables,
        3) Based on dual prices from subproblem, add cut(s) to
            master to better describe how subproblem objective is
            affected by master level decisions.
           )

! The essential, oversimplified, nature of each cut, based on the point YF is:
    SubproblemProfit <= SubproblemProfitAt(YF) + (Y-YF)*RateOfChangeOfProfitAt(YF),
 or           BNDSUB <+ SubproblemProfitAt(YF) + (Y-YF)*DualPriceAt(YF),
    i.e., we create a linear tangent at point YF of the subproblem 
    profit function;
  
! Keywords: Benders decomposition, Cut generation, Decomposition, Inventory
            Scenario, Stochastic Optimization;
SETS:
  ORIGIN: C, UB, Y, YF; ! The stocking points;
  DEST;                 ! Destinations with demands;
  SCENE;                ! Some demand scenarios;
  CUT: RHS;             ! Set of cuts added to master;
 ! Various combinations of above simple sets;
  OXD( ORIGIN, DEST): PC;
  SXD( SCENE, DEST): DEM;
  SXOXD( SCENE, ORIGIN, DEST): X;
  CXO( CUT, ORIGIN): COF;
ENDSETS
DATA: C = 70 90 30; ! Cost/unit to stock; UB= 300 300 300; ! Upper bound on stock; DEST = 1..2; ! Two demand points; SCENE = 1..5; ! Five scenarios; DEM = 100 120 ! Demand scenarios; 130 110 160 140 180 170 190 130; PC = 140 110 ! Marginal profit contribution/unit; 100 160 ! by ORIGIN and DEST; 90 90; BigM = 999999;! A big bound; TolAbs = .5; ! Absolute tolerance for stopping when gap small; CUT = 1..20; ! Cuts allowed; ENDDATA ! Variables: Y(i) = amount stocked at\of origin i, (the linking variables) X(s,i,j) = amount shipped from i to j if scenario is s; SUBMODEL MASTER: ! Maximize bound on profit contribution - purchase cost; Max = OBJMST; ! Maximize LB from SUB + contribution from master; OBJMST = BNDSUB + ACTMST; @FREE(ACTMST); ACTMST = - @SUM( ORIGIN(i): C(i)*Y(i)); ! Stage 0 cost; @FOR( ORIGIN(i): Y(i) <= UB(i)); ! Arbitrary initial upper bound on SUBPROB contribution; BNDSUB <= BigM; ! The profit contribution bounds; @FOR( CUT( k) | k #LE# NCUTS: BNDSUB <= RHS(k) + @SUM( ORIGIN(i): COF(k,i)*Y(i)); ); ENDSUBMODEL
SUBMODEL SUBPROB: ! Maximize profit contribution averaged over all scenarios; MAX= OBJSUB; OBJSUB = @SUM( SXOXD( s,i,j): PC(i,j)*X(s,i,j))/NS; @FOR( SCENE(s): @FOR( ORIGIN(i): ! Cannot sell more than supplied; [DC] @SUM( DEST(j): X(s,i,j)) <= YF(i); ); @FOR( DEST(j): ! Cannot sell more than demanded; @SUM( ORIGIN(i): X(s,i,j)) <= DEM(s,j); ); ); ENDSUBMODEL
SUBMODEL FULL! For reference, here is the full model, not using Benders Decomposition; NSF = @SIZE( SCENE); ! Number of scenarios; Max = OBJMST; ! Maximize - purchase cost + revenue from sales; OBJMST = ACTMST + OBJSUB; @FREE(ACTMST); ACTMST = - @SUM( ORIGIN(i): C(i)*Y(i)); ! Stocking costs from stage 0/master; ! Stage 0, the stocking decisions; @FOR( ORIGIN(i): Y(i) <= UB(i)); ! Cannot stock more than capacity; ! The shipping (stage 1) of the model, averaged over all (equally likely) scenarios; OBJSUB = @SUM( SXOXD( s,i,j): PC(i,j)*X(s,i,j))/NSF; @FOR( SCENE(s): ! For each scenario; @FOR( ORIGIN(i): ! Cannot sell more from source i than supplied/stocked; [DCF] @SUM( DEST(j): X(s,i,j)) <= Y(i); ); @FOR( DEST(j): ! Cannot sell more at destination j than demanded under scenario s; @SUM( ORIGIN(i): X(s,i,j)) <= DEM(s,j); ); ); ENDSUBMODEL
CALC: ! If alternatively, we want to solve the full model w/o Benders...; ! @SOLVE( FULLMODEL); NCUTS = 0; ! No cuts so far; @SET('DUALCO', 1); ! Make sure dual computations are on; @SET('TERSEO', 2); ! Turn off default output; NS = @SIZE( SCENE);! Number scenarios; MORE = 1; ! Keep looping if MORE > 0; @WHILE ( MORE: @WRITE(@NEWLINE(1), 'Solve MASTER at iteration ', NCUTS, @NEWLINE(1)); @SOLVE( MASTER); @WRITE(' Obj/Optimistic bound= ', OBJMST, @NEWLINE(1)); OBJBND = OBJMST; ! Remember optimistic bound; FEASOBJ = ACTMST; ! Remember obj of current feasible soln; ! Variable value rule: Any variable appearing on LHS in a CALC section has its value fixed permanently (until next appearance on a LHS). Any other variable is a decision variable and is liable to change in the next optimization. ! Prepare to send master solution to SUB; @FOR( ORIGIN(i): YF(i) = Y(i)); ! YF(i) is the fixed value of current Y(i); @WRITE( 'Solve SUBPROB ', @NEWLINE(1)); @SOLVE( SUBPROB); @WRITE(' OBJSUB= ', OBJSUB, @NEWLINE(1)); FEASOBJ = FEASOBJ + OBJSUB; @WRITE(' Obj of current feasible solution= ',FEASOBJ,@NEWLINE(1)); ! We are done if the bound and the feasible solution are close enough or we have run out of cut space; @IFC( @ABS(FEASOBJ - OBJBND) #GT# TolAbs #AND# NCUTS #LT# @SIZE(CUT): ! Keep going. Add linear bounds from SUB to MASTER; ! Compute constant term of linear subprob bound; RHST = OBJSUB; NCUTS = NCUTS + 1; ! Build up coefficients in new cut; @FOR( ORIGIN(i): COF(NCUTS,i) = 0); @FOR(SCENE(s): @FOR( ORIGIN(i): TEMP = @DUAL( DC(s,i)); ! Get dual prices; RHST = RHST - YF(i)*TEMP; COF(NCUTS,i) = COF(NCUTS,i) + TEMP; ); ); RHS(NCUTS) = RHST; @WRITE(' Add cut to MASTER:',@NEWLINE(1),' BNDSUB <= ',RHS(NCUTS)); @FOR( ORIGIN(i): @WRITE( ' +',COF(NCUTS,i),'*Y(',i,')'); ); @WRITE(';',@NEWLINE(1)); @ELSE MORE = 0; ! Signal quit; ); ! @IFC ; ); ! @WHILE; @WRITE(@NEWLINE(1),' Bound on optimal obj = ', @FORMAT( OBJBND,'12.2f'), @NEWLINE(1)); @WRITE( ' Best feasible solution obj = ', @FORMAT(FEASOBJ,'12.2f'), @NEWLINE(1)); @WRITE(' Best decisions in master( Amount to stock at each origin:',@NEWLINE(1)); @FOR( ORIGIN(i): @WRITE(' ', ORIGIN(i), ' ')); @WRITE(@NEWLINE(1)); @FOR( ORIGIN(i): @WRITE( ' ', @FORMAT(YF(i),'11.2f')); ); @WRITE( @NEWLINE(2)); @WRITE(' Optimal decision under each scenario:', @NEWLINE(1)); @WRITE(' Scenario Origin Destination Ship', @NEWLINE(1)); @FOR( SXOXD( s,i,j) | X(s,i,j) #GT# 0: @WRITE( ' ',s,' ', ORIGIN(i),' ', DEST(j),' ', @FORMAT(X(s,i,j),'10.2f'), @NEWLINE(1)); ); @IFC( NCUTS #EQ# @SIZE(CUT): @WRITE(' Ran out of cut space.', @NEWLINE(1))); ENDCALC