Benders2Stage.java
1package com.gams.examples.benders;
2 
3 import java.io.File;
4 
5 import com.gams.api.GAMSCheckpoint;
6 import com.gams.api.GAMSException;
7 import com.gams.api.GAMSGlobals;
8 import com.gams.api.GAMSJob;
9 import com.gams.api.GAMSModelInstance;
10 import com.gams.api.GAMSModifier;
11 import com.gams.api.GAMSOptions;
12 import com.gams.api.GAMSParameter;
13 import com.gams.api.GAMSParameterRecord;
14 import com.gams.api.GAMSSetRecord;
15 import com.gams.api.GAMSVariable;
16 import com.gams.api.GAMSVariableRecord;
17 import com.gams.api.GAMSWorkspace;
18 import com.gams.api.GAMSWorkspaceInfo;
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 
101  masteri.solve(GAMSModelInstance.SymbolUpdateType.BASECASE);
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 
123  subi.solve(GAMSModelInstance.SymbolUpdateType.BASECASE);
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 }
com::gams::api::GAMSDatabase::getEquation
GAMSEquation getEquation(String identifier)
com::gams::api::GAMSModelInstance::getModelStatus
GAMSGlobals.ModelStat getModelStatus()
com::gams::api::GAMSDatabase::getParameter
GAMSParameter getParameter(String identifier)
com::gams::api::GAMSOptions::defines
void defines(String defStr, String asStr)
com::gams::api::GAMSWorkspaceInfo::setSystemDirectory
void setSystemDirectory(String directory)
com::gams.examples.benders.Benders2Stage
This example demonstrates a sequential implementation of a simple Benders decomposition method for a ...
Definition: Benders2Stage.java:35
com::gams::api::GAMSDatabase::addParameter
GAMSParameter addParameter(String identifier, int dimension)
com::gams::api::GAMSDatabase::addVariable
GAMSVariable addVariable(String identifier, int dimension, GAMSGlobals.VarType varType)
com::gams::api::GAMSDatabase::getSet
GAMSSet getSet(String identifier)
com::gams::api::GAMSWorkspace::addJobFromString
GAMSJob addJobFromString(String source)
com::gams::api::GAMSModelInstance::instantiate
void instantiate(String modelDefinition, GAMSModifier ... modifiers)
com::gams::api::GAMSWorkspace::addCheckpoint
GAMSCheckpoint addCheckpoint()
com::gams::api::GAMSCheckpoint::addModelInstance
GAMSModelInstance addModelInstance()
com::gams::api::GAMSWorkspaceInfo::setWorkingDirectory
void setWorkingDirectory(String directory)
com::gams::api::GAMSJob::OutDB
GAMSDatabase OutDB()
com::gams::api::GAMSDatabase::getVariable
GAMSVariable getVariable(String identifier)
GAMSSymbol< GAMSParameterRecord >::addRecord
T addRecord(Vector< String > keys)
com::gams::api::GAMSOptions::setAllModelTypes
void setAllModelTypes(String value)