Longest flow shop (2)

Permutations in GAMS

This is a follow up to http://yetanothermathprogrammingconsultant.blogspot.com/2012/04/longest-flow-shop.html.

We consider again the small data set:

table proctime(m,j) 'processing times'             job1  job2  job3  job4  job5  job6  job7  job8  job9  job10   machine1   4     2     1     5     4     3     5     2     1     8   machine2   2     5     8     6     7     4     7     3     6     2   machine3   2     8     6     2     2     4     2     7     1     8 ;

This data set comes from:

The first approach to find the permutation with the longest schedule is to enumerate all permutations. GAMS has a built-in option to generate permutations of a set. It is not completely intuitive, in part because the way GAMS stores sets. All sets have fixed ordering with respect to a pool of set elements (the universe). To create a different ordering, an extra index is used. E.g. for 3 jobs we have:

$set n 3 $eval k fact(%n%) set    j           'jobs'                    /job1*job%n%/    p           'permutations'            /p1*p%k%/    allp(p,j,j) 'store all permutations'    m           'processing stage'        /machine1*machine3/ ; option allp>j; option allp:0:1:1; display allp;

The option > creates the permutations. Of course the syntax is not very readable: one would expect a function with a name somehow related to ‘permutation’. One of the major goals when using a modeling system is to produce readable models, and that is certainly not achieved with this syntax.

The result is:

----     13 SET allp  store all permutations INDEX 1 = p1             job1        job2        job3 job1         YES job2                     YES job3                                 YES INDEX 1 = p2             job1        job2        job3 job1         YES job2                                 YES job3                     YES INDEX 1 = p3             job1        job2        job3 job1                     YES job2         YES job3                                 YES INDEX 1 = p4             job1        job2        job3 job1                     YES job2                                 YES job3         YES INDEX 1 = p5             job1        job2        job3 job1                                 YES job2         YES job3                     YES INDEX 1 = p6             job1        job2        job3 job1                                 YES job2                     YES job3         YES

Although I don’t like the syntax, it is actually quite fast. If we set n to 10 and disable the display statement – it would take a long time to print all these permutations – we can generate all 10!=3,628,800 permutations in less than 2.5 seconds. This compares quite good to the following R session:

We can create and evaluate the schedules as follows:

loop(p,    putclose$(mod(ord(p),1000) = 0) "iteration " ord(p):0:0/;    completion(m,j) = 0;    loop((m,allp(p,jj,j)),       completion(m,jj) = max(completion(m,jj-1),completion(m-1,jj)) + proctime(m,j);    );    obj(p) = smax((m,j),completion(m,j)); );

Unfortunately this loop is very slow, it took almost 7 hours on my machine.

Results

The following picture shows the three schedules we discussed:

• optimal
• initial schedule: formed by sequence job1,job2,job3,…
• a longest schedule

A MIP Formulation

Of course with a little thought we can improve on this solution technique and write a MIP model:

* * assignment constraints * assign1(j).. sum(jj, position(j,jj)) =e= 1; assign2(j).. sum(jj, position(jj,j)) =e= 1; processing(m,j).. completion(m,j) =e= start(m,j) + sum(jj,position(jj,j)*proctime(m,jj)); * * start(m,j) = max(completion(m,j-1),completion(m-1,j)) * max1(m,j)$(ord(j)>1).. start(m,j) =g= completion(m,j-1); max2(m,j)$(ord(m)>1).. start(m,j) =g= completion(m-1,j); max3(m,j)..            start(m,j) =l= completion(m,j-1)+H*b(m,j); max4(m,j)..            start(m,j) =l= completion(m-1,j)+H*(1-b(m,j)); start.fx('machine1','job1')=0; cmaxdef.. cmax =e= completion('machine3','job10'); model flowshop /all/; option optcr=0; solve flowshop maximizing cmax using mip;

It is not 100% trivial as we need to make sure the max is attained exactly. We use the following set of inequalities:

In the case of minimization as formulated in http://yetanothermathprogrammingconsultant.blogspot.com/2012/04/longest-flow-shop.html we just could use bounds as the minimization would make sure that

start(m,j) = max(completion(m,j-1),completion(m-1,j))

would hold exactly (for the places where it matters – if there is slack in some machine this relation may not hold).

This model would look much simpler with a constraint programming formulation. There we could have used constructs like all-different and max to simplify things. Many modern modeling languages support constraint programming (unfortunately GAMS does not). Also note that AMPL using Cplex indicator variables allows this to be formulated quite elegantly (see http://yetanothermathprogrammingconsultant.blogspot.com/2009/07/formulation-ciminaibi.html). GAMS does not directly support indicator variables in the language, but instead some option files must be used.

The above MIP model solves in 0.2 seconds, clearly beating our brute force method.