# ----------------------------------------
# BENDERS DECOMPOSITION FOR
# THE LOCATION-TRANSPORTATION PROBLEM
# (using primal formulation of subproblem)
# ----------------------------------------
reset;
model trnloc1p.mod;
data trnloc1.dat;

# Set options
option presolve 0;

option solver gurobi;
option omit_zero_rows 1;
option display_eps .000001;


# Define suffix for storing dual rays
suffix dunbdd OUT;

problem Master: Build, Max_Ship_Cost, Total_Cost, Cut_Defn;

problem SubProblem: Ship, Ship_Cost, Supply, Demand;

option presolve_eps 1e-12;

let nCUT := 0;
let Max_Ship_Cost := 0;
let {i in ORIG} build[i] := 1;

param GAP default Infinity;

repeat { printf "\nITERATION %d\n\n", nCUT+1;

  solve SubProblem;
  printf "\n";

  if SubProblem.result = "infeasible" then { printf "Infeasible\n";
      let nCUT := nCUT + 1;
      let cut_type[nCUT] := "ray";
      let {i in ORIG} supply_price[i,nCUT] := -Supply[i].dunbdd;
      let {j in DEST} demand_price[j,nCUT] := -Demand[j].dunbdd;
      }
   else {
      if Ship_Cost - Max_Ship_Cost <= 0.000001 then break;

      let GAP := min (GAP, Ship_Cost - Max_Ship_Cost);
      option display_1col 0;
      display GAP, Ship;
      let nCUT := nCUT + 1;
      let cut_type[nCUT] := "point";
      let {i in ORIG} supply_price[i,nCUT] := Supply[i].dual;
      let {j in DEST} demand_price[j,nCUT] := Demand[j].dual;
      }
      
   printf "\nRE-SOLVING MASTER PROBLEM\n\n";
   solve Master;
   
   printf "\n";
   option display_1col 20;
   display Build;

   let {i in ORIG} build[i] := Build[i];    
};

option display_1col 0;
display Ship;
display Build;