Optimization and modeling in Python

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:

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

Google OR-Tools:

  • 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:

Github

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:

Categories: Integer programming, Operations Research [EN], Solvers [EN]
Tags: , ,

You May Also Like

Warehouse receiving process simulation

04/05/2024

Scheduling a CNC job shop machine park

03/19/2024

Simulation for yogurt supply chain optimization

03/15/2024

Leave a Reply

9 comments

Maziyar Khadivi says:

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

Keivan Tafakkori M.Sc. says:

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

Juan Laiglecia says:

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

Keivan Tafakkori M.Sc. says:

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.

Leandro says:

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!

Leave a Reply

Your email address will not be published. Required fields are marked *

This site uses Akismet to reduce spam. Learn how your comment data is processed.