$title Simple GUSS Example (GUSSEX1,SEQ=398)

$onText
GUSS is a procedure that facilitates optimization of a collection of related
problems modifying data in those problems. GUSS allows one to solve the
collection of models in a single pass without needing repeated solves or a
LOOP over multiple solves. GUSS is not really a solver but rather organizes and
passes data to and from another LP, NLP, MIP etc solver. This is all done in a
faster fashion than say when using multiple solves through the GAMS Loop command.

In particular GUSS runs the model repeatedly over user specified data for model
parameters that collectively define alternative scenarios to be run. In doing
this it repeatedly updates the base model with the altered scenario data, then
solves the updated model for that scenario and saves user chosen results for
each scenario.

GUSS was developed by Michael R. Bussieck, Michael C. Ferris, and Timo Lohmann.
It is documented in /modlib/adddocs/gusspaper.pdf and in the
GUSS section of the solver manual.

Use of GUSS for an existing model requires six steps

1. Definition of scenarios to run
2. Definition of parameters holding scenario specific data for the items
   in the model that to be changed
3. Definition of parameters that will hold scenario specific model results
   for the items that the user wished to save
4. Definition of a set that tells GUSS the scenarios to run, data to change
   and results to save
5. Modification of the solve statement to identify that scenarios will be run
6. Development of code to report the scenario results

These 6 steps are demonstrated used the simple transport model

Bussieck, M R, Ferris, M C, and Lohmann, T, GUSS: Solving Collections of Data
Related Models within GAMS. In Kallrath, J, Ed, Algebraic Modeling Systems: Modeling
and Solving Real World Optimization Problems. Springer, Berlin Heidelberg, 2012, pp. 35-56.

Keywords: linear programming, transportation problem, scheduling, GUSS, scenario analysis
$offText

Set
   i 'canning plants' / seattle,  san-diego /
   j 'markets'        / new-york, chicago, topeka /;

Parameter
   a(i) 'capacity of plant i in cases'
        / seattle    350
          san-diego  600 /

   b(j) 'demand at market j in cases'
        / new-york   325
          chicago    300
          topeka     275 /;

Table d(i,j) 'distance in thousands of miles'
              new-york  chicago  topeka
   seattle         2.5      1.7     1.8
   san-diego       2.4              1.4;

Scalar f 'freight in dollars per case per thousand miles' / 90 /;

Parameter c(i,j) 'transport cost in thousands of dollars per case';
c(i,j) = f*d(i,j)/1000;

Variable
   x(i,j) 'shipment quantities in cases'
   z      'total transportation costs in thousands of dollars';

Positive Variable x;

Equation
   cost      'define objective function'
   supply(i) 'observe supply limit at plant i'
   demand(j) 'satisfy demand at market j';

cost..      z =e= sum((i,j), c(i,j)*x(i,j));

supply(i).. sum(j, x(i,j)) =l= a(i);

demand(j).. sum(i, x(i,j)) =g= b(j);

Model transport / all /;

Set ScenariosToRun 'scenarios' / base, run1, run2 /;

Table newsupply(ScenariosToRun,i) 'updater for a'
         seattle  san-diego
   base      350        600
   run1      300        650
   run2      401        549;

Table newdemand(ScenariosToRun,j) 'updater for b'
         new-york  chicago  topeka
   base       325      300     275
   run1       325      300     275
   run2       350      300     250;

Set mattrib / system.GUSSModelAttributes /;

Parameter
   xl(ScenariosToRun,i,j)         'collector for level of x'
   xm(ScenariosToRun,i,j)         'collector for marginal of x'
   demandl(ScenariosToRun,j)      'collector for level of demand'
   demandm(ScenariosToRun,j)      'collector for marginal of demand'
   srep(ScenariosToRun, mattrib)  'model attibutes like modelstat etc'
   o(*)                           'GUSS options' / SkipBaseCase 1 /;

Set dict / ScenariosToRun.scenario.''
           o.             opt     .srep
           a.             param   .newsupply
           b.             param   .newdemand
           x.             level   .xl
           x.             marginal.xm 
           demand.        level   .demandl
           demand.        marginal.demandm /;

solve transport using lp minimizing z scenario dict;

* The solve statement from the line above results can also be expressed
* in traditional GAMS syntax as follows:
$onText
loop(ScenariosToRun,
   a(i) = newsupply(ScenariosToRun,i);
   b(j) = newdemand(ScenariosToRun,j);
   solve transport using lp minimizing z;
   xl(ScenariosToRun,i,j) = x.l(i,j);
   xm(ScenariosToRun,i,j) = x.m(i,j);
   demandl(ScenariosToRun,j) = demand.l(j);
   demandm(ScenariosToRun,j) = demand.m(j);
);
$offText

Set Error(scenariostorun) 'empty solution';
Error(scenariostorun) = sum((i,j), xl(scenariostorun,i,j)) = 0;

abort$card(error) 'Missing solution for some scenarios', error;

* Now rerun GUSS with scaleOpt=1 and compare solution
Parameter
   sxl(ScenariosToRun,i,j)         'collector for level of x'
   sxm(ScenariosToRun,i,j)         'collector for marginal of x'
   sdemandl(ScenariosToRun,j)      'collector for level of demand'
   sdemandm(ScenariosToRun,j)      'collector for marginal of demand';

Set sdict / ScenariosToRun.scenario.''
           o.             opt     .srep
           a.             param   .newsupply
           b.             param   .newdemand
           x.             level   .sxl
           x.             marginal.sxm 
           demand.        level   .sdemandl
           demand.        marginal.sdemandm /;

x.scale(i,j) = 10;
demand.scale(j) = 100;
transport.scaleOpt = 1;

solve transport using lp minimizing z scenario sdict;

abort$(smax((ScenariosToRun,i,j), abs(xl(ScenariosToRun,i,j)-sxl(ScenariosToRun,i,j))) > 1e-3) 'Wrong x.l scaled and unscaled:', xl, sxl;
abort$(smax((ScenariosToRun,i,j), abs(xm(ScenariosToRun,i,j)-sxm(ScenariosToRun,i,j))) > 1e-3) 'Wrong x.m scaled and unscaled:', xm, sxm;
abort$(smax((ScenariosToRun,j), abs(demandl(ScenariosToRun,j)-sdemandl(ScenariosToRun,j))) > 1e-3) 'Wrong demand.l scaled and unscaled:', demandl, sdemandl;
abort$(smax((ScenariosToRun,j), abs(demandm(ScenariosToRun,j)-sdemandm(ScenariosToRun,j))) > 1e-3) 'Wrong demand.m scaled and unscaled:', demandm, sdemandm;