Introduction to Linear Programming with Python - Part 4

Real world examples - Blending Problem

We're going to make some sausages!

We have the following ingredients available to us:

Ingredient Cost (€/kg) Availability (kg)
Pork 4.32 30
Wheat 2.46 20
Starch 1.86 17

We'll make 2 types of sausage:

  • Economy (>40% Pork)
  • Premium (>60% Pork)

One sausage is 50 grams (0.05 kg)

According to government regulations, the most starch we can use in our sausages is 25%

We have a contract with a butcher, and have already purchased 23 kg pork, that will go bad if it's not used.

We have a demand for 350 economy sausages and 500 premium sausages.

We need to figure out how to most cost effectively blend our sausages.

Let's model our problem

pe = Pork in the economy sausages (kg)
we = Wheat in the economy sausages (kg)
se = Starch in the economy sausages (kg)
pp = Pork in the premium sausages (kg)
wp = Wheat in the premium sausages (kg)
sp = Starch in the premium sausages (kg)

We want to minimise costs such that:

Cost = 0.72(pe + pp) + 0.41(we + wp) + 0.31(se + sp)

With the following constraints:

pe + we + se = 350 * 0.05
pp + wp + sp = 500 * 0.05
pe ≥ 0.4(pe + we + se)
pp ≥ 0.6(pp + wp + sp)
se ≤ 0.25(pe + we + se)
sp ≤ 0.25(pp + wp + sp)
pe + pp ≤ 30
we + wp ≤ 20
se + sp ≤ 17
pe + pp ≥ 23


In [1]:
import pulp

In [2]:
# Instantiate our problem class
model = pulp.LpProblem("Cost minimising blending problem", pulp.LpMinimize)

Here we have 6 decision variables, we could name them individually but this wouldn't scale up if we had hundreds/thousands of variables (you don't want to be entering all of these by hand multiple times).

We'll create a couple of lists from which we can create tuple indices.


In [4]:
# Construct our decision variable lists
sausage_types = ['economy', 'premium']
ingredients = ['pork', 'wheat', 'starch']

Each of these decision variables will have similar characteristics (lower bound of 0, continuous variables). Therefore we can use PuLP's LpVariable object's dict functionality, we can provide our tuple indices.

These tuples will be keys for the ing_weight dict of decision variables


In [10]:
ing_weight = pulp.LpVariable.dicts("weight kg",
                                     ((i, j) for i in sausage_types for j in ingredients),
                                     lowBound=0,
                                     cat='Continuous')

In [13]:
type(ing_weight)
ing_weight


Out[13]:
dict
Out[13]:
{('economy', 'pork'): weight_kg_('economy',_'pork'),
 ('economy', 'starch'): weight_kg_('economy',_'starch'),
 ('economy', 'wheat'): weight_kg_('economy',_'wheat'),
 ('premium', 'pork'): weight_kg_('premium',_'pork'),
 ('premium', 'starch'): weight_kg_('premium',_'starch'),
 ('premium', 'wheat'): weight_kg_('premium',_'wheat')}

PuLP provides an lpSum vector calculation for the sum of a list of linear expressions.

Whilst we only have 6 decision variables, I will demonstrate how the problem would be constructed in a way that could be scaled up to many variables using list comprehensions.


In [14]:
# Objective Function
model += (
    pulp.lpSum([
        4.32 * ing_weight[(i, 'pork')]
        + 2.46 * ing_weight[(i, 'wheat')]
        + 1.86 * ing_weight[(i, 'starch')]
        for i in sausage_types])
)

Now we add our constraints, bear in mind again here how the use of list comprehensions allows for scaling up to many ingredients or sausage types


In [15]:
# Constraints
# 350 economy and 500 premium sausages at 0.05 kg
model += pulp.lpSum([ing_weight['economy', j] for j in ingredients]) == 350 * 0.05
model += pulp.lpSum([ing_weight['premium', j] for j in ingredients]) == 500 * 0.05

# Economy has >= 40% pork, premium >= 60% pork
model += ing_weight['economy', 'pork'] >= (
    0.4 * pulp.lpSum([ing_weight['economy', j] for j in ingredients]))

model += ing_weight['premium', 'pork'] >= (
    0.6 * pulp.lpSum([ing_weight['premium', j] for j in ingredients]))

# Sausages must be <= 25% starch
model += ing_weight['economy', 'starch'] <= (
    0.25 * pulp.lpSum([ing_weight['economy', j] for j in ingredients]))

model += ing_weight['premium', 'starch'] <= (
    0.25 * pulp.lpSum([ing_weight['premium', j] for j in ingredients]))

# We have at most 30 kg of pork, 20 kg of wheat and 17 kg of starch available
model += pulp.lpSum([ing_weight[i, 'pork'] for i in sausage_types]) <= 30
model += pulp.lpSum([ing_weight[i, 'wheat'] for i in sausage_types]) <= 20
model += pulp.lpSum([ing_weight[i, 'starch'] for i in sausage_types]) <= 17

# We have at least 23 kg of pork to use up
model += pulp.lpSum([ing_weight[i, 'pork'] for i in sausage_types]) >= 23

In [16]:
# Solve our problem
model.solve()
pulp.LpStatus[model.status]


Out[16]:
1
Out[16]:
'Optimal'

In [18]:
for var in ing_weight:
    var_value = ing_weight[var].varValue
    print("The weight of {0} in {1} sausages is {2} kg".format(var[1], var[0], var_value))


The weight of pork in economy sausages is 7.0 kg
The weight of wheat in premium sausages is 2.75 kg
The weight of starch in economy sausages is 4.375 kg
The weight of pork in premium sausages is 16.0 kg
The weight of wheat in economy sausages is 6.125 kg
The weight of starch in premium sausages is 6.25 kg

In [19]:
total_cost = pulp.value(model.objective)

print("The total cost is €{} for 350 economy sausages and 500 premium sausages".format(round(total_cost, 2)))


The total cost is €140.96 for 350 economy sausages and 500 premium sausages

In [ ]: