Linear Programming with Python and PuLP – Part 6

  • Post author:
  • Post category:Python

Introduction to Linear Programming with Python – Part 6

Mocking conditional statements using binary constraints

In part 5, I mentioned that in some cases it is possible to construct conditional statements using binary constraints.

We will explore not only conditional statements using binary constraints, but combining them with logical operators, ‘and’ and ‘or’.

First we’ll work through some theory, then a real world example as an extension of part 5’s example at the end.

Conditional statement

To start simply, if we have the binary constraint x1 and we want:

if x1 == 1:
    y1 == 1
elif x1 == 0:
    y1 == 0

We can achieve this easily using the following constraint:

y1 == x1

However, if we wanted the opposite:

if x1 == 1:
    y1 == 0
elif x1 == 0:
    y1 == 1

Given that they most both be 1 or 0, we just need the following constraint:

x1 + y1 == 1

Logical ‘AND’ operator

Now for something a little more complex, we can coerce a particular binary constraint to be 1 based on the states of 2 other binary constraints.

If we have binary constraints x1 and x2 and we want to achieve the following:

if x1 == 1 and x2 == 0:
    y1 == 1
else:
    y1 == 0

So that $y_1$ is only 1 in the case that x1 is 1 and x2 is 0. We can use the following 3 constraints to achieve this:

[
y1 >= x1 - x2,
y1 <= x1,
y1 <= (1 - x2)
]

We’ll take a moment to deconstruct this. In our preferred case that x1 = 1 and x2 = 0, the three statments resolve to:

  • y1 ≥ 1
  • y1 ≤ 1
  • y1 ≤ 1

The only value of $y_1$ that fulfils each of these is 1.

In any other case, however, y1 will be zero. Let’s take another example, say x1 = 0 and x2 = 1. This resolves to:

  • y1 ≥ -1
  • y1 ≤ 0
  • y1 ≤ 0

Given that y1 is a binary variable and must be 0 or 1, the only value of y1 that can fulfil each of these is 0.

You can construct 3 constraints so that y1 is equal to 1, only in the case you’re interested in out of the 4 following options:

  • x1 = 1 and x2 = 1
  • x1 = 1 and x2 = 0
  • x1 = 0 and x2 = 1
  • x1 = 0 and x2 = 0

I have created a function for exactly this purpose to cover all cases:

In [1]:
def make_io_and_constraint(y1, x1, x2, target_x1, target_x2):
    """
    Returns a list of constraints for a linear programming model
    that will constrain y1 to 1 when
    x1 = target_x1 and x2 = target_x2; 
    where target_x1 and target_x2 are 1 or 0
    """
    binary = [0,1]
    assert target_x1 in binary
    assert target_x2 in binary
    
    if IOx1 == 1 and IOx2 == 1:
        return [
            y1 >= x1 + x2 - 1,
            y1 <= x1,
            y1 <= x2
        ]
    elif IOx1 == 1 and IOx2 == 0:
        return [
            y1 >= x1 - x2,
            y1 <= x1,
            y1 <= (1 - x2)
        ]
    elif IOx1 == 0 and IOx2 == 1:
        return [
            y1 >= x2 - x1,
            y1 <= (1 - x1),
            y1 <= x2
        ]
    else:
        return [
            y1 >= - (x1 + x2 -1),
            y1 <= (1 - x1),
            y1 <= (1 - x2)
        ]

Logical ‘OR’ operator

This is all well and good for the ‘and’ logical operator. What about the ‘or’ logical operator.

If we would like the following:

if x1 == 1 or x2 == 1:
    y1 == 1
else:
    y1 == 0

We can use the following linear constraints:

y1 <= x1 + x2
y1 * 2 >= x1 + x2

So that:

  • if x1 is 1 and x2 is 1:
    • y1 ≤ 2
    • 2y1 ≥ 2
    • y1 must equal 1
  • if x1 is 1 and x2 is 0:
    • y1 ≤ 1
    • 2y1 ≥ 1
    • y1 must equal 1
  • if x1 is 0 and x2 is 1:
    • y1 ≤ 1
    • 2y1 ≥ 1
    • y1 must equal 1
  • if x1 is 0 and x2 is 0:
    • y1 ≤ 0
    • 2y1 ≥ 0
    • y1 must equal 0

Again, we’ll consider the alternative option:

if x1 == 0 or x2 == 0:
    y1 == 1
else:
    y1 == 0

We can use the following linear constraints:

y1 * 2 <= 2 - (x1 + x2)
y1 >= 1 - (x1 + x2)

An Example – Scheduling Example Extended

In our last example, we explored the scheduling of 2 factories.

Both factories had 2 costs:

  • Fixed Costs – Costs incurred while the factory is running
  • Variable Costs – Cost per unit of production

We’re going to introduce a third cost – Start up cost.

This will be a cost incurred by turning on the machines at one of the factories.

In this example, our start-up costs will be:

  • Factory A – €20,000
  • Factory B – €400,000

Let’s start by reminding ourselves of the input data.

In [2]:
import pandas as pd
import pulp

factories = pd.DataFrame.from_csv('csv/factory_variables.csv', index_col=['Month', 'Factory'])
factories
Out[2]:
Max_Capacity Min_Capacity Variable_Costs Fixed_Costs
Month Factory
1 A 100000 20000 10 500
B 50000 20000 5 600
2 A 110000 20000 11 500
B 55000 20000 4 600
3 A 120000 20000 12 500
B 60000 20000 3 600
4 A 145000 20000 9 500
B 100000 20000 5 600
5 A 160000 20000 8 500
B 0 0 0 0
6 A 140000 20000 8 500
B 70000 20000 6 600
7 A 155000 20000 5 500
B 60000 20000 4 600
8 A 200000 20000 7 500
B 100000 20000 6 600
9 A 210000 20000 9 500
B 100000 20000 8 600
10 A 197000 20000 10 500
B 100000 20000 11 600
11 A 80000 20000 8 500
B 120000 20000 10 600
12 A 150000 20000 8 500
B 150000 20000 12 600
In [3]:
demand = pd.DataFrame.from_csv('csv/monthly_demand.csv', index_col=['Month'])
demand
Out[3]:
Demand
Month
1 120000
2 100000
3 130000
4 130000
5 140000
6 130000
7 150000
8 170000
9 200000
10 190000
11 140000
12 100000

We’ll begin by defining our decision variables, we have an additional binary variable for switching on the factory.

In [4]:
# Production
production = pulp.LpVariable.dicts("production",
                                     ((month, factory) for month, factory in factories.index),
                                     lowBound=0,
                                     cat='Integer')

# Factory Status, On or Off
factory_status = pulp.LpVariable.dicts("factory_status",
                                     ((month, factory) for month, factory in factories.index),
                                     cat='Binary')

# Factory switch on or off
switch_on = pulp.LpVariable.dicts("switch_on",
                                    ((month, factory) for month, factory in factories.index),
                                    cat='Binary')

We instantiate our model and define our objective function, including start up costs

In [5]:
# Instantiate the model
model = pulp.LpProblem("Cost minimising scheduling problem", pulp.LpMinimize)

# Select index on factory A or B
factory_A_index = [tpl for tpl in factories.index if tpl[1] == 'A']
factory_B_index = [tpl for tpl in factories.index if tpl[1] == 'B']

# Define objective function
model += pulp.lpSum(
    [production[m, f] * factories.loc[(m, f), 'Variable_Costs'] for m, f in factories.index]
    + [factory_status[m, f] * factories.loc[(m, f), 'Fixed_Costs'] for m, f in factories.index]
    + [switch_on[m, f] * 20000 for m, f in factory_A_index]
    + [switch_on[m, f] * 400000 for m, f in factory_B_index]
)

Now we begin to build up our constraints as in Part 5

In [6]:
# Production in any month must be equal to demand
months = demand.index
for month in months:
    model += production[(month, 'A')] + production[(month, 'B')] == demand.loc[month, 'Demand']

# Production in any month must be between minimum and maximum capacity, or zero.
for month, factory in factories.index:
    min_production = factories.loc[(month, factory), 'Min_Capacity']
    max_production = factories.loc[(month, factory), 'Max_Capacity']
    model += production[(month, factory)] >= min_production * factory_status[month, factory]
    model += production[(month, factory)] <= max_production * factory_status[month, factory]

# Factory B is off in May
model += factory_status[5, 'B'] == 0
model += production[5, 'B'] == 0

But now we want to add in our constraints for switching on.

A factory switches on if:

  • It is off in the previous month (m-1)
  • AND it on in the current month (m).

As we don’t know if the factory is on before month 0, we’ll assume that the factory has switched on if it is on in month 1.

In [7]:
for month, factory in factories.index:
    # In month 1, if the factory ison, we assume it turned on
    if month == 1:
        model += switch_on[month, factory] == factory_status[month, factory]
    
    # In other months, if the factory is on in the current month AND off in the previous month, switch on = 1
    else:
        model += switch_on[month, factory] >= factory_status[month, factory] - factory_status[month-1, factory]
        model += switch_on[month, factory] <= 1 - factory_status[month-1, factory]
        model += switch_on[month, factory] <= factory_status[month, factory]
        

We’ll then solve our model

In [8]:
model.solve()
pulp.LpStatus[model.status]
Out[8]:
'Optimal'
In [9]:
output = []
for month, factory in production:
    var_output = {
        'Month': month,
        'Factory': factory,
        'Production': production[(month, factory)].varValue,
        'Factory Status': factory_status[(month, factory)].varValue,
        'Switch On': switch_on[(month, factory)].varValue
    }
    output.append(var_output)
output_df = pd.DataFrame.from_records(output).sort_values(['Month', 'Factory'])
output_df.set_index(['Month', 'Factory'], inplace=True)
output_df
Out[9]:
Factory Status Production Switch On
Month Factory
1 A 1 70000 1
B 1 50000 1
2 A 1 45000 0
B 1 55000 0
3 A 1 70000 0
B 1 60000 0
4 A 1 30000 0
B 1 100000 0
5 A 1 140000 0
B 0 0 0
6 A 1 60000 0
B 1 70000 1
7 A 1 90000 0
B 1 60000 0
8 A 1 70000 0
B 1 100000 0
9 A 1 100000 0
B 1 100000 0
10 A 1 170000 0
B 1 20000 0
11 A 1 80000 0
B 1 60000 0
12 A 1 100000 0
B 0 0 0

Interestingly, we see that it now makes economic sense to keep factory B on after it turns off in month 5 up until month 12.

Previously, we had the case that it was not economic to run factory B in month 10, but as there is now a significant cost to switching off and back on, the factory runs through month 10 at its lowest capacity (20,000 units).


For those interested in using my function defined above (make_io_and_constraint). Instead of:

model += switch_on[month, factory] >= factory_status[month, factory] - factory_status[month-1, factory]
model += switch_on[month, factory] <= 1 - factory_status[month-1, factory]
model += switch_on[month, factory] <= factory_status[month, factory]

You could write:

for constraint in make_io_and_constraint(switch_on[month, factory], 
                                        factory_status[month, factory], 
                                        factory_status[month-1, factory], 0, 1):
    model += constriant

Introduction
Part 1 – Introduction to Linear Programming
Part 2 – Introduction to PuLP
Part 3 – Real world examples – Resourcing Problem
Part 4 – Real world examples – Blending Problem
Part 5 – Using PuLP with pandas and binary constraints to solve a scheduling problem
Part 6 – Mocking conditional statements using binary constraints