Resources

Details Examples (omod0e)


/***************************************************************/
/*                                                             */
/*          S A S   S A M P L E   L I B R A R Y                */
/*                                                             */
/*    NAME: omod0e                                             */
/*   TITLE: Details Examples (omod0e)                          */
/* PRODUCT: OR                                                 */
/*  SYSTEM: ALL                                                */
/*    KEYS: OR                                                 */
/*   PROCS: OPTMODEL                                           */
/*    DATA:                                                    */
/*                                                             */
/* SUPPORT:                             UPDATE:                */
/*     REF:                                                    */
/*    MISC: Examples from the Details section of the OPTMODEL  */
/*          chapter of Mathematical Programming.               */
/*                                                             */
/***************************************************************/

/*  Named Parameters  */
data coeff;
   input c_xx c_x c_y c_xy c_yy;
   datalines;
   1 -1  -2 -1 1
   ;
proc optmodel;
   var x, y;
   number c_xx, c_x, c_y, c_xy, c_yy;
   read data coeff into c_xx c_x c_y c_xy c_yy;
   min z=c_xx*x**2 + c_x*x + c_y*y + c_xy*x*y + c_yy*y**2;
   solve;

/*  Indexing  */
proc optmodel;
   number p{1..3};
   p[1]=5;
   p[2]=7;
   p[3]=9;
   put p[*]=;
quit;

proc optmodel;
   number p{i in 1..3} init 3 + 2*i;
   put (sum{i in 1..3} p[i]);
quit;

proc optmodel;
   number n init 100000;
   var x{1..n};
   min z = sum{i in 1..n}(x[i] - log(i))**2;
   solve;
quit;

/*  Parameters  */
proc optmodel;
   number pi=4*atan(1);
   number r;
   number circum=2*pi*r;
   r=1;
   put circum;       /* prints 6.2831853072 */
   r=2;
   put circum;       /* prints 12.566370614 */
quit;

/*  Expressions  */
proc optmodel;
   var x init 0.5 >= 0 <= 1;
   var y init (0.5 >= 0) <= 1;
   print x y;
   expand;
quit;

proc optmodel;
   number original{i in 1..8} = sin(i);
   number sorted{i in 1..8} init original[i];
   call sortn(of sorted[*]);
   print original sorted;

/*  Index Sets  */

proc optmodel;
   put (sum{i in 1..10} i**2);
quit;

/*  AND Aggregation Expression  */
proc optmodel;
   put (and{i in 1..5} i < 10); /* returns 1 */
   put (and{i in 1..5} i NE 3); /* returns 0 */
quit;

/*  CARD Function  */
proc optmodel;
   put (card(1..3));
quit;

/*  CROSS Expression  */
proc optmodel;
   set s1 = 1..2;
   set<string> s2 = {'a', 'b'};
   set<number, string> s3=s1 cross s2;
   put 's3 is ' s3;
   set<number, string, number> s4 = s3 cross 4..5;
   put 's4 is ' s4;
quit;

/*  IN Expression  */
proc optmodel;
   set s = 1..10;
   put (5 in s);         /* outputs 1 */
   put (-1 not in s);    /* outputs 1 */
   set<num, str> t = {<1,'a'>, <2,'b'>, <2,'c'>};
   put (<2, 'b'> in t);   /* outputs 1 */
   put (<1, 'b'> in t);  /* outputs 0 */
quit;

/*  Index Set Expression  */
proc optmodel;
   put ({i in 1..5 : i NE 3}); /* outputs {1,2,4,5} */
quit;

/*  INTER Expression  */

proc optmodel;
   put ({1,3} inter {2,3}); /* outputs {3} */
quit;

/*  INTER Aggregation Expression  */

proc optmodel;
   put (inter{i in 1..3} i..i+3); /* outputs {3,4} */
quit;

/*  MAX Aggregation Expression  */

proc optmodel;
   put (max{i in 2..5} 1/i);
quit;

/*  MIN aggregation Expression  */

proc optmodel;
   put (min{i in 2..5} 1/i);
quit;

/*  OR Aggregation Expression  */

proc optmodel;
   put (or{i in 1..5} i = 2); /* returns 1 */
   put (or{i in 1..5} i = 7); /* returns 0 */
quit;

/*  PROD Aggregation Expression  */

proc optmodel;
   number n = 5;
   put (prod{i in 1..n} i); /* outputs 120 */
quit;

/*  RANGE Expression  */

proc optmodel;
   put (10..30 by 7); /* outputs {10,17,24} */
quit;

/*  Set Constructor Expression  */

proc optmodel;
   put ({1,2,3,2}); /* outputs {1,2,3} */
quit;

proc optmodel;
   number m = 3, n = 4;
   var x{1..4} init 1;
   string y = 'c';
   put ({<'a', x[3]>, <'b', m>, <y, m/n>});
quit;

/*  SETOF Aggregation Expression  */

proc optmodel;
   put (setof{i in 1..3}<i, i*i, i**3>);
quit;

/*  SLICE Expression  */

proc optmodel;
   put (slice(<1,*>, {<1,3>, <1,0>, <3,1>}));
   put (slice(<*,2,*>, {<1,2,3>, <2,4,3>, <2,2,5>}));
quit;

proc optmodel;
   set<str,str> dep = {<'B','A'>, <'C','B'>, <'D','C'>};
   set<str,str> cl;
   set<str> cn;
   cl = dep;
   cn =  (setof{<i,j> in dep} i) inter (setof{<i,j> in dep} j);
   for {node in cn}
       cl = cl union (slice(<*,node>,cl) cross slice(<node,*>,cl));
   put cl;
quit;

/*  SUM Expression   */

proc optmodel;
   put (sum {i in 1..10} i);  /* outputs 55 */
quit;

/*  SYMDIFF Expression  */

proc optmodel;
   put ({1,3} symdiff {2,3}); /* outputs {1,2} */
quit;

/*  Tuple Expression  */

proc optmodel;
   put (<1,2,3> in setof{i in 1..2}<i,i+1,i+2>);
   put ({<1,'a'>, <2,'b'>} cross {<3,'c'>, <4,'d'>});
quit;

/*  UNION Expression  */

proc optmodel;
   put ({1,3} union {2,3}); /* outputs {1,3,2} */
quit;

/*  UNION Aggregation Expression  */

proc optmodel;
   put (union{i in 1..3} i..i+3); /* outputs {1,2,3,4,5,6} */
quit;

/*  WITHIN Expression  */

proc optmodel;
   put ({1,3} within {2,3});     /* outputs 0 */
   put ({1,3} not within {2,3}); /* outputs 1 */
   put ({1,3} within {1,2,3});   /* outputs 1 */
quit;

/*  Data Set Input/Output  */
proc optmodel;
   set LOCS = {'New York', 'Washington', 'Boston'};  /* locations */
   set DOW = 1..7;  /* day of week */
   var s{LOCS, DOW} init 1;
   create data soldata from [location day_of_week]={LOCS, DOW} sale=s;
quit;

data abcd;
   input letter $ @@;
   datalines;
a b c d
;

proc optmodel;
   set<num> subscripts=1..4;
   string  letter{subscripts};
   read data abcd into [_N_] letter[5-_N_];
   print letter;
quit;

data revdata;
   input month rev @@;
   datalines;
1 200 2 345 3 362 4 958
5 659 6 804 7 487 8 146
9 683 10 732 11 652 12 469
;

proc optmodel;
   set m = 1..3;
   var revenue{1..12};
   read data revdata into [_N_] revenue=rev;
   create data qtr1 from [month]=m revenue[month];
   create data qtr2 from [month]=m revenue[month+3];
   create data qtr3 from [month]=m revenue[month+6];
   create data qtr4 from [month]=m revenue[month+9];

/*  Formatted Output  */
proc optmodel;
   number a=1.7, b=2.8;
   set s={a,b};
   put a b;        /* list output */
   put a= b=;      /* named output */
   put 'Value A: ' a 8.1 @30 'Value B: ' b 8.; /* formatted */
   string str='Ratio (A/B) is:';
   put str (a/b);  /* strings and expressions */
   put s=;         /* named set output */
quit;

proc optmodel;
   number square{i in 0..5} = i*i;
   number recip{i in 1..5} = 1/i;
   print square recip;

proc optmodel;
   set R=1..6;
   set C=1..4;
   number a{i in R, j in C} = 10*i+j;
   print a;
quit;

/*  ODS Table and Variable Names  */
proc optmodel printlevel=2;
   ods output PrintTable=expt ProblemSummary=exps DerivMethods=exdm
              SolverOptions=exso SolutionSummary=exss OptStatistics=exos
              Timing=exti;
   performance details;
   var x{1..2} >= 0;
   min z = 2*x[1] + 3 * x[2] + x[1]**2 + 10*x[2]**2
           + 2.5*x[1]*x[2] + x[1]**3;
   con c1: x[1] - x[2] <= 1;
   con c2: x[1] + 2*x[2] >= 100;
   solve;
   print x;
quit;

title2 'PrintTable';
proc print data=expt;
run;

title2 'ProblemSummary';
proc print data=exps;
run;

title2 'SolverOptions';
proc print data=exso;
run;

title2 'DerivMethods';
proc print data=exdm;
run;

title2 'SolutionSummary';
proc print data=exss;
run;

title2 'OptStatistics';
proc print data=exos;
run;

title2 'Timing';
proc print data=exti;
run;

title2;

/*  Constraints  */
proc optmodel;
   var x, y;
   min r = x**2 + y**2;
   con c: x+y >= 1;
   solve;
   print x y;
quit;

proc optmodel;
   var x{1..3};
   con b: sum{i in 1..3}(x[i] - i) = 0;
   expand b;
quit;

proc optmodel;
    var x init 1;
    con c1: x**2 <= 5;
    con c2: 5 >= x**2;
    con c3: -x**2 >= -5;
    con c4: -5 <= -x**2;
    expand;
    print c1.body c2.body c3.body c4.body;
quit;

/*  Suffixes  */

proc optmodel;
     var x, y;
     min z = x + y;
     con c: x + 2*y <= 3;
     solve;
     print x.lb x.ub x.init x.sol;
     print y.lb y.ub y.init y.sol;
     print c.lb c.ub c.body c.dual;

     x.lb=0;
     y.lb=0;
     c.lb=1;
     solve;
     print x.lb x.ub x.init x.sol;
     print y.lb y.ub y.init y.sol;
     print c.lb c.ub c.body c.dual;
quit;

/*  Dual Value  */

proc optmodel;
   var x, y;
   min z = 6*x + 7*y;
   con
      4*x +   y >=  5,
       -x - 3*y <= -4,
        x +   y <=  4;
   solve;
   print x y;
   expand _ACON_ ;
   print _ACON_.dual _ACON_.body;

   _ACON_[1].lb = _ACON_[1].lb - 1;
   solve;
   _ACON_[2].ub = _ACON_[2].ub + 1;
   solve;
quit;

/*  Reduced Costs  */

proc optmodel;
   var x >= 0, y >= 0, z >= 0;
   max cost = 4*x + 3*y - 5*z;
   con
      -x + y + 5*z <= 15,
      3*x - 2*y - z <= 12,
      2*x + 4*y + 2*z <= 16;
   solve;
   print x y z;
   print x.rc y.rc z.rc;
quit;

/*  Model Update  */

proc optmodel;
   var x, y;
   min r = x**2 + y**2;
   con c: x+y >= 1;
   solve;
   print x y;

   drop c;
   solve;
   print x y;

   restore c;
   fix x=0.3;
   solve;
   print x y c.dual;

   unfix x;
   solve;
   print x y c.dual;
   quit;

proc optmodel printlevel=0;
   var x1 >= 0, x2 >= 0, x3 >= 0, x4 >= 0, x5 >= 0;

   minimize z = x1 + x2 + x3 + x4;

   con a1: x1 + x2 + x3          <= 4;
   con a2:               x4 + x5 <= 6;
   con a3: x1 +          x4      >= 5;
   con a4:      x2 +          x5 >= 2;
   con a5:           x3          >= 3;

   solve with lp;

   print _var_.name _var_ _var_.rc _var_.status;
   print _con_.name _con_.lb _con_.body _con_.ub _con_.dual _con_.status;

/*  OPTMODEL Options  */
proc optmodel;
   number n = 1/7;
   for {i in 1..9 by 4}
   do;
      reset options pdigits=(i);
      print n;
   end;
   reset options pdigits; /* reset to default */
quit;

/* Calling FCMP Functions */
proc fcmp outlib=work.funcs.test;
   /* arithmetic geometric mean */
   function agm(a0, b0);
      a=a0; b = b0;
      if a<=0 or b<=0 then return(0);
      do until( a - b < a/1e12 );
         a1 = 0.5*a + 0.5*b;
         b1 = sqrt(a*b);
         a = a1; b = b1;
      end;
      return(a);
   endsub;
run;

/* libraries must be specified with the CMPLIB option */
option cmplib=work.funcs;

proc optmodel;
   title 'AGM(1,2)';
   print (agm(1,2));

proc optmodel;
   title 'AGM optimization';
   /* find x where agm(1,x) == 23 */
   var x init 1;
   num c = 23;
   min z=(agm(1,x)-c)^2;
   solve;
   print x;
quit;

/* FCMP Output Arguments */
proc fcmp outlib=work.funcs.test;
   subroutine do_sqr(x, sq, text $);
      outargs sq, text;
      sq = x*x;
      text = 'This is an example of output arguments';
   endsub;
run;

option cmplib=work.funcs;

proc optmodel;
   string s init repeat(' ', 79); /* reserve 80 bytes */
   number n;
   call do_sqr(7, n, s);
   print s n;
quit;

/* Constant Matrix Argument to an FCMP Function */
proc fcmp outlib=work.funcs.test;
   function evalpoly(x, coeff[*]);
      z = 0;
      do i = dim1(coeff) to 1 by -1;
         z = z * x + coeff[i];
      end;
      return (z);
   endsub;
run;

option cmplib=work.funcs;

proc optmodel;
   num coeff{1..3} = [1, -2, 1];
   var x;
   min z=evalpoly(x, coeff);
   solve;
   print x;
quit;

/* OPTMODEL Declarations for FCMP Matrices */
proc optmodel;
   /* the following arrays can be used as matrices */
   num N;
   num mat1{1..N};         /* OK */
   set S1 = 1..5;
   num mat2{S1};           /* OK */
   set S2 = {S1,S1};
   num mat3{S2};           /* OK */
   num mat4{S2 cross S2};  /* OK */

   /* the following arrays cannot be used as matrices */
   num L init 1;
   num arr1{L..N};         /* lower bound is not literally 1 */
   num arr2{1..10 by 3};   /* step size is not 1 */
   set S3 init S1;
   num arr3{S3};           /* S3 is modifiable */
   S3 = {3, 5, 7};
   num arr4{i in 1..N, j in 1..N: j >= i}; /* selection expression used */
   num arr5{i in 1..N, j in 1..i};         /* index dependency on 'i' */
quit;

/*  More on Index Sets  */
proc optmodel;
   put (setof {i in 1..3, j in 1..3 : j >= i} <i, j>);
quit;

proc optmodel;
   put ({i in 1..3, i..3});
quit;

proc optmodel;
   number N = 3;
   set<num,str> S = {<1,'a'>,<2,'b'>,<3,'a'>,<4,'b'>};
   put ({i in 1..N, <(i),j> in S});
   put ({i in 1..N, j in slice(<i,*>, S)});
quit;

/*  PROC NLP Example 7  */
proc nlp all;
   parms amountx amounty amounta amountb amountc
         pooltox pooltoy ctox ctoy pools = 1;
   bounds 0 <= amountx amounty amounta amountb amountc,
               amountx <= 100,
               amounty <= 200,
          0 <= pooltox pooltoy ctox ctoy,
          1 <= pools <= 3;
   lincon amounta + amountb = pooltox + pooltoy,
          pooltox + ctox = amountx,
          pooltoy + ctoy = amounty,
          ctox + ctoy    = amountc;
   nlincon nlc1-nlc2 >= 0.,
           nlc3 = 0.;
   max f;
   costa = 6; costb = 16; costc = 10;
   costx = 9; costy = 15;
   f = costx * amountx + costy * amounty
       - costa * amounta - costb * amountb - costc * amountc;
   nlc1 = 2.5 * amountx - pools * pooltox - 2. * ctox;
   nlc2 = 1.5 * amounty - pools * pooltoy - 2. * ctoy;
   nlc3 = 3 * amounta + amountb - pools * (amounta + amountb);
run;

proc optmodel;
   var amountx init 1 >= 0 <= 100,
       amounty init 1 >= 0 <= 200;
   var amounta init 1 >= 0,
       amountb init 1 >= 0,
       amountc init 1 >= 0;
   var pooltox init 1 >= 0,
       pooltoy init 1 >= 0;
   var ctox init 1 >= 0,
       ctoy init 1 >= 0;
   var pools init 1 >=1 <= 3;
   con amounta + amountb = pooltox + pooltoy,
       pooltox + ctox = amountx,
       pooltoy + ctoy = amounty,
       ctox + ctoy    = amountc;
   number costa, costb, costc, costx, costy;
   costa = 6; costb = 16; costc = 10;
   costx = 9; costy = 15;
   max f = costx * amountx + costy * amounty
          - costa * amounta - costb * amountb - costc * amountc;
   con nlc1: 2.5 * amountx - pools * pooltox - 2. * ctox >= 0,
       nlc2: 1.5 * amounty - pools * pooltoy - 2. * ctoy >= 0,
       nlc3: 3 * amounta + amountb - pools * (amounta + amountb)
              = 0;
   solve;
   print amountx amounty amounta amountb amountc
         pooltox pooltoy ctox ctoy pools;