* Trim loss minimization problem 
* Westerlund and coworkers
* Note that the integer variables are considered continuous
* and modeled through the introduction of binary variables
* and appropriate constraints.

OPTION ITERLIM=100000;

* Sets
SET i roll type /1*4/;
SET j pattern number /1*4/;
SET k /1*3/;
SET p /1*5/;

ALIAS(j,jj);

* Parameters
SCALAR 	nj   / 3 /;

PARAMETER n(i)  order size for each paper roll
         / '1'     15
           '2'     28
           '3'     21
           '4'     30/;

PARAMETER b(i)  width of paper roll
         / '1'     290
           '2'     315
           '3'     350
           '4'     455/;

PARAMETER mupp(j)  upper bound on number of repeats for pattern j
         / '1'     30
           '2'     30
           '3'     30
           '4'     30/;

PARAMETER rs starting point for r;

TABLE rs(i,j) 
		1	2	3	4
	1 	1 	0 	2 	0
	2 	2	1 	1 	0
	3 	0 	3 	0 	0
	4	2 	1 	2 	0 ;

VARIABLES
	 r(i,j)    number of rolls of type i in pattern j
	 rr(i,j,k) binary variables representing r
	 y(j)      existence of pattern j 
	 m(j)      repeats of pattern j
	 mm(j,p)   binary variables representing m
     objval    objective function variable;

FREE VARIABLES    objval;
BINARY VARIABLE y, rr, mm;
POSITIVE VARIABLE r, m;


EQUATIONS
	 f Objective function
	 numroll(i)  order constraints ensuring sufficient production
	 widthL(j)   width lower bound constraint
	 widthU(j)   width upper bound constraint
	 rL(j)       logical constraint on r
	 sumr(j)     logical constraint on r
	 mL(j)       logical constraint on m
	 mU(j)       logical constraint on m
	 sumbil      lower bound on total number of patterns made
	 yy(j)       ordering of y variables to reduce degeneracy
	 lmm(j)      ordering of m variables to reduce degeneracy
	 grr(i,j)    representation of r through binary variables
	 gmm(j)      representation of r through binary variables;

f  .. objval =e=SUM(j,m(j)) + 0.1 * y('1') + 0.2 * y('2') + 0.3 * y('3') + 0.4 * y('4');
numroll(i)  .. SUM(j,m(j) * r(i,j)) =g= n(i);
widthL(j)  .. SUM(i,-b(i) * r(i,j)) + 1750 * y(j) =l= 0;
widthU(j)  .. SUM(i,b(i) * r(i,j)) - 1850 * y(j) =l= 0;
rL(j)  .. y(j) - SUM(i,r(i,j)) =l= 0;
sumr(j)  .. SUM(i,r(i,j)) - 5 * y(j) =l= 0;
mL(j)  .. y(j) - m(j) =l= 0;
mU(j)  .. m(j) - mupp(j) * y(j) =l= 0;
sumbil  .. SUM(j,m(j)) =g= 19;
yy(j)$(ORD(j) le nj) .. SUM(jj$(ORD(jj) eq ORD(j)+1),y(jj)) - y(j) =l= 0;
lmm(j)$(ORD(j) le nj)  .. SUM(jj$(ORD(jj) eq ORD(j)+1),m(jj)) - m(j) =l= 0;

grr(i,j) .. r(i,j) =e= rr(i,j,'1') + 2*rr(i,j,'2') + 4*rr(i,j,'3');
gmm(j) .. m(j) =e= mm(j,'1') + 2*mm(j,'2') + 4*mm(j,'3') + 8*mm(j,'4') + 16*mm(j,'5');

* Bounds
r.LO(i,j)  = 0;
r.UP(i,j)  = 5;
m.LO(j)  = 0;
m.UP(j)  = mupp(j);

* Starting point (global solution)
* y.L('1') = 1; y.L('2') = 1; y.L('3') = 1; y.L('4') = 0;
* m.L('1') = 9; m.L('2') = 7; m.L('3') = 3; m.L('4') = 0;
* r.L(i,j) = rs(i,j);

MODEL test /ALL/;
SOLVE test USING MINLP MINIMIZING objval;