Source-linked AI summary
Computing in Operations Research using Julia
Miles Lubin, Iain Dunning
TL;DR
Operations-research computing faces a divide between expressive but slower high-level languages and efficient but more cumbersome low-level languages. This paper explores Julia through optimization modeling and numerical algorithm implementations, finding that its performance results extend to realistic operations-research problems. It also develops publicly available optimization codes and evaluates Julia’s modeling and implementation capabilities.
Problem
The paper addresses the divide between expressive high-level languages and efficient low-level languages in operations-research numerical computing.
Method
The authors demonstrate Julia through linear and nonlinear algebraic modeling and a partial practical simplex implementation for numerical optimization.
Results
Julia performed within a factor of two of C on common basic tasks, and the paper tests whether this performance holds for realistic operations-research problems.
Takeaways & Limitations
The paper provides publicly available Julia codes for optimization-related tools and evaluates Julia’s suitability for realistic optimization computation.
Takeaways & Limitations
Retrofitting just-in-time compilation to existing languages has had mixed success because language design can conflict with compiler inference and compatibility with wider package ecosystems.
Abstract
from arXiv · showhide
The state of numerical computing is currently characterized by a divide between highly efficient yet typically cumbersome low-level languages such as C, C++, and Fortran and highly expressive yet typically slow high-level languages such as Python and MATLAB. This paper explores how Julia, a modern programming language for numerical computing which claims to bridge this divide by incorporating recent advances in language and compiler design (such as just-in-time compilation), can be used for implementing software and algorithms fundamental to the field of operations research, with a focus on mathematical optimization. In particular, we demonstrate algebraic modeling for linear and nonlinear optimization and a partial implementation of a practical simplex code. Extensive cross-language benchmarks suggest that Julia is capable of obtaining state-of-the-art performance.
1. Introduction
The paper examines Julia as a language designed to combine expressive programming with efficient numerical code, using operations-research optimization as motivating cases. It presents modeling tools and algorithm implementations, then evaluates whether Julia’s performance extends to realistic optimization problems.
- Motivation and approach: The paper uses linear and nonlinear programming to explore modern programming-language advances affecting operations-research computation.Linear programming is presented as a foundational operations-research problem whose algorithms and extensions occupy a large part of related computation.
- Motivation and approach: Julia is designed from the ground up to be expressive while enabling LLVM-based JIT compilation to generate efficient code.JIT compilation infers information at run time and uses those inferences to optimize generated machine code.
- Performance evaluation: Julia performed within a factor of two of C on common basic tasks in benchmarks reported by its authors.The paper tests whether these performance results hold for realistic problems of interest to operations research.
- Contributions: The authors develop publicly available Julia codes demonstrating features that facilitate implementing optimization-related tools.The paper includes JuMP for mixed-integer algebraic modeling, solver interfaces, nonlinear extensions, and a partial practical simplex implementation.
- Performance evaluation: The paper evaluates Julia’s suitability for low-level numerical optimization through a realistic partial implementation of the simplex algorithm.This complements the algebraic modeling demonstrations with an implementation-level performance assessment.
2. JuMP
JuMP is a Julia-based algebraic modeling language designed to combine expressive mathematical modeling with efficient model construction. It uses metaprogramming to transform expressions into sparse representations without operator overloading, achieving competitive benchmark performance.
- Interpreted-language AMLs can require substantial time to build sparse model representations, especially when models are rebuilt repeatedly during simulations.
- JuMP is a Julia package that combines high-level language benefits with the speed of commercial optimization products.
- Julia macros transform natural mathematical expressions into sparse internal model representations without operator overloading.This approach uses metaprogramming to bypass operator-overloading overhead during model construction.
- JuMP provides sufficient functionality for linear optimization while using a number of code lines similar to AMPL.
- JuMP macros reserve expression storage and fill coefficient arrays directly, avoiding many temporary objects associated with operator overloading.
- In p-median benchmarks, JuMP was within a factor of two of AMPL, comparable to Gurobi’s C++ interface, and an order of magnitude faster than Python-based modeling.
L LP MPS MPS LP MPS LP MPS LP
The benchmark tables report total processing time in seconds for linear-quadratic control models, while the software can produce LP and MPS output files.
- Table 2 reports linear-quadratic control benchmark results using grid size N=M and total processing time in seconds.
- The modeling system can produce output files in LP and MPS formats when available.
- The benchmark comparison includes JuMP/Julia, AMPL, Gurobi/C++, PuLP/PyPy, and Pyomo.
3. Nonlinear Modeling
Nonlinear AMLs can be impractically slow to build large-scale models in high-level languages because expression handling creates temporary objects or requires slow string parsing. Julia’s metaprogramming and JIT compilation generate nonlinear expressions and derivative evaluators with performance comparable to commercial AMLs, while the implementation remains a proof of concept.
- Large-scale nonlinear models can be significantly slower to build in open-source AMLs, sometimes becoming impractical.
- MATLAB and Python lack programmatic access to their optimized expression parsers, forcing operator overloading or manual string parsing for expressions such as y*sin(x).Operator overloading can create many temporary objects, while string parsing may be slow and disconnect expressions from the surrounding language.
- Julia’s first-class expression-parser access and metaprogramming features facilitate generated code with performance comparable to commercial AMLs.
- 3.1. Implementation in Julia: Julia macros transform abstract nonlinear expression trees into runtime code by splicing numeric constants and variable placeholders without expensive eval calls.The macro receives symbols before runtime values exist, then generates code that replicates the input expression with values inserted at runtime.
- 3.1. Implementation in Julia: Symbolic chain-rule differentiation produces sparse Jacobians, and JIT-compiled evaluation is reported as fast as AMPL through the low-level amplsolver library.Equivalent expression trees are differentiated only once, reducing repeated symbolic work.
- 3.2. Computational tests: Julia performs as fast as AMPL, if not faster, while Pyomo is over 50x slower on the largest instances and YALMIP cannot process the largest clnlbeam instance within an hour.The comparison includes the AMPL intermediate nl-file step and uses nonlinear benchmark instances.
4. Implementing Optimization Algorithms
The paper implements selected simplex operations in Julia to test whether high-level code can approach low-level performance while exploiting sparse linear algebra. The partial implementation emphasizes sparse matrix-vector products, ratio testing, and modified vector updates, then benchmarks Julia against other languages.
- Benchmark design: The benchmark compares Julia implementations with C++, MATLAB, and Python on core operations using realistic instances whose sparsity affects execution time.The C++ comparison includes a variant with explicit bounds checking.
- Simplex implementation: Sparse matrix-vector products exploit selected columns, sparse vectors, and matrix formats tailored to each operation.The implementation uses CSC or CSR storage and can skip dot products for unselected columns.
- Simplex implementation: The two-pass minimum ratio test is designed to avoid numerical instability from small α_i values, while permitting infeasibilities up to tolerance ϵ_D.The implementation handles both upper and lower bounds, although the displayed algorithm is simplified.
- Results: Vector-sparse routines significantly reduce execution time except on the smaller greenbea instance, although PyPy performs relatively poorly on them.The results qualitatively corroborate earlier reports that Julia is within a factor of 2 of equivalent low-level implementations.