Solving linear programming problems (LPPs) using Pyomo

In this tutorial, we learn how to solve linear programming problems (LPPs) using Pyomo. Pyomo stands for Python Optimization Modeling Objects. This is an algebraic modeling language.

Pyomo is a Python-based, open-source optimization modeling language with a diverse set of optimization capabilities. — http://www.pyomo.org/

In this tutorial, we consider an optimization problem and then we formulate it as a linear programming problem (LPP). After that, we learn how to write a Pyomo-based model using Python to solve the problem.

Example problem

Consider a chocolate manufacturing company that produces only two types of chocolate — A and B. Both the chocolates require Milk and Choco only. To manufacture each unit of A and B, the following quantities are required:

• Each unit of A requires 1 unit of Milk and 3 units of Choco,
• Each unit of B requires 1 unit of Milk and 2 units of Choco.

The company kitchen has a total of 5 units of Milk and 12 units of Choco. On each sale, the company makes a profit of

• Rs. 6 per unit A sold,
• Rs. 5 per unit B sold.

Now, the company wishes to maximize its profit. How many units of A and B should it produce, respectively?

Mathematical formulation

Now, we formulate the above optimization problem as an LPP. An LPP has mainly three components — decision variables, the objective function, and constraints.

Decision variables:

• x be the total number of units produced by A,
• y be the total number of units produced by B.

Objective function:

Maximize 6x + 5y

Constraints:

x + y ≤ 5

3x + 2y ≤ 12

x ≥ 0, y ≥ 0

Join Medium with my referral link - Dr. Soumen Atta, Ph.D.

Installing Pyomo

We can use the following command to install Pyomo:

pip install pyomo

In this tutorial, we assume that Python3 and IBM ILOG CPLEX Optimization Studio or later are installed on your machine.

Pyomo model

Now, we write down the Pyomo model for the example problem mentioned above.

Import required packages

from pyomo.environ import * 
from pyomo.opt import SolverFactory

Create a model

model = ConcreteModel(name=”(My-Pyomo-Model)”)

Here, we have created a Concrete Pyomo model. Note that the argument name is optional.

Define decision variables of the model

model.x = Var(domain = NonNegativeReals, name = ‘x’) 
model.y = Var(domain = NonNegativeReals, name = ‘y’)

Define the objective function of the model

model.objective = Objective(expr = 6*model.x + 5*model.y, sense = maximize)

Add constraints to the model

model.constraints = ConstraintList()
model.constraints.add(model.x + model.y <= 5); # constraint-1 
model.constraints.add(3*model.x + 2*model.y <= 12); # constraint-2 
# to suppress unnecessary outputs, use a semicolon at the end

Note that we don’t need to include the non-negativity constraints (x≥0, y≥0) because, by default, these are already included in the definition of the decision variables.

Create a solver

solver = SolverFactory(“/opt/ibm/ILOG/CPLEX_Studio1271/cplex/bin/x86–64_linux/cplexamp”)

Here, we have used the CPLEX solver. We can use other solvers too. Note that we may need to specify the absolute solver path if it is not available as a system PATH variable.

Solve the model

sol = solver.solve(model)

Print the solution

print(sol)

The output of the above print statement is shown below:

Problem: 
- Lower bound: -inf
 Upper bound: inf
 Number of objectives: 1
 Number of constraints: 2
 Number of variables: 2
 Sense: unknown
Solver: 
- Status: ok
 Message: CPLEX 12.7.1.0\x3a optimal solution; objective 27; 2 dual simplex iterations (1 in phase I)
 Termination condition: optimal
 Id: 0
 Error rc: 0
 Time: 0.2711672782897949
Solution: 
- number of solutions: 0
 number of solutions displayed: 0

As we can see above, the problem has been solved to optimality and the objective function value is 27.

Retrieving decision variables’ values

print(‘x =’, model.x.value)
print(‘y =’, model.y.value)

The outputs of the above print statements are shown below:

x = 2.0000000000000004
y = 2.9999999999999996

Retrieve the objective value

print(‘The objective value is =’,model.objective())

The output of the above print statement is:

The objective value is = 27.0

Save the model in a file

model.write();

The above command will create a file named ‘(My-Pyomo-Model)’.lp. The contents of this file are shown below:

\* Source Pyomo model name=‘(My-Pyomo-Model)’ *\
max 
x3:
+6 x1
+5 x2
s.t.
c_u_x4_:
+1 x1
+1 x2
<= 5
c_u_x5_:
+3 x1
+2 x2
<= 12
c_e_ONE_VAR_CONSTANT: 
ONE_VAR_CONSTANT = 1.0
bounds
 0 <= x1 <= +inf
 0 <= x2 <= +inf
end

Display model information

model.display()

This will display the following information:

Model ‘(My-Pyomo-Model)’
Variables:
 x : Size=1, Index=None
 Key : Lower : Value : Upper : Fixed : Stale : Domain
 None : 0 : 2.0000000000000004 : None : False : False : NonNegativeReals
 y : Size=1, Index=None
 Key : Lower : Value : Upper : Fixed : Stale : Domain
 None : 0 : 2.9999999999999996 : None : False : False : NonNegativeReals
Objectives:
 objective : Size=1, Index=None, Active=True
 Key : Active : Value
 None : True : 27.0
Constraints:
 constraints : Size=2
 Key : Lower : Body : Upper
 1 : None : 5.0 : 5.0
 2 : None : 12.0 : 12.0

This is the end of this tutorial. I have already written similar tutorials for solving linear programming problems (LPPs) using different algebraic modeling languages. You can find them from the links mentioned below:

Python API of DOcplex for solving linear programming problems

Solving Linear Programming Problems (LPPs) Using PuLP and Python

If you have any queries regarding this tutorial then please feel free to contact me. You can find me at https://www.soumenatta.com/.

Dr. Soumen Atta, Ph.D.