Experimental Analysis of Algebraic Modelling Languages for Mathematical Optimization

1 Introduction

Many real-world problems are routinely solved using modern optimization tools (e.g. Abhishek et al., 2010; Fragniere and Gondzio, 2002; Groër et al., 2011; Paulavičius and Žilinskas, 2014; Paulavičius et al., 2020a, 2020b; Pistikopoulos et al., 2015). Internally, these tools use the combination of a mathematical model with an appropriate solution algorithm (e.g. Cosma et al., 2020; Fernández et al., 2020; Gómez et al., 2019; Lee et al., 2019; Paulavičius and Žilinskas, 2014; Paulavičius et al., 2014; Paulavičius and Adjiman, 2020; Stripinis et al., 2019, 2021) to solve the problem at hand. Thus, the way mathematical models are formulated is critical to the impact of optimization in real life.

Mathematical modelling is the process of translating real-world business problems into mathematical formulations whose theoretical and numerical analysis can provide insight, answers, and guidance beneficial for the originating application (Kallrath, 2004), including the current Covid-19 pandemic (Rothberg, 2020). Algebraic modelling languages (AMLs) are declarative optimization modelling languages, which bridge the gap between model formulation and the proper solution technique (Fragniere and Gondzio, 2002). They enable the formulation of a mathematical model as a human-readable set of equations while not requiring to specify how the described model should be solved or what specific solver should be used.

Models written in an AML are known for the high degree of similarity to the mathematical formulation. This aspect distinguishes AMLs from other types of modelling languages, like object-oriented (e.g. OptimJ), solver specific (e.g. LINGO), or general-purpose (e.g. TOMLAB) modelling languages. Such an algebraic design approach allows practitioners without specific programming or modelling knowledge to be efficient in describing the problems to be solved. It is also important to note that AML is then responsible for creating a problem instance that a solution algorithm can tackle (Kallrath, 2004). Since many AMLs are integral parts of a specific modelling system, it is essential to isolate a modelling language’s responsibilities from the overall system. In general, AMLs are sophisticated software packages that provide a crucial link between an optimization model’s mathematical concept and the complex algorithmic routines that compute optimal solutions. Typically, AML software automatically reads a model and data, generates an instance, and conveys it to a solver in the required form (Fourer, 2013).

From the late 1970s, many AMLs were created (e.g., GAMS, McCarl et al., 2016, AMPL, Fourer, 2003) and are still actively developed and used today. Lately, new open-source competitors to the traditional AMLs started to emerge (e.g., Pyomo, Hart et al., 2017, 2011, JuMP, Dunning et al., 2017; Lubin and Dunning, 2015). Therefore, we feel that a review and comparison of the traditional and emerging AMLs are needed to examine how the current landscape of AMLs looks.

The remainder of the paper is organized as follows. In Section 2, we review the essential characteristics of AMLs and motivate our selection of AMLs for the current review. In Section 3, we investigate how the requirements for a modern AML are met within each of the chosen languages. In Section 4, we examine the characteristics of AMLs using an extensive benchmark. In Section 5, we investigate the presolve impact on solving. Finally, we conclude the paper in Section 6.

2 Algebraic Modelling Languages

The first algebraic modelling languages, developed in the late 1970s, were game-changers. They allowed separating the model formulation from the implementation details (Kallrath, 2004) while keeping the notation close to the problem’s mathematical formulation (Fragniere and Gondzio, 2002). Since the data appears to be more volatile than the problem structure, most modelling language designers insist on the data and model structure being separated (Hürlimann, 1999). Therefore, the central idea in modern AMLs is the differentiation between abstract models and concrete problem instances (Hart et al., 2011). A specific model instance is generated from an abstract model using data. This way, the model and data together specify a particular instance of an optimization problem for which a solution can be sought. This is realized by replicating every entity of an abstract model over the different elements of the data set. Such a feature often is referred to as a set-indexing ability of the AML (Fragniere and Gondzio, 2002).

Essential characteristics of a modern AML could be defined in the following way (Kallrath, 2004):

1.

problems are represented in a declarative1 way;

2.

there is a clear separation between problem definition and the solution process;

3.

there is a clear separation between the problem structure and its data.

Besides, the support for mathematical expressions and operations needed for describing non-linear models is often considered an important feature of an AML (Kallrath, 2004). Moreover, it is worth observing that most interpreters included in today’s AMLs are based on automatic differentiation (Fragniere and Gondzio, 2002), a process in which the modeling language can compute the derivatives of problems from the model description without the assistance of the user (Kallrath, 2004). This motivates us to include automatic differentiation as an additional, important feature of a modern AML.

The algebraic expressions are useful in describing individual models and describing manipulations on models and transformations of data. Thus, almost as soon as AML became available, users started finding ways to adapt model notations to implement sophisticated solution strategies and iterative schemes. These efforts stimulated the evolution within AMLs of scripting features, including statements for looping, testing, and assignment (Fourer, 2013). Therefore, scripting capabilities are an integral part of AMLs.

For this review, we have chosen five AMLs: AIMMS, AMPL, GAMS, JuMP, and Pyomo. The selection was based on the following criteria:

•

AMLs which won 2012 INFORMS Impact Prize award2 dedicated to the originators of the five most important algebraic modelling languages: AIMMS, AMPL, GAMS, LINDO/LINGO, and MPL;

•

the popularity of AMLs based on NEOS Server model input statistics for the year 2020;3

•

open-source options that are attractive for the academic society or in situations where budgets are tight.

We have chosen to include GAMS and AMPL based on NEOS Server popularity, respectively, with 49% and 47% share of jobs executed via the NEOS platform in 2020. AIMMS was added as an example of AML with a graphical application development environment. JuMP and Pyomo were included as the most prominent open-source AMLs. We have decided to exclude MPL since it has not been updated for the last five years. We have also excluded LINGO as a solver specific modelling language.

3 Comparative Analysis of AMLs Characteristics

In the following section, we investigate how each of the chosen languages meets the requirements for a modern AML defined in the previous section. The websites of the AMLs and vendor documentation are used for this comparison. Any support of the identified features and capabilities are validated against the documentation the suppliers of the AMLs provide. Besides, an in-depth survey concluded by Robert Fourer in Linear Programming Software Survey (Fourer, 2017) is also used as a reference. Later on, a more practical comparison of AML characteristics is conducted to identify the potential ease of use of AML in daily work.

3.2 Practical Comparison of AMLs

For the first practical comparison of the selected AMLs, a classical Dantzig Transportation Problem was chosen (Dantzig, 1963). In this problem, we are given the factories’ supplies and the markets’ demands for a single commodity. We have also given the unit costs of shipping the product from factories to the markets. The goal is to find the least costly shipping schedule that meets the requirements at markets and supplies at factories.

The transportation problem formulated as a model in all five considered AML is compared based on the following criteria:

•

model size in bytes;

•

model size in the number of code lines;

•

model size in the number of language primitives used;

•

model instance creation time.

Since the transportation problem is a linear programming (LP) type of problem, we have chosen to measure the model instance creation time as the time needed to export a concrete model instance to MPS 4 format supported by most LP solvers. The following sources provided sample implementations of the transportation problem for the AMLs under consideration: •

AIMMS Wikipedia page;5

•

AMPL model in GNU Linear Programming Kit;6

•

GAMS Model Library;7

•

JuMP Examples;8

•

Pyomo Gallery.9

Transportation problem models in all five AMLs are given in Appendix B, Listings 3–6. It should be noted that the textual representation of an AIMMS model presents the model as a tree of attributed identifier nodes. It reflects how the model is given to the modeller in the AIMMS IDE and is typically generated by the AIMMS IDE. Moreover, it is worth noting that for the sake of simplicity, the problem model samples are concrete models, i.e. data of the model instance is described alongside the model structure.

OPEN IN VIEWER

Table 2

Comparison of transportation problem models.

Criteria	AIMMS	AMPL	GAMS	JuMP	Pyomo
Size in bytes	2229	683	652	632	1235
Lines of code	68	24	31	18	29
Primitives used	9	5	8	4	6

A comparison of the sample Transportation Problem model’s characteristics in all reviewed AMLs is given in Table 2. The simplification of the model implementations provided in the literature sources is made in the following way:

•

all optional comments, explanatory texts, and documentation are removed;

•

all empty lines are excluded;

•

parts of the code responsible for calling the solver and displaying results are omitted;

•

while counting AML primitives generic functions (sum, for), data loading directives (data), arithmetical and logical operators are excluded.

It can be seen from Table 2 that models implemented in AMPL, GAMS, and JuMP are the most compact ones, while the model written in AIMMS is much more verbose, and Pyomo lies in the middle. AIMMS propagates the creation of models using a graphical user interface (GUI) while keeping the model’s source code hidden from the modeller. Naturally, there is not much focus on how the model is stored. We can argue that while the GUI-based approach might be convenient to some modellers, it enforces greater vendor lock-in and makes the model’s extensibility and maintainability harder.

While comparing the number of language primitives required to create a model, JuMP and AMPL showed the best results, which allows us to predict that these modelling languages might have a more gentle learning curve. Therefore, we can conclude that in the context of the reviewed algebraic modelling languages, JuMP allows formulating an optimization problem most concisely.

The creation time of the transportation problem model instance defined in each AMLs was used to measure a model loading. The process was done in the following steps: 1.

loading model instance from a problem definition written in the native AML;

2.

exporting model instance to MPS format;

3.

measuring total execution time;

4.

investigating characteristics of an instance model.

Since the AIMMS system’s creators did not respond to an academic license request, we could not include AIMMS in the benchmark. Generated model instances in MPS format can be found in our GitHub repository’s models directory (Jusevičius and Paulavičius, 2019).

OPEN IN VIEWER

Table 3

Characteristics of the created transportation model instances.

Characteristic	AMPL	GAMS	JuMP	Pyomo
Constraints	6	6	6	6
Non zero elements	13	19	13	13
Variables	7	7	7	7

The characteristics of the created model instances can be seen in Table 3. We can conclude that all modelling languages have created a model instance using the same amount of variables and constraints. However, the definition of nonzero elements is different between GAMS and other modelling systems.

OPEN IN VIEWER

Table 4

Total time of consecutive transportation model instance creation runs.

No. of runs	AMPL	GAMS	JuMP	Pyomo
1 run	30 ms	170 ms	28341 ms	720 ms
10 runs	220 ms	1730 ms	32199 ms	7280 ms
100 runs	2130 ms	16490 ms	58151 ms	79600 ms

In Table 4, the benchmark results of model instance creation time are provided. We have tried to run multiple consecutive model instance creations (10 runs, 100 runs) to identify if the modelling system uses any caching. We can exhibit that AMPL showed significantly better results compared to others. This allows concluding that AMPL is the most optimized from a performance point of view. On the other hand, the poor JuMP performance confirms the Dunning et al. (2017) statement that JuMP has a noticeable start-up cost10 of a few seconds even for the smallest instances. In our case, only the initialization of the JuMP package took around 7 seconds. We also observed a significant speed-up in multiple consecutive model instances creation, which also confirms (Dunning et al., 2017) results. When a family of models is solved multiple times within a single session, this compilation cost is only paid for the first time that an instance is solved.

4 Performance Benchmark of AMLs

All examined AMLs support all types of traditional optimization problems; however, it is unclear how efficiently each AML can handle large model loading and what optimizations are applied during model instance creation. It would also be of great value to analyse how each of the modelling languages performs within an area of the specific type of optimization problems (e.g. linear, quadratic, nonlinear, mixed-integer). To give such a comparison and thoroughly examine characteristics of AMLs, a more extensive benchmark involving much larger optimization problem models is needed. Therefore, a large and extensive library of sample optimization problems for the analysed AMLs has to be used.

4.1 AMLs Testing Library

We have chosen the GAMS Model Library11 as a reference for creating such a sample optimization problem suite against which future research will be done. Automated shell script gamslib-convert.sh was created to build such a library. It can be found in the tools directory of our GitHub repository (Jusevičius and Paulavičius, 2019). A detailed explanation of how the test library creation tool works and the issues identified in the GAMS Library are provided in Appendix A. As a result of the transformation, we compiled a library consisting of 296 sample problems in AMPL, GAMS, JuMP, and Pyomo scalar model formats.

4.2 Model Instance Creation Time

The generated library was used to determine the amount of time each modelling system requires to create a model instance of a particular problem. We wrote load-benchmark.sh shell script available in the tools directory of our GitHub repository (Jusevičius and Paulavičius, 2019), which loads each model into the particular modelling system and then exports it to the format understandable by the solvers. We have chosen the .nl (Gay, 2005) format as the target format acceptable by the solvers, as .nl supports a wide range of optimization problem types. The benchmark measures the time the modelling system takes to perform both model instance creation and export operations.

We have chosen to exclude sample problems with conversion errors from the benchmark (more information about them in Appendix A). Only the models that were successfully processed by all modelling systems were compared. This reduced the scope of our benchmark to 268 models.

OPEN IN VIEWER

Fig. 1

Average model instance creation time.

OPEN IN VIEWER

Fig. 2

Average large model instance creation time

Benchmark methodology, hardware, and software specifications can be found in our GitHub repository (Jusevičius and Paulavičius, 2019). Detailed results are available in the model-loading-times.xlsx workbook in the benchmark section of our GitHub repository (Jusevičius and Paulavičius, 2019). We have provided a summary of the average model instance creation time split by the problem type in Fig. 1. We can see the trend exhibited in the transportation problem model benchmark persists. AMPL is still a definite top performer, while JuMP and Pyomo perform the worst. There are no significant variations between different optimization problem types except for JuMP, where the model instance creation time tends to vary significantly while working with different types of problems. Moreover, as confidence intervals show, the variation between different models of the same type is also more significant once using JuMP We tend to believe this is caused by Julia’s dynamic nature and the mix of run time compilation and caching of similar JuMP models.

We have observed that the average difference between AMPL and other contenders increases when the models become larger. Comparing instance creation times of large models (models having more than 500 equations, 8 such models in the testing library) reveals 11 times the difference between AMPL and GAMS, 38 times the difference between AMPL and Pyomo, and close to 100 times the difference between AMPL and JuMP. The difference between GAMS and Pyomo stayed roughly the same – around 3.5 times. The summary of the large model instance creation time can be seen in Fig. 2.

Thus, we can conclude that out of the reviewed AMLs, AMPL is a clear top-performing AML when it comes to the model instance creation time.

4.3 JuMP Benchmark

A similar time benchmark of the model instance creation has already been conducted (Dunning et al., 2017), where a smaller set of large models is used. While some of the trends exhibited in our benchmark persist (AMPL is the fastest, GAMS comes second), JuMP performance in our and Dunning et al. (2017) benchmarks differs significantly. This leads us to compare the benchmark methodology and results by conducting the benchmark described by Dunning et al. (2017).

First of all, our and their time benchmark methodologies differ. In comparison, we are trying to be solver independent and instruct AML to export the generated model instance to NL file, Dunning et al. (2017) attempt to solve the model using Gurobi solver and measure the time until Gurobi reports model characteristics. We believe that while our approach can be impacted by the system’s input/output performance, using a specific solver heavily depends on how the solver interface is implemented for a particular AML.

In the following, we have conducted two benchmarks – one as described in the original article and the second one using our method of exporting to a NL file. Results of the benchmarks can be seen in Tables 5, 6.

OPEN IN VIEWER

Table 5

JuMP benchmark using Dunning et al. (2017) method (in milliseconds).

Model	AMPL	GAMS	Pyomo	JuMP (DIRECT)	JuMP (CACHE)
lqcp-500	2093	2271	17000	17388	37317
lqcp-1000	8075	11995	139201	24590	44575
lqcp-1500	18222	38813	322604	39370	66566
lqcp-2000	32615	93586	575406	57597	88833
fac-25	407	480	7442	17517	39245
fac-50	2732	2884	43106	21331	47735
fac-75	9052	12422	150550	31582	57432
fac-100	20998	29144	393200	61326	93129

Before running the benchmarks, we had to rewrite some parts of the sample lqcp and facility JuMP models since syntax changes were introduced between JuMP v0.12 (used by Dunning et al. (2017)) and JuMP v0.21.5 (used by us, the latest version at the time of writing). Our benchmark was also conducted using newer versions of other AMLs–AMPL Version 20190207, GAMS v.32, Pyomo 5.7 (Python 3.8.3), Gurobi 9.0.

Additionally, we wanted to test JuMP’s new abstraction layer’s performance for working with solvers called MathOptInterface.jl (MOI). Therefore, we have tried both CachingOptimizer and DIRECT modes. As seen in Table 5, the DIRECT mode performed much better than the CachingOptimizer mode for both lqcp and facility models. An average difference in the instance creation time is close to two times, which leads us to suggest modellers evaluate MOI type choice based on specific use cases carefully.

OPEN IN VIEWER

Table 6

JuMP benchmark using export to NL method (in milliseconds).

Model	AMPL	GAMS	Pyomo	JuMP
lqcp-500	2716	3265	39988	20424
lqcp-1000	10503	14394	161404	80578
lqcp-1500	25402	49822	307121	483268
lqcp-2000	42780	125564	>10 min	>10 min
fac-25	409	502	9420	8163
fac-50	2837	2993	43087	31799
fac-75	10879	13457	143286	219548
fac-100	23474	32128	328170	>10 min

Overall, both benchmarks confirmed our observation that JuMP suffers from the long warm-up time required to pre-compile JuMP libraries. Results were also consistent with the patterns exhibited during the full gamslib benchmark performed earlier. We could not reproduce the JuMP performance metrics reported by Dunning et al. (2017), where JuMP always outperforms Pyomo. Using the original benchmark method, JuMP outperformed Pyomo only once the model size increased. However, while using our export to NL file method, JuMP, on the contrary, started to fall behind Pyomo once model size increased.

The reported differences between our and the original benchmark (Dunning et al., 2017) might be caused by multiple factors such as different JuMP versions used, improved Pyomo performance, or different Gurobi versions. It is important to stress that JuMP is a very actively developed AML that underwent significant changes during the last years. We think that it could be valuable to explore why the performance could have degraded and the reasons for such slow I/O operations performance revealed during writing to a NL file benchmark.