Benders2StageMT.java
1package com.gams.examples.benders;
2
3import java.io.File;
4import java.util.HashMap;
5import java.util.LinkedList;
6import java.util.List;
7import java.util.Map;
8import java.util.Map.Entry;
9
12import com.gams.api.GAMSGlobals;
13import com.gams.api.GAMSJob;
16import com.gams.api.GAMSOptions;
23
39@SuppressWarnings("unchecked")
40public class Benders2StageMT {
41
42 public static void main(String[] args) {
44 if (args.length > 0)
45 wsInfo.setSystemDirectory( args[0] );
46
47 File workingDirectory = new File(System.getProperty("user.dir"), "Benders2StageMT");
48 workingDirectory.mkdir();
49 wsInfo.setWorkingDirectory(workingDirectory.getAbsolutePath());
50
51 GAMSWorkspace ws = new GAMSWorkspace( wsInfo );
52
53 GAMSJob dataJob = ws.addJobFromString(data);
54
55 GAMSOptions optData = ws.addOptions();
56 optData.defines("useBig", "1");
57 optData.defines("nrScen", "100");
58
59 dataJob.run(optData);
60
61 optData.dispose();
62
63 GAMSParameter scenarioData = dataJob.OutDB().getParameter("ScenarioData");
64
65 GAMSOptions opt = ws.addOptions();
66 opt.defines("datain", dataJob.OutDB().getName());
67 int maxiter = 40;
68 opt.defines("maxiter", Integer.toString(maxiter));
69 opt.setAllModelTypes( "cplexd" );
70
71 GAMSCheckpoint cpMaster = ws.addCheckpoint();
72 GAMSCheckpoint cpSub = ws.addCheckpoint();
73
74 ws.addJobFromString(masterModel).run(opt, cpMaster, dataJob.OutDB());
75
76 GAMSModelInstance masteri = cpMaster.addModelInstance();
77 GAMSParameter cutconst = masteri.SyncDB().addParameter("cutconst", 1, "Benders optimality cut constant");
78 GAMSParameter cutcoeff = masteri.SyncDB().addParameter("cutcoeff", 2, "Benders optimality coefficients");
79 GAMSVariable theta = masteri.SyncDB().addVariable("theta", 0, GAMSGlobals.VarType.FREE, "Future profit function variable");
80 GAMSParameter thetaFix = masteri.SyncDB().addParameter("thetaFix", 0, "");
81 masteri.instantiate("masterproblem max zmaster using lp", opt,
82 new GAMSModifier[] { new GAMSModifier(cutconst), new GAMSModifier(cutcoeff), new GAMSModifier(theta,GAMSGlobals.UpdateAction.FIXED,thetaFix) }
83 );
84
85 ws.addJobFromString(subModel).run(opt, cpSub, dataJob.OutDB());
86
87
88 int numThreads = 4;
89 GAMSModelInstance[] subi = new GAMSModelInstance[numThreads];
90 LinkedList<Tuple> demQueue = new LinkedList<Tuple>();
91
92 subi[0] = cpSub.addModelInstance();
93 GAMSParameter received = subi[0].SyncDB().addParameter("received", 1, "units received from first stage solution");
94 GAMSParameter demand = subi[0].SyncDB().addParameter("demand", 1, "stochastic demand");
95 subi[0].instantiate("subproblem max zsub using lp", opt, new GAMSModifier[] { new GAMSModifier(received), new GAMSModifier(demand) } );
96
97 for (int i = 1; i < numThreads; i++)
98 {
99 subi[i] = subi[0].copyModelInstance();
100 }
101
102 opt.dispose();
103
104 double lowerbound = Double.NEGATIVE_INFINITY, upperbound = Double.POSITIVE_INFINITY, objmaster = Double.POSITIVE_INFINITY;
105 int iter = 1;
106 do {
107 System.out.println("Iteration: " + iter);
108 // Solve master
109 if (iter == 1) // fix theta for first iteration
110 thetaFix.addRecord().setValue( 0.0 );
111 else
112 thetaFix.clear();
113
115 System.out.println(" Master " + masteri.getModelStatus() + " : obj=" + masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel());
116
117 if (iter > 1)
118 upperbound = masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel();
119
120 objmaster = masteri.SyncDB().getVariable("zmaster").getFirstRecord().getLevel() - theta.getFirstRecord().getLevel();
121
122 for (GAMSSetRecord s : dataJob.OutDB().getSet("s"))
123 {
124 Map<String,Double> demDict = new HashMap<String, Double>();
125 for (GAMSSetRecord j : dataJob.OutDB().getSet("j")) {
126 String[] keys = new String[] { s.getKey(0), j.getKey(0) };
127 demDict.put( j.getKey(0), new Double( scenarioData.findRecord( keys ).getValue()) );
128 }
129 String item1 = s.getKey(0);
130 Double item2 = new Double( scenarioData.findRecord(new String[] { s.getKey(0), "prob" } ).getValue() );
131 Tuple items = new Tuple( item1, item2, demDict );
132 demQueue.add(items);
133 }
134
135 for (int i = 0; i < numThreads; i++)
136 subi[i].SyncDB().getParameter("received").clear();
137
138 for (GAMSVariableRecord r : masteri.SyncDB().getVariable("received"))
139 {
140 cutcoeff.addRecord(new String[] {Integer.toString(iter), r.getKey(0)});
141 for (int i = 0; i < numThreads; i++)
142 {
143 subi[i].SyncDB().getParameter("received").addRecord(r.getKeys()).setValue( r.getLevel() );
144 }
145 }
146
147 cutconst.addRecord(Integer.toString(iter));
148 double objsubsum = 0.0;
149
150 // solve multiple model instances in parallel
151 Object queueMutex = new Object();
152 Object ioMutex = new Object();
153 Wrapper<Double>[] objsub = new Wrapper[numThreads];
154 Object[] coef = new Object[numThreads];
155 Wrapper<Double>[] cons = new Wrapper[numThreads] ;
156
157 for (int i = 0; i < numThreads; i++)
158 {
159 objsub[i] = new Wrapper<Double>(new Double(0.0));
160 cons[i] = new Wrapper<Double>(new Double(0.0));
161 coef[i] = new HashMap<String, Double>();
162 for (GAMSSetRecord j : dataJob.OutDB().getSet("j")) {
163 Map<String, Double> cmap = (Map<String, Double>) coef[i];
164 cmap.put( j.getKey(0), new Double(0.0) );
165 }
166 }
167
168 // Solve of sub problems in parallel
169 Scenario[] sc = new Scenario[numThreads];
170 for (int i = 0; i < numThreads; i++)
171 {
172 sc[i] = new Scenario( i, subi[i], cons[i], (Map<String, Double>) coef[i], demQueue, objsub[i], queueMutex, ioMutex );
173 sc[i].start();
174 }
175
176
177 // Synchronize all sub problems
178 for (int i = 0; i < numThreads; i++)
179 {
180 try {
181 sc[i].join();
182 } catch (InterruptedException e) {
183 e.printStackTrace();
184 }
185 }
186
187 for (int i = 0; i < numThreads; i++)
188 {
189 objsubsum += objsub[i].get().doubleValue();
190 double new_consValue = cutconst.findRecord( Integer.toString(iter) ).getValue() + cons[i].get().doubleValue();
191 cutconst.findRecord( Integer.toString(iter) ).setValue( new_consValue );
192
193 for (GAMSSetRecord j : dataJob.OutDB().getSet("j"))
194 {
195 Map<String, Double> map = (Map<String, Double>) coef[i];
196 String[] keys = new String[] { Integer.toString(iter), j.getKey(0) };
197 double newvalue = cutcoeff.findRecord( keys ).getValue( ) + map.get( j.getKey(0) ).doubleValue();
198 cutcoeff.findRecord( keys ).setValue( newvalue );
199 }
200 }
201
202 lowerbound = Math.max(lowerbound, objmaster + objsubsum);
203
204 iter++;
205 if (iter == maxiter + 1)
206 throw new GAMSException("Benders out of iterations");
207
208 System.out.println(" lowerbound: " + lowerbound + " upperbound: " + upperbound + " objmaster: " + objmaster);
209
210 } while ((upperbound - lowerbound) >= 0.001 * (1 + Math.abs(upperbound)));
211
212 masteri.dispose();
213
214 for(GAMSModelInstance inst : subi)
215 {
216 inst.dispose();
217 }
218 }
219
221 static class Scenario extends Thread {
222 int _i;
223 GAMSModelInstance _submi;
224 Wrapper<Double> _cutconst;
225 List<Tuple> _demQueue;
226 Map<String, Double> _cutcoeff;
227 Wrapper<Double> _objsub;
228 Object _queueMutex;
229 Object _ioMutex;
230
241 public Scenario(int i , GAMSModelInstance submi, Wrapper<Double> cutconst, Map<String, Double> cutcoeff, LinkedList<Tuple> demQueue,
242 Wrapper<Double> objsub, Object queueMutex, Object ioMutex) {
243 _i = i;
244 _submi = submi;
245 _cutconst = cutconst;
246 _cutcoeff = cutcoeff;
247 _demQueue = demQueue;
248 _objsub = objsub;
249 _queueMutex = queueMutex;
250 _ioMutex = ioMutex;
251 }
252
254 public void run() {
255 while (true)
256 {
257 Tuple demDict;
258
259 synchronized (_queueMutex)
260 {
261 if (_demQueue.size() == 0)
262 break;
263 else
264 demDict = _demQueue.remove(0); // dequeue
265 }
266
267 _submi.SyncDB().getParameter("demand").clear();
268 for (Entry<String, Double> kv : demDict.getItem3().entrySet())
269 _submi.SyncDB().getParameter("demand").addRecord(kv.getKey()).setValue( kv.getValue() );
270
272
273
274 synchronized (_ioMutex)
275 {
276 System.out.println(" Thread "+_i+" : Sub " + _submi.getModelStatus().toString() + " : obj=" + _submi.SyncDB().getVariable("zsub").getFirstRecord().getLevel());
277 }
278
279 double probability = demDict.getItem2().doubleValue();
280 double new_objsubValue = _objsub.get().doubleValue() + probability * _submi.SyncDB().getVariable("zsub").getFirstRecord().getLevel();
281 _objsub.set( new Double( new_objsubValue ) );
282
283 for (Entry<String, Double> kv : demDict.getItem3().entrySet())
284 {
285 double new_custconstValue = _cutconst.get().doubleValue() + probability * _submi.SyncDB().getEquation("market").findRecord( kv.getKey() ).getMarginal() * kv.getValue().doubleValue();
286 _cutconst.set( new Double( new_custconstValue ) );
287 double new_cutcoeffValue = _cutcoeff.get( kv.getKey() ).doubleValue() + probability * _submi.SyncDB().getEquation("selling").findRecord(kv.getKey()).getMarginal();
288 _cutcoeff.put( kv.getKey(), new_cutcoeffValue );
289 }
290 }
291 }
292 }
293
294 static String data = "Sets \n"+
295 "i factories /f1*f3/ \n"+
296 "j distribution centers /d1*d5/ \n"+
297 " \n"+
298 "Parameter \n"+
299 "capacity(i) unit capacity at factories \n"+
300 " /f1 500, f2 450, f3 650/ \n"+
301 "demand(j) unit demand at distribution centers \n"+
302 " /d1 160, d2 120, d3 270, d4 325, d5 700 / \n"+
303 "prodcost unit production cost /14/ \n"+
304 "price sales price /24/ \n"+
305 "wastecost cost of removal of overstocked products /4/ \n"+
306 " \n"+
307 "Table transcost(i,j) unit transportation cost \n"+
308 " d1 d2 d3 d4 d5 \n"+
309 " f1 2.49 5.21 3.76 4.85 2.07 \n"+
310 " f2 1.46 2.54 1.83 1.86 4.76 \n"+
311 " f3 3.26 3.08 2.60 3.76 4.45; \n"+
312 " \n"+
313 "$ifthen not set useBig \n"+
314 "Set \n"+
315 " s scenarios /lo,mid,hi/ \n"+
316 " \n"+
317 "Table ScenarioData(s,*) possible outcomes for demand plus probabilities \n"+
318 " d1 d2 d3 d4 d5 prob \n"+
319 " lo 150 100 250 300 600 0.25 \n"+
320 " mid 160 120 270 325 700 0.50 \n"+
321 " hi 170 135 300 350 800 0.25; \n"+
322 "$else \n"+
323 "$if not set nrScen $set nrScen 10 \n"+
324 "Set s scenarios /s1*s%nrScen%/;\n"+
325 "parameter ScenarioData(s,*) possible outcomes for demand plus probabilities;\n"+
326 "option seed=1234; \n"+
327 "ScenarioData(s,'prob') = 1/card(s); \n"+
328 "ScenarioData(s,j) = demand(j)*uniform(0.6,1.4); \n"+
329 "$endif \n"+
330 " \n";
331
332 static String masterModel = "Sets \n"+
333 "i factories \n"+
334 "j distribution centers \n"+
335 " \n"+
336 "Parameter \n"+
337 "capacity(i) unit capacity at factories \n"+
338 "prodcost unit production cost \n"+
339 "transcost(i,j) unit transportation cost \n"+
340 " \n"+
341 "$if not set datain $abort 'datain not set' \n"+
342 "$gdxin %datain% \n"+
343 "$load i j capacity prodcost transcost \n"+
344 " \n"+
345 "* Benders master problem \n"+
346 "$if not set maxiter $set maxiter 25 \n"+
347 "Set \n"+
348 " iter max Benders iterations /1*%maxiter%/ \n"+
349 " \n"+
350 "Parameter \n"+
351 " cutconst(iter) constants in optimality cuts \n"+
352 " cutcoeff(iter,j) coefficients in optimality cuts \n"+
353 " \n"+
354 "Variables \n"+
355 " ship(i,j) shipments \n"+
356 " product(i) production \n"+
357 " received(j) quantity sent to market \n"+
358 " zmaster objective variable of master problem \n"+
359 " theta future profit \n"+
360 "Positive Variables ship; \n"+
361 " \n"+
362 "Equations \n"+
363 " masterobj master objective function \n"+
364 " production(i) calculate production in each factory \n"+
365 " receive(j) calculate quantity to be send to markets \n"+
366 " optcut(iter) Benders optimality cuts; \n"+
367 " \n"+
368 "masterobj.. \n"+
369 " zmaster =e= theta -sum((i,j), transcost(i,j)*ship(i,j)) \n"+
370 " - sum(i,prodcost*product(i)); \n"+
371 " \n"+
372 "receive(j).. received(j) =e= sum(i, ship(i,j)); \n"+
373 " \n"+
374 "production(i).. product(i) =e= sum(j, ship(i,j)); \n"+
375 "product.up(i) = capacity(i); \n"+
376 " \n"+
377 "optcut(iter).. theta =l= cutconst(iter) + \n"+
378 " sum(j, cutcoeff(iter,j)*received(j)); \n"+
379 " \n"+
380 "model masterproblem /all/; \n"+
381 " \n"+
382 "* Initialize cut to be non-binding \n"+
383 "cutconst(iter) = 1e15; \n"+
384 "cutcoeff(iter,j) = eps; \n"+
385 " \n";
386
387 static String subModel = "Sets \n"+
388 " i factories \n"+
389 " j distribution centers \n"+
390 " \n"+
391 "Parameter \n"+
392 " demand(j) unit demand at distribution centers \n"+
393 " price sales price \n"+
394 " wastecost cost of removal of overstocked products \n"+
395 " received(j) first stage decision units received \n"+
396 " \n"+
397 "$if not set datain $abort 'datain not set' \n"+
398 "$gdxin %datain% \n"+
399 "$load i j demand price wastecost \n"+
400 " \n"+
401 "* Benders' subproblem \n"+
402 " \n"+
403 "Variables \n"+
404 " sales(j) sales (actually sold) \n"+
405 " waste(j) overstocked products \n"+
406 " zsub objective variable of sub problem \n"+
407 "Positive variables sales, waste \n"+
408 " \n"+
409 "Equations \n"+
410 " subobj subproblem objective function \n"+
411 " selling(j) part of received is sold \n"+
412 " market(j) upperbound on sales \n"+
413 "; \n"+
414 " \n"+
415 "subobj.. \n"+
416 " zsub =e= sum(j, price*sales(j)) - sum(j, wastecost*waste(j)); \n"+
417 " \n"+
418 "selling(j).. sales(j) + waste(j) =e= received(j); \n"+
419 " \n"+
420 "market(j).. sales(j) =l= demand(j); \n"+
421 " \n"+
422 "model subproblem /subobj,selling,market/; \n"+
423 " \n"+
424 "* Initialize received \n"+
425 "received(j) = demand(j); \n"+
426 " \n";
427}
428
430class Wrapper<T> {
431 T _value;
432 public Wrapper(T value) { _value = value; }
433 public T get() { return _value; }
434 public void set(T anotherValue) { _value = anotherValue; }
435}
436
440class Tuple {
441 String _item1;
442 Double _item2;
443 Map<String, Double> _item3;
444 public Tuple(String item1, Double item2, Map<String, Double> item3) {
445 _item1 = item1;
446 _item2 = item2;
447 _item3 = item3;
448 }
449 public String getItem1() { return _item1; }
450 public Double getItem2() { return _item2; }
451 public Map<String, Double> getItem3() { return _item3; }
452}
GAMSCheckpoint addCheckpoint()
Provides package namespace for Java interface and examples to General Algebraic Model System (GAMS)...
GAMSModelInstance addModelInstance()
GAMSVariable getVariable(String identifier)
GAMSParameter addParameter(String identifier, int dimension)
void instantiate(String modelDefinition, GAMSModifier ... modifiers)
&#160;
GAMSEquation getEquation(String identifier)
GAMSParameter getParameter(String identifier)
void defines(String defStr, String asStr)
void setSystemDirectory(String directory)
T addRecord(Vector< String > keys)
GAMSDatabase OutDB()
GAMSJob addJobFromString(String source)
void setAllModelTypes(String value)
GAMSVariable addVariable(String identifier, int dimension, GAMSGlobals.VarType varType)
GAMSSet getSet(String identifier)
This example demonstrates a parallel implementation of a simple Benders decomposition method for a st...
void setWorkingDirectory(String directory)
GAMSGlobals.ModelStat getModelStatus()