Linear Programming With Python

a linear programming problem may be defined as the problem of maximizing or minimizing a linear function subject to linear constraints. the constraints may be equalities or inequalities.

here is a simple example: find numbers x1 and x2 that maximize the sum x1 + x2 subject to the constraints x1 ≥ 0, x2 ≥ 0, and

x1 + 2x2 ≤ 4
4x1 + 2x2 ≤ 12
−x1 + x2 ≤ 1

the first two constraints, x1 ≥ 0 and x2 ≥ 0 are called nonnegativity constraints . the other constraints are then called the main constraints . the function to be maximized (or minimized) is called the objective function . here, the objective function is x1 + x2.

two classes of problems, called here the standard maximum problem and the standard minimum problem , play a special role. in these problems, all variables are constrained to be nonnegative, and all main constraints are inequalities.

standard maximum problem

we are given an m-vector, an n-vector, and an m x n matrix of real numbers

find an n-vector, to maximize

subject to the constraints

and

standard minimum problem

find an n-vector, to minimize subject to the constraints and

example

the example at the beginning of this post corresponds to the standard maximum problem, where:

the solutions are:

from pulp import *
from fractions import fraction

prob = lpproblem("example of standard maximum problem",lpmaximize)

# nonnegativity constraints
x1=lpvariable("x1",0)
x2=lpvariable("x2",0)

# objective function
prob += x1 + x2, "maximum value of x1 + x2"

# main constraints
prob += x1 + 2*x2 <= 4, "constraint 1"
prob += 4*x1 + 2*x2 <= 12, "constraint 2"
prob += -x1 + x2 <= 1, "constraint 3"

# the problem is solved using pulp's choice of solver
prob.solve()

# status of the solution
print(f"status: {lpstatus[prob.status]}")

for v in prob.variables():
print(f"{v.name} = {str(fraction(v.varvalue).limit_denominator())}")

# maximum value of the objective function
print(f"max (x1 + x2) = {str(fraction(value(prob.objective)).limit_denominator())}")

status: optimal
x1 = 8/3
x2 = 2/3
max (x1 + x2) = 10/3

the importance of the standard problem derives from the fact that all linear programming problems can be converted to standard form.

what is more, a minimum problem can be changed to a maximum problem by:

• multiplying the objective function by -1
• multiplying the constraint inequalities by -1 and reversing the inequalities

let's prove it by changing the previous example to a minimum problem:

from pulp import *
from fractions import fraction

prob = lpproblem("example of standard minimum problem",lpminimize)

# nonnegativity constraints
x1=lpvariable("x1",0)
x2=lpvariable("x2",0)

# objective function
prob += -x1 - x2, "minimum value of -x1 - x2"

# main constraints
prob += -x1 - 2*x2 >= -4, "constraint 1"
prob += -4*x1 - 2*x2 >= -12, "constraint 2"
prob += x1 - x2 >= -1, "constraint 3"

# the problem is solved using pulp's choice of solver
prob.solve()

# status of the solution
print(f"status: {lpstatus[prob.status]}")

for v in prob.variables():
    print(f"{v.name} = {str(fraction(v.varvalue).limit_denominator())}")

# maximum value of the objective function
print(f"min (-x1 - x2) = {str(fraction(value(prob.objective)).limit_denominator())}")

status: optimal
x1 = 8/3
x2 = 2/3
min (-x1 - x2) = -10/3

duality

for every linear program there is a dual linear program. the dual of the standard maximum problem (as defined above) is defined to be the standard minimum problem:

find an m-vector, to minimize subject to the constraints and

therefore, the dual problem of the initial example is:

from pulp import *
from fractions import fraction

prob = lpproblem("dual problem",lpminimize)

# nonnegativity constraints
y1=lpvariable("y1",0)
y2=lpvariable("y2",0)
y3=lpvariable("y3",0)

# objective function
prob += 4*y1 + 12*y2 + y3, "minimum value of 4*y1 + 12*y2 + y3"

# main constraints
prob += y1 + 4*y2 - y3 >= 1, "constraint 1"
prob += 2*y1 + 2*y2 + y3 >= 1, "constraint 2"

# the problem is solved using pulp's choice of solver
prob.solve()

# status of the solution
print(f"status: {lpstatus[prob.status]}")

for v in prob.variables():
    print(f"{v.name} = {str(fraction(v.varvalue).limit_denominator())}")

# maximum value of the objective function
print(f"min (4*y1 + 12*y2 + y3) = {str(fraction(value(prob.objective)).limit_denominator())}")

status: optimal
y1 = 1/3
y2 = 1/6
y3 = 0
min (4*y1 + 12*y2 + y3) = 10/3

as the reader must have noticed, the number of main constraints of the standard problem equals the number of nonnegative constraints of its dual. according to the equilibrium theorem , strict inequality in a constraint in a standard problem implies that the complementary constraint in the dual is satisfied with equality, and viceversa.

let's verify that assertion:

print([x1 + 2*x2 == 4, 4*x1 + 2*x2 == 12, -x1 + x2 == 1] == [not x for x in [y1.varvalue == 0, y2.varvalue == 0, y3.varvalue == 0]])
print([y1 + 4*y2 - y3 == 1, 2*y1 + 2*y2 + y3 == 1] == [not x for x in [x1.varvalue == 0, x2.varvalue == 0]])

true
true

although, in this case, the value of the objective function corresponding to the vectors x and y is the same, this is not always true. it can be proved that, when this is true, the solution is optimal for both problems.

more generally, the values of the standard objective function (compatible with its constraints) are always ≤ than the values of the dual objective function (compatible with its constraints). in other words, the dual problem provides an upper bound to the optimal value of the original problem.