multi3.run

Warning

The original AMPL book does not reflect many of the latest features available in AMPL. To programmatically interact with your models you should use APIs such as our popular Python API.

To explore examples showcasing these recent features, please visit:

• AMPL Streamlit Apps

• AMPL+Python Book

• AMPL Colab Notebooks

multi3.run

# ----------------------------------------
# DANTZIG-WOLFE DECOMPOSITION FOR
# MULTI-COMMODITY TRANSPORTATION
# ----------------------------------------

model multi3.mod;
data multi3.dat;

param nIter default 0;
param nNegRC;

let {p in PROD} nPROP[p] := 0;
let {p in PROD} price_convex[p] := 1;
let {i in ORIG, j in DEST} price[i,j] := 0;

option omit_zero_rows 1;
option display_1col 0;
option display_eps .000001;

# ----------------------------------------------------------

problem MasterI: Artificial, Weight, Excess, Multi, Convex;
problem SubI: Artif_Reduced_Cost, Trans, Supply, Demand;

repeat {

   let nIter := nIter + 1;
   let nNegRC := 0;
   printf "\nPHASE I -- ITERATION %d\n", nIter;
   problem SubI;

   for {p in PROD} { printf "\nPRODUCT %s\n\n", p;

      let P := p;
      fix Trans; unfix {i in ORIG, j in DEST} Trans[i,j,p];

      solve SubI;
      printf "\n";
      display {i in ORIG, j in DEST} Trans[i,j,p];

      if Artif_Reduced_Cost < - 0.00001 then {
         let nNegRC := nNegRC + 1;
         let nPROP[p] := nPROP[p] + 1;
         let {i in ORIG, j in DEST}
            prop_ship[i,j,p,nPROP[p]] := Trans[i,j,p];
         let prop_cost[p,nPROP[p]] :=
            sum {i in ORIG, j in DEST} cost[i,j,p] * Trans[i,j,p];
      };
   };

   if nNegRC = 0 then {
      printf "\n*** NO FEASIBLE SOLUTION ***\n";
      break;
   };

   solve MasterI;
   printf "\n";
   display Weight; display Multi.dual;
   display {i in ORIG, j in DEST}
      limit[i,j] - sum {p in PROD, k in 1..nPROP[p]}
         prop_ship[i,j,p,k] * Weight[p,k];

   if Excess <= 0.00001 then break;
   else {
      let {i in ORIG, j in DEST} price[i,j] := Multi[i,j].dual;
      let {p in PROD} price_convex[p] := Convex[p].dual;
   };
};

# ----------------------------------------------------------

printf "\nSETTING UP FOR PHASE II\n";

problem MasterII: Total_Cost, Weight, Multi, Convex;
problem SubII: Reduced_Cost, Trans, Supply, Demand;

solve MasterII;
printf "\n";
display Weight; display Multi.dual; display Multi.slack;

let {i in ORIG, j in DEST} price[i,j] := Multi[i,j].dual;
let {p in PROD} price_convex[p] := Convex[p].dual;

repeat {

   let nIter := nIter + 1;
   let nNegRC := 0;
   printf "\nPHASE II -- ITERATION %d\n\n", nIter;
   problem SubII;

   for {p in PROD} { printf "\nPRODUCT %s\n\n", p;

      let P := p;
      fix Trans; unfix {i in ORIG, j in DEST} Trans[i,j,p];

      solve SubII;
      printf "\n";
      display {i in ORIG, j in DEST} Trans[i,j,p];

      if Reduced_Cost < - 0.00001 then  {
         let nNegRC := nNegRC + 1;
         let nPROP[p] := nPROP[p] + 1;
         let {i in ORIG, j in DEST}
            prop_ship[i,j,p,nPROP[p]] := Trans[i,j,p];
         let prop_cost[p,nPROP[p]] :=
            sum {i in ORIG, j in DEST} cost[i,j,p] * Trans[i,j,p];
      };
   };

   if nNegRC = 0 then break;

   solve MasterII;

   printf "\n";
   display Weight;

   let {i in ORIG, j in DEST} price[i,j] := Multi[i,j].dual;
   let {p in PROD} price_convex[p] := Convex[p].dual;
};

# ----------------------------------------------------------

printf "\nPHASE III\n";

let {i in ORIG, j in DEST, p in PROD}
   Trans[i,j,p] := sum {k in 1..nPROP[p]} prop_ship[i,j,p,k] * Weight[p,k];

param true_Total_Cost
   := sum {i in ORIG, j in DEST, p in PROD} cost[i,j,p] * Trans[i,j,p].val;

printf "\n";
display true_Total_Cost;
display Trans;