Experimental Analysis of Algebraic Modelling Languages for Mathematical Optimization

3.1Comparative Analysis of the Features

We start by analysing how selected AMLs satisfy the three essential characteristics defined in the previous Section 2. In all reviewed AMLs, optimization problems are represented in a declarative way. Furthermore, since all of them are part of a specific modelling system, a clear separation between problem definition and the solution process in the context of the modelling system exists. The separation between the problem structure and its data is supported in all reviewed languages. It should be noted that GAMS, JuMP, and Pyomo also allow initiating data structures during their declaration, while AIMMS and AMPL only support it as a separate step in the model instance building process. However, while it might be convenient for building a simple model, we do not consider the lack of direct data structure initiation as an advantage since, in real-world cases, it is rarely needed. Therefore, we can conclude that all reviewed languages fulfill the essential characteristics of modern AMLs.

Next, in Table 1, we provide an overview of the key features each AML supports. For creating such a summary, we used the information provided by the AML vendors on their websites. All reviewed AMLs allow modelling problems in a solver independent manner. Additionally, AIMMS, JuMP, and Pyomo provide a more powerful way to define advanced algorithms using R, Julia, or Python programming languages. The ease of data input for the model differs among AMLs. While all of them support input from a flat file, some more advanced scenarios such as reading data from relational databases are more straightforward in AIMMS, JuMP, or Pyomo. AMPL and GAMS require a complicated setup instead (e.g. using ODBC drivers) to access the database. Wherein JuMP or Pyomo, a standard Julia or Python driver could be used to get data from relational and any other type of database supported by Python or Julia. Manipulation (e.g. transformation) of data is only supported by AIMMS, JuMP, and Pyomo.

When it comes to solver support, AMPL is the one supporting the most. However, it should be noticed that the categorization of solvers by supported problem types is different among vendors. Thus, in this comparison, we have reflected the information available from vendors harmonizing it across all of them. Solvers supported by JuMP and Pyomo require additional explanation. First, both AMLs support solvers compatible with AMPL (via AmplNLWriter package or ASL interface). Therefore, any solver that is equipped with an AMPL interface can be used by JuMP or Pyomo. This could allow us to state that JuMP and Pyomo support all AMPL solvers. However, we have excluded solvers supported via the AMPL interface. It might be needed for some commercial solvers to request a particular version from the solver’s vendor that comes with the AMPL interface. Second, since both AMLs are open-source, multiple third-party packages add support for specific solvers for each of AMLs. In Section 1, we counted only the solvers mentioned on the official JuMP and Pyomo websites.

Presolving capabilities are only available in AIMMS and AMPL. JuMP and Pyomo have programming interfaces for creating custom presolvers, however, none of them are provided out of the box. Only AIMMS provides a visualization of the solver results out of the box. Using Python or Julia libraries, it is possible to visualize the results produced by Pyomo and JuMP. However, it requires custom development, and none of the standard JuMP or Pyomo libraries are supporting that.

It is important to conclude that JuMP and Pyomo are open-source AMLs built on top of general-purpose programming languages, making them fundamentally different from the competitors. This allows researchers familiar with Julia or Python to learn, improve, and use JuMP or Pyomo much more comfortably. At the same time, it is practically impossible to introduce improvements to commercial counterparts.

3.2Practical 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:

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 MPS44 format supported by most LP solvers. The following sources provided sample implementations of the transportation problem for the AMLs under consideration:

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:

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:

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.

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 cost1010 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.1AMLs Testing Library

We have chosen the GAMS Model Library1111 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.2Model 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.

Fig. 1

Average model instance creation time.

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.3JuMP 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.

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.

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.

5.1Presolving in AMPL

To assess AMPL’s presolving performance, we gathered presolving characteristics while performing the model instance creation time benchmark. We have used 286 models that were successfully converted from GAMS original model to the AMPL scalar model.

A detailed report of the presolving applied to the specific model can be seen in the benchmark section of our GitHub repository (Jusevičius and Paulavičius, 2019), while the summary of it can be found in Table 7. We observed that AMPL presolver managed to simplify the models in 52.8% of the cases, out of which 5 times it could determine that the problem solution is not feasible, thus not even requiring to call the solver. On average, once applied, the AMPL presolver managed to reduce the model size by removing 18.42% of constraints and 10.73% of variables.

We can conclude that AMPL presolver is an efficient way to simplify larger problems, leading to improved solution finding performance once invoking a solver with an already reduced problem model instance. Moreover, determining not feasible models can help modellers debug the problem definition process and find errors in the model definition. This allows us to argue that presolving is an important capability of any modern AML.

5.2Presolve Impact on Solving

To evaluate if AMPL presolving has a positive impact on problem-solving, an additional benchmark was conducted. The benchmark included 146 out of 151 models to which AMPL has applied presolve in the model instance creation benchmark. Five models that AMPL presolve determined to be not feasible were excluded from the benchmark. Shell script solve-benchmark.sh provided in the tools directory of our GitHub repository (Jusevičius and Paulavičius, 2019) was created for executing such a benchmark. The script solves each model using one of the solvers and gathers output statistics to a report file.

We have chosen to solve the models using Gurobi and BARON solvers. Gurobi Optimizer (v.8.1.0) was chosen for solving LP, MIP, QCP, and MIQCP type of problems. At the same time, BARON (v.18.11.12) global solver was chosen for solving NLP, MINLP, MCP, MPEC, CNS, and DNLP problems. The solvers’ choice was motivated by the support for particular problem types.1313 and the popularity of solvers based on NEOS Server statistics.1414 Two attempts to solve each model were made. One with AMPL presolver turned on (default setting), and the second one with AMPL presolver turned off. After each run, solvers statistics, including iterations count, solve time (pure solve phase execution time), and objective, were gathered.

It is important to note that both BARON and Gurobi solvers have their presolve mechanisms (Puranik and Sahinidis, 2017; Achterberg et al., 2019). Thus the provided model is simplified by the solver too. This might result in very similar models being solved by the solver despite the AMPL presolve being turned on or off. However, the focus was on estimating AMPL presolve impact in real-life situations; therefore, a full benchmark was executed without changing the default solver behaviour. Later on, an additional benchmark was made to estimate the impact of AMPL presolve once solver presolve functionality is turned off.

Detailed AMPL presolve impact on solving report can be found in our GitHub repository’s (Jusevičius and Paulavičius, 2019) file ampl-solving-times.xlsx sheet Benchmark 1. Table 8 summarizes the positive and negative impact AMPL presolve had on solving the problems iteration and time-wise. Positive impact means fewer iterations or time was needed to solve a problem once the presolve was turned on. A negative impact means the opposite that more iterations or time was required. The impact is considered neutral if the number of iterations did not change or the required time was within the one-second tolerance level.

During this benchmark, 6 models failed to be solved due to solver limitations. Two models were deemed to be not feasible, and two were solved during the AMPL presolve phase. Solvers were capable of solving 41 models during the solver’s presolve phase. Moreover, for six models, the mismatching objective was found with AMPL presolve turned on and off. Overall, AMPL presolve positively impacted 26.43% of the cases iteration-wise and 47.86% time-wise. However, it hurt 20.71% of cases iteration-wise and 23.57% time-wise.

As mentioned earlier, both BARON and Gurobi solvers have their presolve mechanisms. An additional benchmark was made to test the impact of AMPL presolve with disabled solvers presolving. Since only Gurobi allows the user to disable presolve functionality, a subset of models previously solved with Gurobi was chosen. Detailed benchmark results can be seen in our GitHub repository’s (Jusevičius and Paulavičius, 2019) file ampl-solving-times.xlsx sheet Benchmark 2. The summary of the benchmark is provided in Tables 9, 10. Gurobi could not solve two MIP problems (clad and mws) in a reasonable time once Gurobi’s presolve functionality was turned off. Those models were excluded from the benchmark.

As seen once comparing these results in Tables 9 and 10, the AMPL presolve had a greater positive effect both iteration-wise (+25.5%) and time-wise (+10.2%) once Gurobi presolve was turned off. AMPL presolve also had a less neutral impact once the solver presolving was off, thus leading to the conclusion that during the first benchmark, some models were simplified to very similar ones before actually solving them.

As we can see from the benchmarks, presolving done by AML has inconclusive effects on the actual problem solving both iterations and time-wise. However, a positive impact is always more significant than the negative one, and it especially becomes evident once the solver does not have or use its problem presolving mechanisms. This allows us to conclude that the presolving capability of AML is an important feature of a modern algebraic modelling language. We can also advise choosing AML having presolving capabilities in cases the solver used to solve the problem does not have its presolving mechanism.