Wednesday, January 14, 2015

Pyomo and GAMS

Nice Pyomo introduction: https://www.youtube.com/watch?v=cjMkVHjhSBI

A few small mistakes wrt GAMS:

image

GAMS does not have a complements keyword. Instead it is using a simple notation to match variables to equations in the model statement:

MODEL TRNSP / PROFIT.X, SUPPLY.W, DEMAND.P/ ;

Secondly it is mentioned that GAMS has a presolve facility. Amazingly it does not. It leaves it to the solver to do that. For LP and MIP problems this may be not so bad, although a built-in presolver has some advantages (smaller problems and solutions files to exchange with the solver and most likely better error messages – we have more context available – the importance of this point is often underestimated IMHO). For NLP/MINLP models a built-in presolver can be more important: many nonlinear solvers do not have a presolver built-in. In general a nonlinear presolver is difficult to do inside the solver as the solver doesn’t know the algebraic form – the math – of the nonlinear functions they are executing (of course they can work on the linear part of the model). The most popular NLP solver with GAMS is CONOPT which actually has very good nonlinear presolving facilities (tailored for use with GAMS) and provides excellent diagnostics. But models to be solved with say MINOS or SNOPT can benefit hugely from a presolver (I have seen models being solved much more easily by AMPL with these solvers thanks to the presolving capabilities of AMPL; it really can make a huge difference). Having said this, careful modeling can lead to models where there is little room for the presolver to make a large difference (I try to implement my models this way – I always feel guilty of not doing a good job when the presolver if taking out large chunks of the model. Actually cleaning up a model to get rid of parts that will be presolved away is often a good exercise.).

Sometimes small presolve-like operations need to be performed by GAMS. Often this has to do with the fact that GAMS has an objective variable while solvers have an objective function. For an LP/MIP we would not care (just add an objective row with a single nonzero element at the location of the objective variable). But in case of QPs it is important to substitute out the objective variable (for NLPs this can also be important: linearly constrained problems may be much easier to solve). Even those simple transformation steps are only done in a somewhat half-baked fashion. Here is an LP file for a small QP problem fed into Cplex:

\ENCODING=ISO-8859-1
\Problem name: gamsmodel

Maximize
_obj: zobj + [ 2 z(bin1) ^2 + 2 z(bin2) ^2 + 2 z(bin3) ^2 + 2 z(bin4) ^2
       + 2 z(bin5) ^2 ] / 2
Subject To
objective:         zobj  = 0
bincapacity(bin1): 12 x(job1.bin1) + 11 x(job2.bin1) + 9 x(job3.bin1)
                    + 10 x(job4.bin1) + z(bin1)  = 10
bincapacity(bin2): 12 x(job1.bin2) + 11 x(job2.bin2) + 9 x(job3.bin2)
                    + 10 x(job4.bin2) + z(bin2)  = 12
bincapacity(bin3): 12 x(job1.bin3) + 11 x(job2.bin3) + 9 x(job3.bin3)
                    + 10 x(job4.bin3) + z(bin3)  = 10
bincapacity(bin4): 12 x(job1.bin4) + 11 x(job2.bin4) + 9 x(job3.bin4)
                    + 10 x(job4.bin4) + z(bin4)  = 14
bincapacity(bin5): 12 x(job1.bin5) + 11 x(job2.bin5) + 9 x(job3.bin5)
                    + 10 x(job4.bin5) + z(bin5)  = 16
assign(job1):      x(job1.bin1) + x(job1.bin2) + x(job1.bin3) + x(job1.bin4)
                    + x(job1.bin5)  = 1
assign(job2):      x(job2.bin1) + x(job2.bin2) + x(job2.bin3) + x(job2.bin4)
                    + x(job2.bin5)  = 1
assign(job3):      x(job3.bin1) + x(job3.bin2) + x(job3.bin3) + x(job3.bin4)
                    + x(job3.bin5)  = 1
assign(job4):      x(job4.bin1) + x(job4.bin2) + x(job4.bin3) + x(job4.bin4)
                    + x(job4.bin5)  = 1
Bounds
      zobj Free
0 <= x(job1.bin1) <= 1
0 <= x(job1.bin2) <= 1
0 <= x(job1.bin3) <= 1
0 <= x(job1.bin4) <= 1
0 <= x(job1.bin5) <= 1
0 <= x(job2.bin1) <= 1
0 <= x(job2.bin2) <= 1
0 <= x(job2.bin3) <= 1
0 <= x(job2.bin4) <= 1
0 <= x(job2.bin5) <= 1
0 <= x(job3.bin1) <= 1
0 <= x(job3.bin2) <= 1
0 <= x(job3.bin3) <= 1
0 <= x(job3.bin4) <= 1
0 <= x(job3.bin5) <= 1
0 <= x(job4.bin1) <= 1
0 <= x(job4.bin2) <= 1
0 <= x(job4.bin3) <= 1
0 <= x(job4.bin4) <= 1
0 <= x(job4.bin5) <= 1
Binaries
x(job1.bin1)  x(job1.bin2)  x(job1.bin3)  x(job1.bin4)  x(job1.bin5)
x(job2.bin1)  x(job2.bin2)  x(job2.bin3)  x(job2.bin4)  x(job2.bin5)
x(job3.bin1)  x(job3.bin2)  x(job3.bin3)  x(job3.bin4)  x(job3.bin5)
x(job4.bin1)  x(job4.bin2)  x(job4.bin3)  x(job4.bin4)  x(job4.bin5)
End

Obviously there is some junk left around from what one could say is an incomplete attempt to remove the objective variable. Of course in practice the Cplex presolver will take these parts out.