Benders2Stage.java
1package com.gams.examples.benders;
2 
3 import java.io.File;
4 
7 import com.gams.api.GAMSGlobals;
8 import com.gams.api.GAMSJob;
10 import com.gams.api.GAMSModifier;
11 import com.gams.api.GAMSOptions;
12 import com.gams.api.GAMSParameter;
14 import com.gams.api.GAMSSetRecord;
15 import com.gams.api.GAMSVariable;
17 import com.gams.api.GAMSWorkspace;
19 
35 public class Benders2Stage {
36 
37  public static void main(String[] args) {
38  GAMSWorkspaceInfo wsInfo = new GAMSWorkspaceInfo();
39  if (args.length > 0)
40  wsInfo.setSystemDirectory( args[0] );
41 
42  File workingDirectory = new File(System.getProperty("user.dir"), "Benders2Stage");
43  workingDirectory.mkdir();
44  wsInfo.setWorkingDirectory(workingDirectory.getAbsolutePath());
45 
46  GAMSWorkspace ws = new GAMSWorkspace( wsInfo );
47 
48  GAMSJob dataJob = ws.addJobFromString(data);
49 
50  GAMSOptions optData = ws.addOptions();
51  optData.defines("useBig", "1");
52  optData.defines("nrScen", "100");
53 
54  dataJob.run(optData);
55 
56  optData.dispose();
57 
58  GAMSParameter scenarioData = dataJob.OutDB().getParameter("ScenarioData");
59 
60  GAMSOptions opt = ws.addOptions();
61  opt.defines("datain", dataJob.OutDB().getName());
62  int maxiter = 40;
63  opt.defines("maxiter", Integer.toString(maxiter));
64  opt.setAllModelTypes( "cplex" );
65 
66  GAMSCheckpoint cpMaster = ws.addCheckpoint();
67  GAMSCheckpoint cpSub = ws.addCheckpoint();
68 
69  ws.addJobFromString(masterModel).run(opt, cpMaster, dataJob.OutDB());
70 
71  GAMSModelInstance masteri = cpMaster.addModelInstance();
72  GAMSParameter cutconst = masteri.SyncDB().addParameter("cutconst", 1, "Benders optimality cut constant");
73  GAMSParameter cutcoeff = masteri.SyncDB().addParameter("cutcoeff", 2, "Benders optimality coefficients");
74  GAMSVariable theta = masteri.SyncDB().addVariable("theta", 0, GAMSGlobals.VarType.FREE, "Future profit function variable");
75  GAMSParameter thetaFix = masteri.SyncDB().addParameter("thetaFix", 0, "");
76  masteri.instantiate("masterproblem max zmaster using lp", opt,
77  new GAMSModifier[] { new GAMSModifier(cutconst), new GAMSModifier(cutcoeff), new GAMSModifier(theta,GAMSGlobals.UpdateAction.FIXED,thetaFix) }
78  );
79 
80  ws.addJobFromString(subModel).run(opt, cpSub, dataJob.OutDB());
81 
82  GAMSModelInstance subi = cpSub.addModelInstance();
83  GAMSParameter received = subi.SyncDB().addParameter("received", 1, "units received from master");
84  GAMSParameter demand = subi.SyncDB().addParameter("demand", 1, "stochastic demand");
85  subi.instantiate("subproblem max zsub using lp", opt,
86  new GAMSModifier[] { new GAMSModifier(received), new GAMSModifier(demand) }
87  );
88 
89  opt.dispose();
90 
91  double lowerbound = Double.NEGATIVE_INFINITY, upperbound = Double.POSITIVE_INFINITY, objmaster = Double.POSITIVE_INFINITY;
92  int iter = 1;
93  do {
94  System.out.println("Iteration: " + iter);
95  // Solve master
96  if (1 == iter) // fix theta for first iteration
97  thetaFix.addRecord().setValue( 0 );
98  else
99  thetaFix.clear();
100 
102 
103  System.out.println(" Master " + masteri.getModelStatus() + " : obj=" + masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel());
104  if (1 < iter)
105  upperbound = masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel();
106  objmaster = masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel() - theta.getFirstRecord().getLevel();
107 
108  // Set received from master
109  received.clear();
110  for (GAMSVariableRecord r : masteri.SyncDB().getVariable("received")) {
111  received.addRecord( r.getKeys()).setValue( r.getLevel() );
112  cutcoeff.addRecord(new String[] { Integer.toString(iter), r.getKey(0) });
113  }
114 
115  cutconst.addRecord(Integer.toString(iter));
116  double objsub = 0.0;
117  for(GAMSSetRecord s : dataJob.OutDB().getSet("s")) {
118  demand.clear();
119  for(GAMSSetRecord j : dataJob.OutDB().getSet("j")) {
120  demand.addRecord(j.getKeys()).setValue( scenarioData.findRecord(new String[] {s.getKey(0), j.getKey(0)}).getValue() );
121  }
122 
124  System.out.println(" Sub " + subi.getModelStatus() + " : obj=" + subi.SyncDB().getVariable("zsub").getFirstRecord().getLevel());
125 
126  double probability = scenarioData.findRecord(new String[] {s.getKey(0), "prob"}).getValue();
127  objsub += probability * subi.SyncDB().getVariable("zsub").getFirstRecord().getLevel();
128  for (GAMSSetRecord j : dataJob.OutDB().getSet("j")) {
129  GAMSParameterRecord record = cutconst.findRecord(Integer.toString(iter));
130  double newValue = record.getValue() + probability * subi.SyncDB().getEquation("market").findRecord(j.getKeys()).getMarginal() * demand.findRecord(j.getKeys()).getValue();
131  record.setValue( newValue );
132  record = cutcoeff.findRecord(new String[] {Integer.toString(iter), j.getKey(0)});
133  newValue = record.getValue() + probability * subi.SyncDB().getEquation("selling").findRecord(j.getKeys()).getMarginal();
134  record.setValue( newValue );
135  }
136  }
137 
138  lowerbound = Math.max(lowerbound, objmaster + objsub);
139 
140  iter++;
141  if (iter == maxiter + 1)
142  throw new GAMSException("Benders out of iterations");
143 
144  System.out.println(" lowerbound: " + lowerbound + " upperbound: " + upperbound + " objmaster: " + objmaster);
145  } while ((upperbound - lowerbound) >= 0.001 * (1 + Math.abs(upperbound)));
146 
147  masteri.dispose();
148  subi.dispose();
149  }
150 
151  static String data = "Sets \n"+
152  "i factories /f1*f3/ \n"+
153  "j distribution centers /d1*d5/ \n"+
154  " \n"+
155  "Parameter \n"+
156  "capacity(i) unit capacity at factories \n"+
157  " /f1 500, f2 450, f3 650/ \n"+
158  "demand(j) unit demand at distribution centers \n"+
159  " /d1 160, d2 120, d3 270, d4 325, d5 700 / \n"+
160  "prodcost unit production cost /14/ \n"+
161  "price sales price /24/ \n"+
162  "wastecost cost of removal of overstocked products /4/ \n"+
163  " \n"+
164  "Table transcost(i,j) unit transportation cost \n"+
165  " d1 d2 d3 d4 d5 \n"+
166  " f1 2.49 5.21 3.76 4.85 2.07 \n"+
167  " f2 1.46 2.54 1.83 1.86 4.76 \n"+
168  " f3 3.26 3.08 2.60 3.76 4.45; \n"+
169  " \n"+
170  "$ifthen not set useBig \n"+
171  "Set \n"+
172  " s scenarios /lo,mid,hi/ \n"+
173  " \n"+
174  "Table ScenarioData(s,*) possible outcomes for demand plus probabilities \n"+
175  " d1 d2 d3 d4 d5 prob \n"+
176  " lo 150 100 250 300 600 0.25 \n"+
177  " mid 160 120 270 325 700 0.50 \n"+
178  " hi 170 135 300 350 800 0.25; \n"+
179  "$else \n"+
180  "$if not set nrScen $set nrScen 10 \n"+
181  "Set s scenarios /s1*s%nrScen%/;\n"+
182  "parameter ScenarioData(s,*) possible outcomes for demand plus probabilities;\n"+
183  "option seed=1234; \n"+
184  "ScenarioData(s,'prob') = 1/card(s); \n"+
185  "ScenarioData(s,j) = demand(j)*uniform(0.6,1.4); \n"+
186  "$endif \n"+
187  " \n";
188 
189  static String masterModel = "Sets \n"+
190  "i factories \n"+
191  "j distribution centers \n"+
192  " \n"+
193  "Parameter \n"+
194  " capacity(i) unit capacity at factories \n"+
195  " prodcost unit production cost \n"+
196  " transcost(i,j) unit transportation cost \n"+
197  " \n"+
198  "$if not set datain $abort 'datain not set' \n"+
199  "$gdxin %datain% \n"+
200  "$load i j capacity prodcost transcost \n"+
201  " \n"+
202  "* Benders master problem \n"+
203  "$if not set maxiter $set maxiter 25 \n"+
204  "Set \n"+
205  " iter max Benders iterations /1*%maxiter%/ \n"+
206  " \n"+
207  "Parameter \n"+
208  " cutconst(iter) constants in optimality cuts \n"+
209  " cutcoeff(iter,j) coefficients in optimality cuts \n"+
210  " \n"+
211  "Variables \n"+
212  " ship(i,j) shipments \n"+
213  " product(i) production \n"+
214  " received(j) quantity sent to market \n"+
215  " zmaster objective variable of master problem \n"+
216  " theta future profit \n"+
217  "Positive Variables ship; \n"+
218  " \n"+
219  "Equations \n"+
220  " masterobj master objective function \n"+
221  " production(i) calculate production in each factory \n"+
222  " receive(j) calculate quantity to be send to markets \n"+
223  " optcut(iter) Benders optimality cuts; \n"+
224  " \n"+
225  "masterobj.. \n"+
226  " zmaster =e= theta -sum((i,j), transcost(i,j)*ship(i,j)) \n"+
227  " - sum(i,prodcost*product(i)); \n"+
228  " \n"+
229  "receive(j).. received(j) =e= sum(i, ship(i,j)); \n"+
230  " \n"+
231  "production(i).. product(i) =e= sum(j, ship(i,j)); \n"+
232  "product.up(i) = capacity(i); \n"+
233  " \n"+
234  "optcut(iter).. theta =l= cutconst(iter) + \n"+
235  " sum(j, cutcoeff(iter,j)*received(j)); \n"+
236  " \n"+
237  "model masterproblem /all/; \n"+
238  " \n"+
239  "* Initialize cut to be non-binding \n"+
240  "cutconst(iter) = 1e15; \n"+
241  "cutcoeff(iter,j) = eps; \n"+
242  " \n";
243 
244  static String subModel = "Sets \n"+
245  " i factories \n"+
246  " j distribution centers \n"+
247  " \n"+
248  "Parameter \n"+
249  " demand(j) unit demand at distribution centers \n"+
250  " price sales price \n"+
251  " wastecost cost of removal of overstocked products \n"+
252  " received(j) first stage decision units received \n"+
253  " \n"+
254  "$if not set datain $abort 'datain not set' \n"+
255  "$gdxin %datain% \n"+
256  "$load i j demand price wastecost \n"+
257  " \n"+
258  "* Benders' subproblem \n"+
259  " \n"+
260  "Variables \n"+
261  " sales(j) sales (actually sold) \n"+
262  " waste(j) overstocked products \n"+
263  " zsub objective variable of sub problem \n"+
264  "Positive variables sales, waste \n"+
265  " \n"+
266  "Equations \n"+
267  " subobj subproblem objective function \n"+
268  " selling(j) part of received is sold \n"+
269  " market(j) upperbound on sales \n"+
270  "; \n"+
271  " \n"+
272  "subobj.. \n"+
273  " zsub =e= sum(j, price*sales(j)) - sum(j, wastecost*waste(j)); \n"+
274  " \n"+
275  "selling(j).. sales(j) + waste(j) =e= received(j); \n"+
276  " \n"+
277  "market(j).. sales(j) =l= demand(j); \n"+
278  " \n"+
279  "model subproblem /subobj,selling,market/; \n"+
280  " \n"+
281  "* Initialize received \n"+
282  "received(j) = demand(j); \n"+
283  " \n";
284 }
GAMSEquation getEquation(String identifier)
GAMSGlobals.ModelStat getModelStatus()
GAMSParameter getParameter(String identifier)
void defines(String defStr, String asStr)
void setSystemDirectory(String directory)
This example demonstrates a sequential implementation of a simple Benders decomposition method for a ...
GAMSParameter addParameter(String identifier, int dimension)
GAMSVariable addVariable(String identifier, int dimension, GAMSGlobals.VarType varType)
GAMSSet getSet(String identifier)
GAMSJob addJobFromString(String source)
void instantiate(String modelDefinition, GAMSModifier ... modifiers)
GAMSCheckpoint addCheckpoint()
GAMSModelInstance addModelInstance()
void setWorkingDirectory(String directory)
GAMSDatabase OutDB()
GAMSVariable getVariable(String identifier)
T addRecord(Vector< String > keys)
void setAllModelTypes(String value)