$ontext
Filename: LMP1.gms
Author: Nick Sahinidis, August 2002

Purpose:
Generate and solve random linear multiplicative models of
"Type 1."  Problem instances are generated as proposed by:
   H. Konno and T. Kuno, "Linear multiplicative programming,"
   Mathematical Programming, 56(51-64), 1992.
$offtext

options optcr=0, optca=1.e-6,
    limrow=0, limcol=0,
    solprint=off,
    reslim = 10000;

sets m constraints /1*220/
     n variables   /1*200/
     p products    /1*5/
     c cases       /1*10/
     i instances   /1*5/ ;

*for each case to be solved, we have a
*different (m,n,p) triplet
table cases(c,*)
   1   2   3
1  20  30  2
2  120 100 2
3  220 200 2
4  20  30  3
5  120 120 3
6  200 180 3
7  20  30  4
8  100 100 4
9  200 200 4
10 200 200 5 ;

parameters cc(p,n)  cost coefficients
           A(m,n)   constraint coefficients
           b(m)     left-hand-side
           rep(c,*) summary report ;

parameters mactual, nactual, pactual,
           ResMin, Resmax, NodMin, Nodmax;

variables y(p), x(n), obj ;
x.lo(n) = 0;

equations objective, constraints(m), products(p);

objective .. obj =E= prod(p $ (ord(p) le pactual), y(p));

products(p) $ (ord(p) le pactual) ..
         y(p) =E= sum(n $ (ord(n) le nactual), cc(p,n)*x(n));

constraints(m) $ (ord(m) le mactual) ..
         b(m) =L= sum(n $ (ord(n) le nactual), A(m,n)*x(n)) ;

model lmp1 /all/;
lmp1.workspace = 32;

rep(c,'AvgResUsd') = 0;
rep(c,'AvgNodUsd')= 0;
loop (c,

   mactual = cases(c,'1');
   nactual = cases(c,'2');
   pactual = cases(c,'3');

   ResMin  = inf;
   Resmax  = 0;
   NodMin = inf;
   Nodmax = 0;

   loop(i,

      cc(p,n) = uniform(0,100) $ (ord(p) le pactual and ord(n) le nactual);
      A(m,n)  = uniform(0,100) $ (ord(m) le mactual and ord(n) le nactual);
      b(m)    = uniform(0,100) $ (ord(m) le mactual);

*     make sure all problems have the same zero starting point
      x.l(n)=0; y.l(p)=0;
      solve lmp1 minimizing obj using nlp;
      rep(c,'AvgResUsd')  = rep(c,'AvgResUsd') + lmp1.resusd;
      rep(c,'AvgNodUsd') = rep(c,'AvgNodUsd') + lmp1.nodusd;
      ResMin = min(ResMin, lmp1.resusd);
      NodMin = min(NodMin, lmp1.nodusd);
      ResMax = max(ResMax, lmp1.resusd);
      NodMax = max(NodMax, lmp1.nodusd);

   );

   rep(c,'MinResUsd') = ResMin;
   rep(c,'MaxResUsd') = ResMax;
   rep(c,'MinNodUsd')= NodMin;
   rep(c,'MaxNodUsd')= NodMax;

);

rep(c,'AvgResUsd') = rep(c,'AvgResUsd')/card(i);
rep(c,'AvgNodUsd')= rep(c,'AvgNodUsd')/card(i);

display rep;