Operations Research (OR) involves experiments with optimization models. The aim is to find the best design, plan, or decision for a system or a human. Accordingly, these models consist of objectives and constraints. However, most of the available packages or software for OR are not free or open-source. Thus, the pace of knowledge transfer or even reproducibility of results generated by the optimization models has been reduced. Besides, the available optimization packages or software do not provide the same level of flexibility as we know it from e.g. a programming language such as Python*. *Additionally, most available software or packages are not easy to use.

Overall, building on previous results in OR has always been difficult and time-consuming, compared with other fields in computer science such as machine learning (ML) and simulation (SM). Accordingly, in this blog, I am motivated to introduce and review popular high-level interfaces for modeling optimization problems in Python. As a result, SCM and OR analysts can leverage better and more flexible tools.

Notably, I focus on packages that can model and solve ** mixed-integer linear programs**. A high-level programming language such as Python is used to formulate a problem. The problem is then solved by optimizers (solvers) written in low-level languages.

## Difference between modeling interfaces and solvers in Python

At first, the difference between modeling interfaces (modeling languages) and solvers for optimization problems should be emphasized and explained:

- An
:*interface*- Understands a specific language for commands (syntax).
- Depends on the programming language.
- Returns an intermediate file that a solver reads as its input (e.g., in .mps, .lp, .nl, .osil, .gms formats).
- Can select one from various solvers when solving a problem.
- Is freely available in any programming language.
- Provides model-building tools and basic data structures.
- Sets solver parameters.

- A
:*solver*- Can be commercial (is better) or open-source (might be inferior when solving large-scale optimization problems).
- May have a specific interface for itself.
- Provides libraries that programming languages understand.
- Does factorization, file parsing, parameter parsing, pre-solving, sparse matrix, array storage, memory management, and timing.
- Returns information describing the solution and status of the optimization problem.

## A test problem

Let us consider the following optimization model:

The solution space (feasible region) (as x is a discrete and y is a continuous decision variable) can be graphically derived as follows:

According to the gradient of the objective (the red arrow) and direction of the optimization (maximization), the green point is the optimal solution, in which x=1, y=1 and the optimal value of the objective is 7. I solve this simple mixed-integer linear programming model with eleven interface packages in Python and analyze their syntax and similarities.

## 1. Optimization with MIP in Python

** Python-MIP**:

- Modeling language for linear programming and mixed-integer linear programming in Python.
- Supported solvers: CLP, CBC, and Gurobi.

The following commented code aims at solving the proposed mixed-integer linear programming model with “mip” (the name of the package) in Python:

# Installation (Uncomment the Line Below) #!pip install mip # Import package import mip as op # Define environment prob = op.Model("MyOptProblem") # Define decision variables x = [prob.add_var(var_type=op.INTEGER) for i in range(1)] y = [prob.add_var() for i in range(1)] # Add objective function to the environment prob.objective = op.maximize(2*x[0]+5*y[0]) # Add constraints to the environment prob += 5*x[0]+3*y[0]<=10 prob += 2*x[0]+7*y[0]<=9 # The status of the solution print(prob.optimize()) # To display optimal decision variables print('x: ', x[0].x) print('y: ', y[0].x) # To display optimal value of objective function print("Optimal Value of Objective Is = ",prob.objective_value)

## 2. Optimization with PuLP in Python

** PuLP**:

- Modeling language for linear programming and mixed-integer linear programming in Python.
- Supported solvers: GLPK, COIN-OR CLP/CBC, CPLEX, GUROBI, MOSEK, XPRESS, CHOCO, MIPCL, and SCIP.

The following commented code aims at solving the proposed mixed-integer linear programming model with “pulp” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install PuLP # Import package import PuLP as op # Define environment & direction of optimization prob = op.LpProblem("MyOptProblem", op.LpMaximize) # Define decision variables x = op.LpVariable("x", lowBound = 0, upBound = None, cat='Continuous') y = op.LpVariable("y", lowBound = 0, upBound = None, cat='Integer') # Add objective function to the environment prob += 2*x+5*y, "Objective" # Add constraints to the environment prob += 5*x+3*y<=10, "Constraint1" prob += 2*x+7*y<=9, "Constraint2" # Solve the problem (other solvers: prob.solve(SOLVERNAME())) prob.solve() # The status of the solution print("Status:", op.LpStatus[prob.status]) # To display optimal decision variables for variables in prob.variables(): print(variables.name, "=", variables.varValue) # To display optimal value of objective function print("Optimal Value of Objective Is = ", op.value(prob.objective))

## 3. Optimization with Pyomo in Python

** Pyomo**:

- Modeling language for linear programming, quadratic programming, nonlinear programming, mixed-integer linear programming, mixed-integer quadratic programming, mixed-integer nonlinear programming, stochastic programming, generalized disjunctive programming, differential-algebraic equations, bilevel programming, and mathematical programs with equilibrium constraints in Python.
- Supported solvers: Visit this link.

The following commented code aims at solving the proposed mixed-integer linear programming model with “pyomo” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install pyomo # Import package import pyomo.environ as op # Define environment prob = op.ConcreteModel("MyOptProblem") # Define decision variables prob.x = op.Var([1],domain=op.NonNegativeReals) prob.y = op.Var([1],domain=op.PositiveIntegers) # Add objective function to the environment prob.OBJ = op.Objective(expr=2*prob.x[1]+5*prob.y[1]) # Add constraints to the environment prob.Constraint1 = op.Constraint(expr=5*prob.x[1]+3*prob.y[1]<=10) prob.Constraint2 = op.Constraint(expr=2*prob.x[1]+7*prob.y[1]<=9) # Solve the problem solver = op.SolverFactory('SOLVERNAME') results = solver.solve(prob) # To display results print(results)

## 4. Optimization with Google OR-Tools in Python

- Modeling language for constraint programming, linear programming, and mixed-Integer linear programming.
- Supported solvers: Gurobi, CPLEX, SCIP, GLPK, GLOP, and CP-SAT.

The following commented code aims at solving the proposed mixed-integer linear programming model with “ortools” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install ortools # Import package from ortools.linear_solver import pywraplp # Define environment def prob(): #(Other solvers: pywraplp.Solver.CreateSolver('SOLVERNAME')) op = pywraplp.Solver.CreateSolver('SCIP') # Define decision variables x = op.NumVar(0.0, op.infinity(), 'x') y = op.IntVar(0.0, op.infinity(), 'y') # Add objective function to the environment op.Maximize(2*x+5*y) # Add constraints to the environment op.Add(5*x+3*y<=10) op.Add(2*x+7*y<=9) # Solve the problem status = op.Solve() # The status of the solution if status == pywraplp.Solver.OPTIMAL: print('Solution:') print('Objective value =', op.Objective().Value()) print('x =', x.solution_value()) print('y =', y.solution_value()) else: print('The problem does not have an optimal solution.') if __name__ == '__main__': prob()

## 5. Optimization with Pymoo in Python

* pymoo*:

- Modeling language for unconstrained and constrained single- and multi-objective optimization problems.
- Supported solvers: Meta-heuristic algorithms that are introduced in this link.

The following commented code aims at solving the proposed mixed-integer linear programming model with “pymoo” (the name of the package) in Python:

# Installation (uncomment the line below) #!pip install pymoo # Import packages import numpy as np from pymoo.core.problem import Problem from pymoo.algorithms.moo.nsga2 import NSGA2 from pymoo.factory import get_sampling, get_crossover, get_mutation from pymoo.factory import get_termination from pymoo.optimize import minimize termination = get_termination("n_gen", 40) # Solver of the problem algorithm = NSGA2( pop_size=400, n_offsprings=10, sampling=get_sampling("real_random"), crossover=get_crossover("real_sbx", prob=0.9, eta=15), mutation=get_mutation("real_pm", eta=20), eliminate_duplicates=True ) # Define environment class MyOptProblem(Problem): def __init__(self): # Define decision variables largenumber = 10^10 super().__init__(n_var=2, n_obj=1, n_constr=2, xl=np.array([0,0]), xu=np.array([largenumber,largenumber])) def _evaluate(self, x, out, *args, **kwargs): # Add objective function to the environment f1 = -2*x[:,0] + -5*np.round(x[:,1]) out["F"] = np.column_stack([f1]) # Add constraints to the environment g1 = 5*x[:,0] + 3*np.round(x[:,1]) - 10 g2 = 2*x[:,0] + 7*np.round(x[:,1]) - 9 out["G"] = np.column_stack([g1, g2]) prob = MyOptProblem() # Solve the problem res = minimize(prob, algorithm, termination, seed=2, save_history=True, verbose=False) # Display results print("Optimal Value of Objective Is = ", -res.F[0])

## 6. Optimization with GEKKO in Python

** GEKKO**:

- Modeling language for large-scale linear, quadratic, nonlinear, and mixed-integer programming.
- Supported solvers: APOPT, BPOPT, IPOPT, SNOPT, and MINOS.

The following commented code aims at solving the proposed mixed-integer linear programming model with “gekko” (the name of the package) in Python:

# Installation (uncomment the line below) #!pip install gekko # Import packages from gekko import GEKKO # Define environment prob = GEKKO(remote=False) # Define decision variables x = prob.Var(lb=0,ub=None,integer=True) y = prob.Var(lb=0,ub=None) # Add objective function to the environment prob.Obj(-(2*x+5*y)) # Objective # Add constraints to the environment prob.Equation(5*x+3*y<=10) prob.Equation(2*x+7*y<=9) # Solve the problem (1: MINLP solver, 2,3: Other Solvers) prob.options.SOLVER=1 prob.solve(disp=False) # Display results print('Results') print('x: ' + str(x.value)) print('y: ' + str(y.value)) print('Optimal Value of Objective Is = ' + str(-prob.options.objfcnval))

## 7. Optimization with PICOS in Python

*PICOS*:

- Modeling language for several conic and integer programming problems.
- Supported solvers: Visit this link.

The following commented code aims at solving the proposed mixed-integer linear programming model with “picos” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install picos # Import package import picos as op # Define environment prob = op.Problem("MyOptProblem") # Define decision variables x = op.RealVariable("x", lower = 0) y = op.IntegerVariable("y", lower = 0) # Add objective function to the environment prob.set_objective('max', 2*x+5*y) # Add constraints to the environment prob += 5*x+3*y<=10 prob += 2*x+7*y<=9 # Solve the problem prob.solve(solver='glpk') # To display results print('x: ', x.value) print('y: ', y.value) print("Optimal Value of Objective Is = ", prob.obj_value())

## 8. Optimization with CVXPY in Python

*CVXPY*

- Modeling language for convex programming, geometric programming, mixed-integer convex programming, and quasiconvex programming.
- Supported solvers: GLPK_MI, CBC, SCIP, CPLEX, GUROBI, XPRESS, MOSEK, OSQP, SCS, and ECOS.

The following commented code aims at solving the proposed mixed-integer linear programming model with “cvxpy” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install cvxopt # !pip install cvxpy # Import package import cvxpy as op # Define decision variables x = op.Variable(1,integer=True) y = op.Variable(1) bound_x = [0 <= x] bound_y = [0 <= y] # Add objective function to the environment objective = op.Maximize(2*x+5*y) # Add constraints to the environment cons1 = [5*x+3*y<=10] cons2 = [2*x+7*y<=9] # Define environment prob = op.Problem(objective, cons1+cons2+bound_x+bound_y) # Solve the problem prob.solve(solver='GLPK_MI') # To display results print('x: ', x.value) print('y: ', y.value) print("Optimal Value of Objective Is = ", objective.value)

## 9. Optimization with Drake in Python

*DRAKE*

- Modeling language for linear programming, quadratic programming, nonconvex nonlinear programming, convex optimization, mixed-integer convex optimization, and other non-convex mathematical programs.
- Supported solvers: CLP, CSDP, DREAL, GUROBI, IPOPT, MOSEK, OSQP, and SNOPT.

(Note: It currently lacks full support of integer variables, better open-source solvers, and Anaconda-related coding environments).

The following commented code aims at solving the proposed mixed-integer linear programming model with “pydrake” (the name of the package) in Python:

# Installation (Uncomment the Line Below) #!pip install drake # Import package import pydrake as op from pydrake.solvers.gurobi import GurobiSolver # Define environment prob = op.solvers.mathematicalprogram.MathematicalProgram() # Define decision variables x = prob.NewBinaryVariables(4) y = prob.NewContinuousVariables(1) # Add objective function to the environment prob.AddCost(2*(x[0]+2*x[1]+4*x[2]+6*x[3])+5*y[0]) # Add constraints to the environment prob.AddConstraint(5*(x[0]+2*x[1]+4*x[2]+6*x[3])+3*y[0]<=10) prob.AddConstraint(2*(x[0]+2*x[1]+4*x[2]+6*x[3])+7*y[0]<=9) # The status of the solution result = GurobiSolver().Solve(prob) print(result.is_success()) # To display optimal decision variables print('x: ', result.GetSolution(x)) print('y: ', result.GetSolution(y)) # To display optimal value of objective function print("Optimal Value of Objective Is = ", result.get_optimal_cost())

## 10. Optimization with CyLP in Python

*CyLP*

- Modeling language for linear and mixed-integer programs.
- Supported solvers: CLP, CBC, and CGL.

The following commented code aims at solving the proposed mixed-integer linear programming model with “cylp” (the name of the package) in Python:

# Installation (Uncomment the Line Below) #!pip install cylp # Import package from cylp.cy import CyClpSimplex import cylp as op # Define environment prob = op.py.modeling.CyLPModel() # Add variables x = prob.addVariable('x', 1, isInt=True) y = prob.addVariable('y', 1) # Add constraints prob += 5*x[0]+3*y[0]<=10 prob += 2*x[0]+7*y[0]<=9 prob += y[0] >= 0 prob += x[0] >= 0 # Add constraints to the environment prob.objective = -1*(2*x[0]+5*y[0]) # The status of the solution cbcModel = op.cy.CyClpSimplex(prob).getCbcModel() print(cbcModel.solve()) print (cbcModel.status) # To display optimal decision variables print("x: ", cbcModel.primalVariableSolution['x'][0]) print("y: ", cbcModel.primalVariableSolution['x'][0]) # To display optimal value of objective function print("Optimal Value of Objective Is = ", -cbcModel.objectiveValue)

## 11. Optimization with PyMathProg in Python

*PyMathProg*

- Modeling language for linear and mixed-integer programs
- Supported solvers: GLPK

The following commented code aims at solving the proposed mixed-integer linear programming model with “pymprog” (the name of the package) in Python:

# Installation (uncomment the line below) # !pip install pymprog # Import package import pymprog as op # Define environment op.begin('MyOptProb') # Define decision variables x = op.var('x',bounds=(0,None),kind=int) y = op.var('y',bounds=(0,None)) # Add objective function to the environment op.maximize(2 * x + 5 * y, 'objective') # Add constraints to the environment 5*x + 3*y <= 10 2*x + 7*y <= 9 # Solve the problem op.solve() # To display optimal decision variables print("x: ", x.primal) print("y: ", y.primal) # To display optimal value of objective function print("Optimal Value of Objective Is = ", op.vobj()) op.end()

## 12. Optimization with Optlang in Python

*Optlang*

- Modeling language for linear, mixed-integer, and quadratic programs.
- Supported solvers: GLPK, CPLEX, GUROBI, and INSPYRED.

The following commented code aims at solving the proposed mixed-integer linear programming model with “optlang” (the name of the package) in Python:

# Installation (uncomment the line below) #!pip install optlang from __future__ import print_function import optlang as op # Define environment prob = op.Model(name='Simple model') # Define decision variables x = op.Variable('x', lb=0, type="integer") y = op.Variable('y', lb=0) # Add constraints to the environment c1 = op.Constraint(5*x+3*y, ub=10) c2 = op.Constraint(2*x+7*y, ub=9) prob.add([c1, c2]) # Add objective function to the environment obj = op.Objective(2*x+5*y, direction='max') prob.objective = obj # Solve the problem and report status status = prob.optimize() print("status:", prob.status) # To display optimal decision variables for var_name, var in prob.variables.iteritems(): print(var_name, "=", var.primal) # To display optimal value of objective function print("Optimal Value of Objective Is = ", prob.objective.value)

Click the button below to download the complete pack of the above coding examples:

If this article is going to be used in research or other publishing methods, you can cite it as Tafakkori (2021) (in text) and refer to it as follows: Tafakkori, K. (2021). Modeling and solving optimization problems in Python. *Supply Chain Data Analytics*. url: https://www.supplychaindataanalytics.com/modeling-and-solving-optimization-problems-in-python/

## Related content

If you are interested in mathematical programming and optimization in Python you can consider the following courses and programs:

Industrial engineer focusing on leveraging optimization methods and artificial intelligence technologies using multiple programming languages to empower a business to achieve its goals!

## 9 comments

Very insightful post. It seems to be much easier to model the problem in Python compared to commercial software.

That is exactly one of the advantages of some of these packages.

Thanks for your comment, and I am glad you found the article insightful!

Very interesting post. I hope in the near future, more optimization examples will be presented. May be compared with JuLia. Excelent work!!

Thanks for your comment! We try to cover more about optimization and machine learning and their integration in Python and other programming languages such as Julia.

We have started to cover Julia implementations for OR. See under “Blog” and go to the section “Julia” for a updated list of articles related to Julia programming.

I tried pip install all the libraries with python 3.11 and only 3 of them worked easily without any error:

Gekko

OptLang

PyMathProg

If it is a Python-version thing then be aware that you can have several versions of Python on the same machine. You should also use virtual environments.

Another thing that causes problems for me is the combination of anaconda, pip and virtual environments. I stopped using anaconda because of this.

If you have error messages you can post them here for feedback and documentation. Thanks!