Model generation in Julia

In http://www.optimization-online.org/DB_HTML/2015/04/4891.html we see some interesting numbers on generation times of large models in different modeling environments. AMPL and Julia/JuMP are doing very well on these benchmarks. E.g.

As https://plus.google.com/+VictorZverovich/posts/RtHfqgwh5QX indicates: AMPL is faster than Gurobi/C++. This does not come as a surprise: AMPL is around quite a while and is known to be quite fast. GAMS is slow compared to AMPL especially on the lqcp and acpower instances. Note that in general solution times are much larger than the modeling system needs to generate the problem. Only very easy LPs can be relatively expensive to generate.

Why is GAMS so slow on these instances? In https://github.com/mlubin/JuMPSupplement we can find much information so we can reproduce things (very praiseworthy to make this available).

lqcp-2000

Erwin@EKlap2 ~ $ gams lqcp.gms lo=3  | ts -s 00:00:00 --- Job lqcp.gms Start 05/05/15 15:13:48 24.4.1 r50296 WEX-WEI x86 64bit/MS Windows 00:00:00 GAMS 24.4.1   Copyright (C) 1987-2014 GAMS Development. All rights reserved 00:00:00 Licensee: Erwin Kalvelagen                               G140808/0001CV-GEN 00:00:00           Amsterdam Optimization Modeling Group                     DC10455 00:00:00 --- Starting compilation 00:00:00 --- lqcp.gms(65) 3 Mb 00:00:00 --- Starting execution: elapsed 0:00:00.003 00:00:00 --- lqcp.gms(63) 260 Mb 00:00:00 --- Generating MIQCP model cont5 --- lqcp.gms(65) 2406 Mb) 284 Mb 00:01:03 ---   4,004,002 rows  4,006,002 columns  23,994,008 non-zeroes 00:01:03 ---   18,010 nl-code  4,001 nl-non-zeroes --- lqcp.gms(65) 2404 Mb) 2406 Mb 00:01:03 --- Executing GUROBI (Solvelink=5): elapsed 0:01:02.947 00:01:03 00:01:03 Gurobi           24.4.1 r50296 Released Dec 20, 2014 WEI x86 64bit/MS Windows 00:01:03 00:01:03 Gurobi full license. 00:01:03 Gurobi library version 6.0.0 00:01:06 Space for names approximately 174.82 Mb 00:01:06 Use option 'names no' to turn use of names off 00:01:07 Starting Gurobi... 00:01:07 Optimize a model with 4004001 rows, 4006001 columns and 23990006 nonzeros

We do this in 67 seconds: a little bit better than reported. Not much I can do to speed this up: the structure of the offending equation pde has leads and lags which is not something GAMS is always doing extremely fast. Usually my client models are large and sparse and do not have this structure.

I actually think the equation pde is not 100% correctly formulated. With these discretizations it is important to pay extra attention to the border. This is especially the case as GAMS will set all terms with leads and lags that exceed the boundary are implicitly set to zero. Always inspect the equation listing in GAMS when there are leads and lags involved. However these issues are probably not important for these measurements.

fac-100

Erwin@EKlap2 ~ $ gams facility.gms lo=3 | ts -s 00:00:00 --- Job facility.gms Start 05/05/15 16:01:12 24.4.1 r50296 WEX-WEI x86 64bit/MS Windows 00:00:00 GAMS 24.4.1   Copyright (C) 1987-2014 GAMS Development. All rights reserved 00:00:00 Licensee: Erwin Kalvelagen                               G140808/0001CV-GEN 00:00:00           Amsterdam Optimization Modeling Group                     DC10455 00:00:00 --- Starting compilation 00:00:00 --- facility.gms(41) 3 Mb 00:00:00 --- Starting execution: elapsed 0:00:00.008 00:00:00 --- facility.gms(39) 4 Mb 00:00:01 --- Generating MIQCP model facility --- facility.gms(41) 1324 Mb) 235 Mb 00:00:11 ---   4,090,601 rows  4,080,601 columns  11,221,100 non-zeroes 00:00:11 ---   10,201,001 nl-code  3,060,300 nl-non-zeroes 00:00:12 ---   1,020,100 discrete-columns --- facility.gms(41) 1271 Mb) 1323 Mb 00:00:13 --- Executing GUROBI: elapsed 0:00:11.390 00:00:16 00:00:16 Gurobi           24.4.1 r50296 Released Dec 20, 2014 WEI x86 64bit/MS Windows 00:00:16 00:00:16 Gurobi full license. 00:00:16 Gurobi library version 6.0.0 00:00:20 Space for names approximately 165.81 Mb 00:00:20 Use option 'names no' to turn use of names off 00:01:06 Starting Gurobi... 00:01:06 Optimize a model with 3070501 rows, 4080601 columns and 8160800 nonzeros

Here we are actually slower than reported. The pure generation time is not that bad: 11 seconds. The main problem we see here is that the GAMS/Gurobi link is inefficient in handling the setup of the quadratic terms for these large models. This turned out to be a problem with Gurobi: adding individual quadratic constraints is expensive. A fix is on its way. Update: fixed in Gurobi 6.0.4.

acpower-100

Erwin@EKlap2 ~ $ gams opf_66200bus.gms lo=3 | ts -s 00:00:00 --- Job opf_66200bus.gms Start 05/06/15 10:29:26 24.4.1 r50296 WEX-WEI x86 64bit/MS Windows 00:00:00 GAMS 24.4.1   Copyright (C) 1987-2014 GAMS Development. All rights reserved 00:00:00 Licensee: Erwin Kalvelagen                               G140808/0001CV-GEN 00:00:00           Amsterdam Optimization Modeling Group                     DC10455 00:00:00 --- Starting compilation 00:00:00 --- opf_66200bus.gms(17) 2 Mb --- .IEEE66200.bus.gms(1059233) 15 Mb3) 15 Mb 00:00:01 --- opf_66200bus.gms(46) 15 Mb --- .IEEE66200.branch.gms(1322139) 59 Mb2) 41 Mb 00:00:01 --- opf_66200bus.gms(239) 59 Mb 00:00:02 --- Starting execution: elapsed 0:00:01.581 00:00:02 --- opf_66200bus.gms(2381609) 93 Mb 00:00:02 --- Generating NLP model acpower --- opf_66200bus.gms(2381610) 313 Mb8) 107 Mb 00:04:54 ---   132,401 rows  416,601 columns  1,975,401 non-zeroes 00:04:54 ---   14,576,003 nl-code  1,960,800 nl-non-zeroes 00:04:54 --- opf_66200bus.gms(2381610) 239 Mb 00:04:54 --- Executing IPOPT (Solvelink=5): elapsed 0:04:54.060 00:04:55 00:04:55 COIN-OR Ipopt    24.4.1 r50296 Released Dec 20, 2014 WEI x86 64bit/MS Windows 00:04:55 00:04:55 COIN-OR Interior Point Optimizer (Ipopt Library 3.11) 00:04:55 written by A. Waechter. 00:04:55 00:04:55 Note: This is the free version IPOPT, but you could also use the commercially supported and potentially higher performance version IPOPTH. 00:05:51 00:05:51 ****************************************************************************** 00:05:51 This program contains Ipopt, a library for large-scale nonlinear optimization. 00:05:51  Ipopt is released as open source code under the Eclipse Public License (EPL). 00:05:51          For more information visit http://projects.coin-or.org/Ipopt 00:05:51 ****************************************************************************** 00:05:51 00:05:51 This is Ipopt version 3.11, running with linear solver mumps. 00:05:51 00:05:51 Space for names approximately 14.45 MB. 00:05:51 Use statement '<modelname>.dictfile=0;' to turn dictionary off. 00:05:51 Number of nonzeros in equality constraint Jacobian...:  1068500 00:05:51 Number of nonzeros in inequality constraint Jacobian.:        0 00:05:51 Number of nonzeros in Lagrangian Hessian.............:   805300 00:05:51 00:05:54 Total number of variables............................:   161800 00:05:54                      variables with only lower bounds:        0 00:05:54                 variables with lower and upper bounds:    95700 00:05:54                      variables with only upper bounds:        0 00:05:54 Total number of equality constraints.................:   132400 00:05:54 Total number of inequality constraints...............:        0 00:05:54         inequality constraints with only lower bounds:        0 00:05:54    inequality constraints with lower and upper bounds:        0 00:05:54         inequality constraints with only upper bounds:        0 00:05:54

This is pretty bad for a somewhat small and somewhat dense model. We can make the model easily larger and hopefully sparser and “less” non-linear. After adding some extra options to reduce the size of the listing file we see:

Erwin@EKlap2 ~ $ gams acpower_ek.gms lo=3 limrow=0 limcol=0 | ts -s 00:00:00 --- Job acpower_ek.gms Start 05/06/15 10:57:38 24.4.1 r50296 WEX-WEI x86 64bit/MS Windows 00:00:00 GAMS 24.4.1   Copyright (C) 1987-2014 GAMS Development. All rights reserved 00:00:00 Licensee: Erwin Kalvelagen                               G140808/0001CV-GEN 00:00:00           Amsterdam Optimization Modeling Group                     DC10455 00:00:01 --- Starting compilation 00:00:01 --- acpower_ek.gms(17) 2 Mb --- .IEEE66200.bus.gms(1059233) 15 Mb2) 15 Mb 00:00:01 --- acpower_ek.gms(46) 15 Mb --- .IEEE66200.branch.gms(1322139) 59 Mb7) 43 Mb 00:00:02 --- acpower_ek.gms(263) 59 Mb 00:00:04 --- Starting execution: elapsed 0:00:01.566 --- acpower_ek.gms(2381632) 108 Mb4) 100 Mb 00:00:07 --- Generating NLP model acpower --- acpower_ek.gms(2381634) 445 Mb4) 118 Mb 00:00:57 ---   671,601 rows  955,801 columns  3,683,001 non-zeroes 00:00:57 ---   11,999,503 nl-code  3,063,000 nl-non-zeroes 00:00:57 --- acpower_ek.gms(2381634) 384 Mb 00:00:57 --- Executing IPOPT (Solvelink=5): elapsed 0:00:56.518 00:00:58 00:00:58 COIN-OR Ipopt    24.4.1 r50296 Released Dec 20, 2014 WEI x86 64bit/MS Windows 00:00:58 00:00:58 COIN-OR Interior Point Optimizer (Ipopt Library 3.11) 00:00:58 written by A. Waechter. 00:00:58 00:00:58 Note: This is the free version IPOPT, but you could also use the commercially supported and potentially higher performance version IPOPTH. 00:03:56 00:03:56 ****************************************************************************** 00:03:56 This program contains Ipopt, a library for large-scale nonlinear optimization. 00:03:56  Ipopt is released as open source code under the Eclipse Public License (EPL). 00:03:56          For more information visit http://projects.coin-or.org/Ipopt 00:03:56 ****************************************************************************** 00:03:56 00:03:56 This is Ipopt version 3.11, running with linear solver mumps. 00:03:56 00:03:56 Number of nonzeros in equality constraint Jacobian...:  2586800 00:03:56 Number of nonzeros in inequality constraint Jacobian.:        0 00:03:56 Number of nonzeros in Lagrangian Hessian.............:  3461700 00:03:56 00:04:08 Total number of variables............................:   701000 00:04:08                      variables with only lower bounds:        0 00:04:08                 variables with lower and upper bounds:    95700 00:04:08                      variables with only upper bounds:        0 00:04:08 Total number of equality constraints.................:   671600 00:04:08 Total number of inequality constraints...............:        0 00:04:08         inequality constraints with only lower bounds:        0 00:04:08    inequality constraints with lower and upper bounds:        0 00:04:08         inequality constraints with only upper bounds:        0

We see again that GAMS generates the model quickly, but that inside the solver link we have some inefficiencies. I am not sure what is happening in the 3 lost minutes, but it could be related in setting up the second derivatives. Further improvements can be obtained by telling GAMS to consider fixed variables as parameters.

Some general points:

• Generation time is often a minor part of the total turn around time
• Solvers are in most cases much slower to solve the problem than the time the modeling system needs to generate the model
• In many models we spend more time in data manipulation than in pure model generation. Many practical models have most of the code devoted to data manipulation (data import, data conversion, checks, aggregation, consolidation etc.)
• I believe some models will behave better with upcoming version of GAMS and Gurobi. These very large examples are excellent to find weak spots in the modeling software and will lead to improvements.
• I am really impressed with the JuMP performance.
• The models are very large. Lots of meaningful models are way smaller than these (but may be more complex). Systems like Pyomo should not be dismissed just on these numbers.

The paper is an interesting read. It touches upon many concepts that are not widely known.