Easter 98  (April 16, 1998)  Richard Weber

D2: OPTIMIZATION

Contents

Schedules iii

Notation iii

Index iv

1 Preliminaries 1 1.1 Optimization under constraints . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Representation of constraints . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.3 Geometry of linear programming . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 1.4 Extreme points and optimality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3 1.5 *Efficiency of algorithms* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4

2 The Solution of LP Problems 5 2.1 Basic solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2.2 Primal-dual relationships . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6

3 The Simplex Method 10 3.1 Preview of the simplex algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10 3.2 The simplex algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12

4 The Simplex Tableau 14 4.1 Choice of pivot column . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 4.2 Initialization: the two-phase method . . . . . . . . . . . . . . . . . . . . . . . . . . . 15 4.3 Sensitivity: shadow prices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17

5 Lagrangian Methods 19 5.1 The Lagrangian . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 5.2 The Lagrangian sufficiency theorem . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 5.3 Strategy to solve problems with the Lagrangian sufficiency theorem . . . . . . . . . . 21 5.4 Example: use of the Lagrangian sufficiency theorem . . . . . . . . . . . . . . . . . . . 21

6 The Lagrangian Dual 23 6.1 Example: further use of the Lagrangian sufficiency theorem . . . . . . . . . . . . . . . 23 6.2 The Lagrangian dual problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 6.3 The dual problem for LP . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25 6.4 The weak duality theorem in the case of LP . . . . . . . . . . . . . . . . . . . . . . . 26

7 Shadow Prices and Lagrangian Necessity 27 7.1 Sufficient conditions for optimality . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 7.2 Shadow prices . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 7.3 Lagrangian necessity . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 7.4 Convexity of the feasible set in the case of LP . . . . . . . . . . . . . . . . . . . . . . 30 8 Algebra of Linear Programming 31

8.1 Basic feasible solutions and extreme points . . . . . . . . . . . . . . . . . . . . . . . . 31 8.2 Algebra of the simplex method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

9 Two Person Zero-Sum Games 36 9.1 Games with a saddle-point . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 9.2 Example: Two-finger Morra, a game without a saddle-point . . . . . . . . . . . . . . 37 9.3 Determination of an optimal strategy . . . . . . . . . . . . . . . . . . . . . . . . . . . 37 9.4 Example: Colonel Blotto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

10 Maximal Flow in a Network 40 10.1 Max-flow/min-cut theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 10.2 Ford-Fulkerson algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 42 10.3 Minimal cost circulations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43

11 Minimum Cost Circulation Problems 44 11.1 Sufficient conditions for a minimal cost circulation . . . . . . . . . . . . . . . . . . . . 44 11.2 Max-flow as a minimum cost circulation problem . . . . . . . . . . . . . . . . . . . . 45 11.3 The transportation problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46

12 Transportation and Transshipment Problems 48 12.1 The transportation algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 12.2 *Simplex-on-a-graph* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 12.3 Example: optimal power generation and distribution . . . . . . . . . . . . . . . . . . 52

Books

Bazaraa, M., Jarvis, J. and Sherali, H. Linear Programming and Network Flows, second edition, 1990, Wiley.

Luenberger, D. Introduction to Linear and Non-Linear Programming, second edition, 1984, Addison-Wesley.

BJS is a paperback and full of lots of examples and reasonably comprehensive discussions of the history, recent results, etc. It is a little heavy-handed about some of the theory. Luenberger is nicer on the theory, but may be harder going.

Notes, examples sheet and overheads

An up-to-date postscript copy of these notes, the examples sheet and any overhead slides used in lectures can be downloaded from the WWW site:

http://www.statslab.cam.ac.uk/~rrw1/opt/index.html.

Acknowledgement

The original version of these notes were written by Colin Sparrow in 1994.

Schedules

General formulation of constrained problems; the Lagrangian sufficiency theorem. Interpretation of Lagrange multipliers as shadow prices. Examples. Convexity of feasible region; sufficiency of extreme points. Standardization of problems, slack variables, equivalence of extreme points and basic solutions. The primal simplex algorithm, artificial variables, the two phase method. Practical use of the algorithm; the tableau. Examples. The dual linear problem, duality theorem in a standard case, complementary slackness, dual variables and their interpretation as shadow prices. Relationship of the primal simplex algorithm to dual problem. Two person zero-sum games. The Ford-Fulkerson algorithm and the max-flow min-cut theorems in the rational case. Network flows with costs, the transportation algorithm, relationship of dual variables with nodes. Examples. Conditions for optimality in more general networks; *the simplex-on-a-graph algorithm*. *Efficiency of algorithms*. The formulation of simple practical and combinatorial problems as linear programming or network problems.

Notation

n, m  numbers of decision variables and functional constraints
x  decision variable in a primal problem, x ∈ ℝ^{n}
λ  decision variable in the dual problem, λ ∈ ℝ^{m}
f(x)  objective function in a primal problem
x ∈ X  regional constraints, X ⊆ ℝ^{n}
g(x) ≤ b, g(x) = b functional constraints, g : ℝ^{n} → ℝ^{m}
z  slack variable, z ∈ ℝ^{m}
Ax ≤ b, Ax + z = b linear constraints
c^{T} x  linear objective function
L(x, λ)  Lagrangian, L(x, λ) = f(x) − λ^{T} (g(x) − b)
λ ∈ Y  Y = {λ : min_{x∈X} L(x, λ) > −∞}.
B, N  sets of indices of basic and non-basic components of x.
A, p, q  pay-off matrix and decision variables in a matrix game
x_{ij}  flow on arc (i, j)
v, C(S, S̄)  value of a flow, value of cut (S, S̄)
c_{ij}^{−}, c_{ij}^{+}  minimum/maximum allowed flows on arc (i, j)
d_{ij}  costs per unit flow on arc (i, j)
s_{i}, d_{j}  source and demands amounts in transportation problem
λ_{i}, μ_{j}  node numbers in transportation algorithm

Index

anti-symmetric matrix, 39 artificial variables, 15 basic feasible solution, 6 basic solution, 6 basic, 6 basis, 6 capacity, 40 choice of pivot column, 14 circulation problem, minimal cost, 43 circulation, 43 closed, 43 complementary slackness, 7, 20 computational complexity, 4 convex function, 3 convex set, 2 cut, 40 dual feasibility, 20 dual problem, 25 extreme point, 3 feasible circulation, 43 feasible set, 1 Ford-Fulkerson algorithm, 42 functional constraints, 1 integer linear programming (ILP), 4 Lagrange multipliers, 19 Lagrangian dual problem, 24 Lagrangian sufficiency theorem, 20 Lagrangian, 19 linear programming (LP), 2 max-flow/min-cut, 40 minimal cost circulation, 43 mixed strategy, 37 node numbers, 44 non-basic, 6 non-degeneracy assumptions, 32 non-degenerate, 6 pay-off matrix, 36 pivoting, 13 potentials, 44 primal problem, 25 primal/dual theory, 6 regional constraints, 1 revised simplex algorithm, 35 saddle-point, 36 shadow prices, 18, 27 simplex algorithm, 12 simplex tableau, 12 simplex-on-a-graph algorithm, 51 slack variable, 1 spanning tree, 48 surplus variables, 15 tableau, 12 tension, 44 tight, 7 transportation algorithm, 48 transportation problem, 46 tree, 48 two person zero-sum games, 36 two-phase method, 15 value of the game, 38 value of the flow, 40 weak duality, 26

1 Preliminaries

1.1 Optimization under constraints

The general problem we study in this course takes the form

maximize f(x)
subject to x ∈ X
g(x) = b

where

x ∈ ℝ^{n} (n decision variables)
f : ℝ^{n} → ℝ (objective function)
X ⊆ ℝ^{n} (regional constraints)
g : ℝ^{n} → ℝ^{m} (m functional equations)
b ⊆ ℝ^{m}

Note that minimizing f(x) is the same as maximizing −f(x). In lectures we will discuss various examples of constrained optimization problems. We will also talk briefly about the way in which our methods are applied to real-world problems.

1.2 Representation of constraints

Slack variables Sometimes we shall let the constraint be g(x) ≤ b. This form of constraint can be turned into an equality constraint by the addition of a slack variable z. We write

g(x) + z = b, z ≥ 0.

It will frequently be mathematically useful to turn all our constraints into equalities. Regional and functional constraints There is freedom in choosing how to divide the constraints between regional constraints of the form

x ∈ X

and functional constraints of the form

g(x) = b.

Together, these define the feasible set for x. Typically, ‘obvious’ constraints like x ≥ 0 are catered for by defining X in an appropriate way and more complicated constraints, that may change from instance to instance of the problem, are expressed by functional constraints g(x) = b or g(x) ≤ b.

Sometimes the choice is made for mathematical convenience. Methods of solution typically treat region and functional constraints differently.

1.3 Geometry of linear programming

Consider the two problems P and D.

P: maximize x_{1} + x_{2}
subject to x_{1} + 2x_{2} ≤ 6
x_{1} − x_{2} ≤ 3
x_{1}, x_{2} ≥ 0
D: minimize 6λ_{1} + 3λ_{2}
subject to λ_{1} + λ_{2} ≥ 1
2λ_{1} − λ_{2} ≥ 1
λ_{1}, λ_{2} ≥ 0

These are completely linear problems – the objective function and all constraints are linear. The two problems are related, as we shall see. Note that in matrix/vector notation we can write

P: maximize c^{T} x s.t. Ax ≤ b, x ≥ 0,
D: minimize λ^{T} b s.t. A^{T} λ ≥ c, λ ≥ 0.

Concentrate on problem P for the moment. The shaded region is the feasible set defined by the constraints for P.

PSfrag replacements A B C D E F x_{1} = 0

x_{2} = 0
x_{1} + 2x_{2} = 6
x_{1} − x_{2} = 3

x_{1} x_{2}

We now state some definitions and theorems that are at the heart of the theory of linear programming (LP) problems. The theorems are more or less obvious for our example. They will be proved in future lectures.

Definition 1.1 A set S ⊆ ℝ^{n} is a convex set if x, y ∈ S ⟹ λx + (1 − λ)y ∈ S for all x, y ∈ S and 0 ≤ λ ≤ 1.

In other words, the line segment joining x and y lies in S.

Theorem 1.1 The feasible set of a LP problem is convex.

PSfrag replacements

convex not convex

x
y
S
λ = 0
λ = 1
λ = \frac{1}{2}

For functions defined on convex sets we make the following definition.

Definition 1.2 A function f : S → ℝ is a convex function if the set above its graph is convex. Equivalently, λf(x) + (1 − λ)f(y) ≥ f(λx + (1 − λ)y), for all 0 ≤ λ ≤ 1.

PSfrag replacements

convex concave

x
f

For general S and f there may be local minima (maxima) of f over S which are not global minima (maxima). However, if S is a convex set and f is a convex (concave) function then you can be sure that any local minimum (maximum) of f is also a global minimum (maximum). Note that a linear function is both convex and concave, and so all local optima of a linear objective function are also global optima.

1.4 Extreme points and optimality

Notice that for problem P the optimum occurs (regardless of what the linear objective function is) at a ‘corner’ of the feasible set. (In our case at C.)

There may be optima elsewhere on the edge if the objective function is parallel to an edge, but there is always still an optimum at a ‘corner’. This motivates the following definition.

Definition 1.3 We say x is an extreme point of a convex set S if whenever x = λy + (1 − λ)z, for y, z ∈ S, 0 < λ < 1, then x = y = z.

That is, x is not in the interior of any line segment.

PSfrag replacements

feasible set objective function

all solutions on this edge are optimal including the two endpoints

Examples of extreme points of two convex sets PSfrag replacements all boundary points

are extreme

corners are the only extreme points

not extreme S = {x^{2} + y^{2} ≤ 1}

Theorem 1.2 If an LP has a finite optimum it has an optimum at an extreme point of the feasible set.

For LP problems the feasible set will always have a finite number of extreme points (vertices). The feasible set is ‘polyhedral’, though it may be bounded or unbounded. This suggests the following algorithm for solving LPs.

Algorithm:  1. Find all the vertices of the feasible set. 2. Pick the best one.

This will work, but there may be very many vertices. In fact, for Ax ≤ b, x ≥ 0, there can be \binom{n+m}{m} vertices. So if m = n, say, then the number of vertices is of order (2n)^{n}, which increases exponentially in n. This is not a good algorithm!

1.5 *Efficiency of algorithms*

There is an important distinction between those algorithms whose running times (in the worst cases) are exponential functions of ‘problem size’, e.g., (2n)^{n}, and those algorithms whose running times are polynomial functions of problem size, e.g., n^{k}. For example, the problem of finding the smallest number in a list of n numbers is solvable in polynomial-time n by simply scanning the numbers. There is a beautiful theory concerning the computational complexity of algorithms and one of its main messages is that problems solvable in polynomial-time are the ‘easy’ ones.

In fact, there exists a polynomial-time algorithm for general LP problems, but this was not known until 1981. In contrast, no polynomial-time algorithm is known for general integer LP, in which x is restricted to be integer-valued. ILP includes important problems such as bin packing, job-shop scheduling, traveling salesman and many other essentially equivalent problems of a combinatorial nature.

2 The Solution of LP Problems

2.1 Basic solutions

Returning to P. We need a more algebraic (less geometric) characterisation of extreme points. Let us rewrite P with equality constraints, using slack variables.

P: maximize  x_{1} + x_{2}
subject to x_{1} + 2x_{2} + z_{1} = 6
x_{1} − x_{2} + z_{2} = 3
x_{1}, x_{2}, z_{1}, z_{2} ≥ 0

Let us calculate the value of the variables at each of the 6 points marked A–F in our picture of the feasible set for P. The values are:

x_{1} x_{2} z_{1} z_{2} f
A 0 0 6 3 0
B 3 0 3 0 3
C 4 1 0 0 5
D 0 3 0 6 3
E 6 0 0 −3 6
F 0 −3 12 0 −3

At each point there are two zero and two non-zero variables. This is not surprising.

Geometrically: The 4 lines defining the feasible set can be written x_{1} = 0; x_{2} = 0; z_{1} = 0; z_{2} = 0. At the intersection of each pair of lines, two variables are zero.

Algebraically: Constraints Ax + z = b are 2 equations in 4 unknowns. If we choose 2 variables (which can be done in \binom{4}{2} = 6 ways) and set them equal to zero we will be left with two equations in the other two variables. So (provided A and b are ‘nice’) there will be a unique solution for the two non-zero variables.

Instead of calling the slack variables z_{1} and z_{2}, let us call them x_{3} and x_{4} so that we can write P as

P: maximize x_{1} + x_{2}
subject to  Ax = b
x ≥ 0

Note we have to extend A to ( \begin{matrix}1 & 2 & 1 & 0 \\ 1 & −1 & 0 & 1 \end{matrix} ) and x to ( \begin{matrix}x_{1} \\ x_{2} \\ x_{3} \\ x_{4} \end{matrix} ). In this version

of P, A is (m × n) with n > m and there are m equations in n > m unknowns. We can therefore choose (n − m) variables in \binom{n}{m} ways. Set them to zero. There will be a unique solution to Ax = b for the remaining m variables (provided A and b are ‘nice’).

Definition 2.1 A basic solution to Ax = b is a solution with at least n − m zero variables. The solution is a non-degenerate basic solution if exactly n−m variables are zero. The choice of the m non-zero variables is called the basis. Variables in the basis are called basic; the others are called non-basic.

Definition 2.2 If a basic solution satisfies x ≥ 0 then it is called a basic feasible solution.

So A–F are basic solutions (and non-degenerate) and A–D are basic feasible solutions.

Theorem 2.1 Basic feasible solutions ≡ extreme points of the feasible set.

Taking Theorems 1.2 and 2.1 together,

Theorem 2.2 If an LP has a finite solution then there is an optimal basic feasible solution.

So we can do algebra instead of drawing a picture (which is good for computer and good for us if there are many variables and constraints.) Our previous (and foolish) algorithm can be restated:

Algorithm:

1. Find all the basic solutions.

2. Test to see which are feasible.

3. Choose the best basic feasible solution.

Unfortunately it is not usually easy to know which basic solutions will turn out to be feasible before calculating them. Hence, even though there are often considerably fewer basic feasible solutions we will still need to calculate all \binom{n}{m} basic solutions.

2.2 Primal-dual relationships

Now let us look at problem D which can be written, after introducing slack variables v_{1} and v_{2} as

D: minimize  6λ_{1} + 3λ_{2}
subject to λ_{1} + λ_{2} − v_{1} = 1
2λ_{1} − λ_{2} − v_{2} = 1
λ_{1}, λ_{2}, v_{1}, v_{2} ≥ 0

PSfrag replacements A B C D E F λ_{2} λ_{1}

λ_{1} = 0
λ_{2} = 0
v_{1} = 0
v_{2} = 0

The value of the variables, etc., at the points A–F in P (as above) and D are:

P: \begin{aligned}& x_{1} & x_{2} & z_{1} & z_{2} & f \\ A & 0 & 0 & 6 & 3 & 0 \\ B & 3 & 0 & 3 & 0 & 3 \\ C & 4 & 1 & 0 & 0 & 5 \\ D & 0 & 3 & 0 & 6 & 3 \\ E & 6 & 0 & 0 & −3 & 6 \\ F & 0 & −3 & 12 & 0 & −3 \end{aligned}  D: \begin{aligned}& v_{1} & v_{2} & λ_{1} & λ_{2} & f \\ A & −1 & −1 & 0 & 0 & 0 \\ B & 0 & −2 & 0 & 1 & 3 \\ C & 0 & 0 & \frac{2}{3} & \frac{1}{3} & 5 \\ D & −\frac{1}{2} & 0 & \frac{1}{2} & 0 & 3 \\ E & 0 & 1 & 1 & 0 & 6 \\ F & −2 & 0 & 0 & −1 & −3 \end{aligned}

Observe, that for D, as for P above, there are two zero and two non-zero variables at each intersection (basic solution). C and E are feasible for D. The optimum is at C with optimum value 5 (assuming we are minimizing and the other basic solutions are not feasible.)

We make the following observations by comparing lists of basic solutions for P and D. (Greater explanation is given later.)

1. For each solution for P there is a corresponding basic solution for D. [Labels A–F have been chosen so that corresponding solutions have the same labels.] Each pair

(a) has the same value of the objective function.

(b) satisfies complementary slackness, i.e., x_{i}v_{i} = 0, λ_{i}z_{i} = 0,

so for each corresponding pair,

P  D
variables x  constraints
x_{i} basic (x_{i} ≠ 0)  ⟺ constraint: tight (v_{i} = 0)
x_{i} non-basic (x_{i} = 0) ⟺ constraint: slack (v_{i} ≠ 0)
constraints  variables λ
constraint: tight (z_{i} = 0) ⟺  λ_{i} basic (λ_{i} ≠ 0)
constraint: slack (z_{i} ≠ 0) ⟺ λ_{i} non-basic (λ_{i} = 0)

(These conditions determine which basic solutions in P and D are paired; the implications go both ways because in this example all basic solutions are nondegenerate.)

2. There is only one pair that is feasible for both P and D, and that solution is C, which is optimal, with value 5, for both.

3. For any x feasible for P and λ feasible for D we have c^{⊤}x ≤ b^{⊤}λ with equality if and only if x, λ are optima for P and D.

This correspondence between P and D is so symmetrical (and pretty) that it feels as though it ought to be obvious why it works. We will do the necessary proofs for most of it later, but be reassured that it is not all obvious. In particular, it is hard to prove that P and D have the same solution. But let us write down the theorems. Let P (primal) and D (dual) be the LPs

P: max c^{⊤}x s.t. Ax ≤ b, x ≥ 0,
D: min b^{⊤}λ s.t. A^{⊤}λ ≥ c, λ ≥ 0.

Then

Theorem 2.3 (weak duality) If x is feasible for P and λ is feasible for D then c^{⊤}x ≤ b^{⊤}λ. (In particular, if one problem is feasible then the other is bounded.)

Theorem 2.4 (sufficient conditions for optimality) If x is feasible for P and λ is feasible for D, and x, λ satisfy complementary slackness, then x is optimal for P and λ is optimal for D. Furthermore c^{⊤}x = λ^{⊤}b.

Theorem 2.5 (strong duality: necessary conditions for optimality) If both P and D are feasible (each has at least one feasible solution), then ∃ x, λ satisfying the conditions of Theorem 2.4 above.

Remarks.

1. We will prove Theorems 2.3 and 2.4 in quite a general way. Theorem 2.5 is harder and we will only sketch the general proof.

2. Note that if we write D as

max − b^{⊤}λ s.t. − A^{⊤}λ ≤ −c, λ ≥ 0,

and now construct the dual of D (in the same way as we obtained D from P), we recover P. For LP problems the dual of the dual is always the primal.

3. Theorems 2.3–2.5, and the whole idea of constructing the dual problem, work in quite general convex optimization problems (not just LP problems). We shall see some of this. It is a pretty theory.

4. Why do we care about the dual problem instead of just getting on with the solution of P?

(a) The statement of optimality conditions in terms of P and D in Theorem 2.4 is powerful and natural.

(b) It is sometimes easier to solve D than P (and they have the same solution).

(c) For some problems it is natural to consider both P and D together (e.g., two person zero-sum games).

(d) The theory is helpful in constructing algorithms. Theorem 2.4 says that for optimality we need three things: primal feasibility, dual feasibility and complementary slackness. Some algorithms start with solutions that are primal feasible and work towards dual feasibility. Others start with dual feasibility. Yet others (in particular network problems) alternately look at the primal and dual variables and move towards feasibility for both at once.

3 The Simplex Method

3.1 Preview of the simplex algorithm

Let us look very closely at problem P and see if we can construct an algorithm that behaves as follows.

Simplex algorithm

1. Start with a basic feasible solution.

2. Test — is it optimal?

3. If YES — stop.

4. If NO, move to ‘adjacent’ and better b.f.s. Return to 2.

We need to pick a b.f.s. to start. Let us take vertex A.

x_{1} = x_{2} = 0; z_{1} = 6, z_{2} = 3.

[Even for very large problems it is easy to pick a b.f.s. provided the original constraints are m constraints in n variables, Ax ≤ b with b ≥ 0. Once we add slack variables Ax + z = b we have n + m variables and m constraints. If we pick x = 0, z = b this is a b.f.s. More about picking the initial b.f. solutions in other cases later.]

Now we can write problem P as:

x_{i}, z_{i} ≥ 0
x_{1} +2x_{2} +z_{1}  = 6 (1)
x_{1} − x_{2}  +z_{2} = 3 (2)
max  x_{1} + x_{2}  = f (0)

The peculiar arrangement on the page is deliberate. Now it is obvious that A is not optimal because,

1. At A, x_{1} = x_{2} = 0; z_{1} = 6, z_{2} = 3.

2. From the form of the objective function we see that increasing either x_{1} or x_{2} will improve the solution.

3. From (1) we see that it is possible to increase x_{1} to 6 and decrease z_{1} to 0 without violating this equation or making any variable infeasible.

4. From (2) we see that it is possible to increase x_{1} to 3 and decrease z_{2} to 0 before any variable becomes infeasible.

5. Taking (1) and (2) together, then, we can increase x_{1} to 3, decrease z_{1} to 3 and decrease z_{2} to 0, preserving equality in the two constraints, preserving feasibility and increasing the objective function.

That is, we should move to the b.f.s.

x_{1} = 3, x_{2} = 0, z_{1} = 3, z_{2} = 0, (f = 3),

which is vertex B. Note that one variable (x_{1}) has entered the basis and one (z_{2}) has left; i.e., we have moved to an ‘adjacent’ vertex. Why was this easy to see?

1. The objective function f was written in terms of the non-basic variables (x_{1} = x_{2} = 0), so it was easy to see that increasing one of them would improve the solution.

2. Each basic variable (z_{1}, z_{2}) appeared just once in one constraint, so we could consider the effect of increasing x_{1} on each basic variable separately when deciding how much we could increase x_{1}.

This suggests we try and write the problem so that the conditions above hold at B, our new b.f.s.. We can do this by adding multiples of the second equation to the others (which is obviously allowed as we are only interested in variables satisfying the constraints.)

So P can be written,

x_{i}, z_{i} ≥ 0
(1) − (2)  3x_{2} +z_{1} −z_{2} = 3  (1)′
(2)  x_{1} − x_{2}  +z_{2} = 3  (2)′
(0) − (2)  2x_{2}  −z_{2} = f − 3 (0)′

This form of P (equivalent to the original) is what we wanted.

1. The objective function is written in terms of the non-basic variables x_{2}, z_{2}.

2. Basic variables x_{1}, z_{1} each appear just once in one constraint.

The next step is now easy. Remember we are at B:

x_{1} = 3, x_{2} = 0, z_{1} = 3, z_{2} = 0, (f = 3).

1. From the objective function it is obvious that increasing x_{2} is good, whereas increasing z_{2} would be bad.

2. Equation (1)′ shows that we can increase x_{2} to 1 and decrease z_{1} to 0.

3. Equation (2)′ doesn’t impose any restriction on how much we can increase x_{2} (we just would need to increase x_{1} also.)

4. Thus we can increase x_{2} to 1, while decreasing z_{1} to 0 and increasing x_{1} to 4.

So we move to vertex C:

x_{1} = 4, x_{2} = 1, z_{1} = 0, z_{2} = 0, (f = 5).

Now, rewriting the problem again into the desired form for vertex C we obtain

x_{i}, z_{i} ≥ 0
\frac{1}{3}(1)′  x_{2} +\frac{1}{3}z_{1} −\frac{1}{3}z_{2} = 1
(2)′ + \frac{1}{3}(1)′ x_{1} +\frac{1}{3}z_{1} +\frac{2}{3}z_{2} = 4
(0)′ − \frac{2}{3}(1)′  −\frac{2}{3}z_{1} −\frac{1}{3}z_{2} = f − 5

Now it is clear that we have reached the optimum since

1. We know x_{1} = 4, x_{2} = 1, z_{1} = 0, z_{2} = 0 is feasible.

2. We know that for any feasible x, z we have f = 5 − \frac{2}{3}z_{1} − \frac{1}{3}z_{2} ≤ 5. So clearly our solution (with z_{1} = z_{2} = 0) is the best possible.

3.2 The simplex algorithm

The procedure just described for problem P can be formalised. Rather than writing out all the equations each time we write just the coefficients in a table known as the simplex tableau. To repeat what we have just done, we would write:–

x_{1} x_{2} z_{1} z_{2} a_{i0}
z_{1} basic 1 2 1 0 6
z_{2} basic 1 −1 0 1 3
a_{0j}  1 1 0 0 0

If we label the coefficients in the body of the table (a_{ij}), the right hand sides of the equations (a_{i0}), the coefficients in the expression for the objective function as (a_{0j}) and the value of the objective function −a_{00}, so the tableau contains

\begin{matrix}(a_{ij}) & a_{i0} \\ a_{0j} & a_{00} \end{matrix}

The algorithm is

1. Choose a variable to enter the basis. Pick a j such that a_{0j} > 0. (The variable corresponding to column j will enter the basis.) If all a_{0j} ≤ 0, j ≥ 1, then the current solution is optimal. We picked j = 1, so x_{1} is entering the basis.

2. Find the variable to leave the basis. Choose i to minimize a_{i0}/a_{ij} from the set {i : a_{ij} > 0}. If a_{ij} ≤ 0 for all i then the problem is unbounded (see examples sheet) and the objective function can be increased without limit. If there is more than one i minimizing a_{i0}/a_{ij} the problem has a degenerate basic feasible solution (see example sheet.) For small problems you will be OK if you just choose any one of them and carry on regardless. We choose i = 2 since 3/1 < 6/1, so the variable corresponding to equation 2 leaves the basis.

3. Pivot on the element a_{ij}. (i.e., get the equations into the appropriate form for the new basic solution.)

(a) multiply row i by 1/a_{ij}.

(b) add −(a_{kj}/a_{ij})×(row i) to each row k ≠ i, including the objective function row.

We obtain as before

x_{1} x_{2} z_{1} z_{2} a_{i0}
z_{1} basic 0 3 1 −1 3
x_{1} basic 1 −1 0 1 3
a_{0j}  0 2 0 −1 −3

which is the appropriate form for the tableau for vertex B. Check that repeating these instructions on the new tableau, by pivoting on a_{12}, produces the appropriate tableau for vertex C.

x_{1} x_{2} z_{1} z_{2} a_{i0}
x_{2} basic 0 1 \frac{1}{3} −\frac{1}{3}  1
x_{1} basic 1 0 \frac{1}{3} \frac{2}{3} 4
a_{0j}  0 0 −\frac{2}{3} −\frac{1}{3} −5

Note that the columns of the tableau corresponding to the basic variables always make up the columns of an identity matrix. Since the bottom row is now all ≤ 0 the algorithm stops. Don’t hesitate to look back at Subsection 3.1 to see why we take these steps.

4 The Simplex Tableau

4.1 Choice of pivot column

We might have chosen the first pivot as a_{12} which would have resulted in

x_{1} x_{2} z_{1} z_{2} a_{i0}
z_{1} basic 1 2 1 0 6
z_{2} basic 1 −1 0 1 3
a_{0j}  1 1 0 0 0

⟶ x_{1} x_{2} z_{1} z_{2} a_{i0} x_{2} basic \frac{1}{2} 1 \frac{1}{2}  0 3 z_{2} basic \frac{3}{2} 0 \frac{1}{2}  1 6 a_{0j}  \frac{1}{2} 0 −\frac{1}{2}  0 −3

This is the tableau for vertex D. A further iteration, with pivot a_{21} takes us to the optimal solution at vertex C. Therefore both choices of initial pivot column resulted in it requiring two steps to reach the optimum. Remarks.

1. In general, there is no way to tell in advance which choice of pivot column will result in the smallest number of iterations. We may choose any column where a_{0j} > 0. A common rule-of-thumb is to choose the column for which a_{0j} is greatest, since the objective function increases by the greatest amount per unit increase in the variable corresponding to that column.

2. At each stage of the simplex algorithm we have two things in mind. First — a particular choice of basis and basic solution. Second — a rewriting of the problem in a convenient form.

There is always an identity matrix embedded in the tableau corresponding to the basic variables. Hence, when the non-basic variables are set to zero the equations are trivial to solve for the values of the basic variables. They are just given by the right-hand column. Check that provided we start with the equations written in this form in the initial tableau, the simplex algorithm rules ensure that we obtain an identity matrix at each stage.

3. The tableau obviously contains some redundant information. For example, provided we keep track of which equation corresponds to a basic variable, we could omit the columns corresponding to the identity matrix (and zeros in the objective row). This is good for computer programs, but it is probably better to keep the whole thing for hand calculation.

4.2 Initialization: the two-phase method

In our example there was an obvious basic feasible solution with which to start the simplex algorithm. This is not always the case. For example, suppose we have a problem like

maximize −6x_{1} − 3x_{2}
subject to  x_{1} + x_{2} ≥ 1
2x_{1} − x_{2} ≥ 1
3x_{2} ≤ 2
x_{1}, x_{2} ≥ 0

which we wish to solve using the simplex algorithm. We can add slack variables (sometimes called surplus variables when they appear in ≥ constraints) to obtain

maximize −6x_{1} − 3x_{2}
subject to x_{1} + x_{2} − z_{1} = 1
2x_{1} − x_{2} − z_{2} = 1
3x_{2} + z_{3} = 2
x_{i}, z_{i} ≥ 0

but there is no obvious b.f.s. since z_{1} = −1, z_{2} = −1, z_{3} = 2 is not feasible.

The trick is to add extra variables called artificial variables, y_{1}, y_{2} so that the constraints are

\begin{aligned}x_{1} + x_{2} − z_{1} + y_{1} & = & 1 \\ 2x_{1} − x_{2} − z_{2} + y_{2} & = & 1 \\ 3x_{2} + z_{3} & = & 2 \\ x_{i}, z_{i}, y_{i} & ≥ & 0 \end{aligned}

Phase I is to minimize y_{1} + y_{2} and we can start this phase with y_{1} = 1, y_{2} = 1 and z_{3} = 2. (Notice we did not need an artificial variable in the third equation.) Provided the original problem is feasible we should be able to obtain a minimum of 0 with y_{1} = y_{2} = 0 (since y_{1} and y_{2} are not needed to satisfy the constraints if the original problem is feasible). The point of this is that at the end of Phase I the simplex algorithm will have found a b.f.s. for the original problem. Phase II is then to proceed with the solution of the original problem, starting from this b.f.s.

Note: the original objective function doesn’t enter into Phase I, but it is useful to carry it along as an extra row in the tableau since the algorithm will then arrange for it to be in the appropriate form to start Phase II.

Note also: the Phase I objective must be written in terms of the non-basic variables to begin Phase I. This can also be accomplished in the tableau. We start with

x_{1} x_{2} z_{1} z_{2} z_{3} y_{1} y_{2}
y_{1} 1 1 −1 0 0 1 0 1
y_{2} 2 −1 0 −1 0 0 1 1
z_{3} 0 3 0 0 1 0 0 2
Phase II −6 −3 0 0 0 0 0 0
Phase I  0 0 0 0 0 −1 −1 0

Preliminary step. Add rows 1 and 2 to the Phase I objective so that it is written in terms of non-basic variables.

x_{1} x_{2} z_{1} z_{2} z_{3} y_{1} y_{2}
y_{1} 1 1 −1 0 0 1 0 1
y_{2} 2 −1 0 −1 0 0 1 1
z_{3} 0 3 0 0 1 0 0 2
Phase II −6 −3 0 0 0 0 0 0
Phase I  3 0 −1 −1 0 0 0 2

Begin Phase I.

Pivot on ; a_{21} to get \begin{aligned}& x_{1} & x_{2} & z_{1} & z_{2} & z_{3} & y_{1} & y_{2} \\ y_{1} & 0 & \frac{3}{2} & −1 & \frac{1}{2} & 0 & 1 & −\frac{1}{2} & \frac{1}{2} \\ x_{1} & 1 & −\frac{1}{2} & 0 & −\frac{1}{2} & 0 & 0 & \frac{1}{2} & \frac{1}{2} \\ z_{3} & 0 & 3 & 0 & 0 & 1 & 0 & 0 & 2 \\ & 0 & −6 & 0 & −3 & 0 & 0 & 3 & 3 \\ & 0 & \frac{3}{2} & −1 & \frac{1}{2} & 0 & 0 & −\frac{3}{2} & \frac{1}{2} \end{aligned}

Pivot on ; a_{14} to get \begin{aligned}& x_{1} & x_{2} & z_{1} & z_{2} & z_{3} & y_{1} & y_{2} \\ z_{2} & 0 & 3 & −2 & 1 & 0 & 2 & −1 & 1 \\ x_{1} & 1 & 1 & −1 & 0 & 0 & 1 & 0 & 1 \\ z_{3} & 0 & 3 & 0 & 0 & 1 & 0 & 0 & 2 \\ & 0 & 3 & −6 & 0 & 0 & 6 & 0 & 6 \\ & 0 & 0 & 0 & 0 & 0 & −1 & −1 & 0 \end{aligned}

End of Phase I. y_{1} = y_{2} = 0 and we no longer need these variables (so we drop the last two columns and Phase I objective row.) But we do have a b.f.s. to start Phase II with (x_{1} = 1, z_{2} = 1, z_{3} = 2) and the rest of the tableau is already in appropriate form. So we rewrite the last tableau without the y_{1}, y_{2} columns.

Begin Phase II.

x_{1} x_{2} z_{1} z_{2} z_{3}
z_{2} 0 3 −2 1 0 1
x_{1} 1 1 −1 0 0 1
z_{3} 0 3 0 0 1 2
0 3 −6 0 0 6

In one more step we reach the optimum, by pivoting on a_{12}.

x_{1} x_{2} z_{1} z_{2} z_{3}
x_{2} 0 1 −\frac{2}{3}  \frac{1}{3}  0 \frac{1}{3}
x_{1} 1 0 −\frac{1}{3} −\frac{1}{3}  0 \frac{2}{3}
z_{3} 0 0 2 −1 1 1
0 0 −4 −1 0 5

Notice that in fact, the problem we have considered is the same as the problem D, except that x replaces λ and we have added a constraint 2x_{2} ≤ 3 (which is not tight at the optimum). It is interesting to compare the final tableau with the final tableau for problem P (shown again below).

In general, artificial variables are needed when there are constraints like

≤ −1, or ≥ 1, or = 1,

unless constraints happen to be of a special form where it is easy to spot a b.f.s.

If the Phase I objective does not reach zero then the original problem is infeasible.

4.3 Sensitivity: shadow prices

Remarks.

1. Each row of each tableau merely consists of sums of multiples of rows of the original tableau. The objective row = original objective row + scalar multiples of other rows.

Consider the initial and final tableau for problem P.

initial  \begin{matrix}1 & 2 & 1 & 0 & 6 \\ 1 & −1 & 0 & 1 & 3 \\ 1 & 1 & 0 & 0 & 0 \end{matrix}  final  \begin{matrix}0 & 1 & \frac{1}{3} & −\frac{1}{3} & 1 \\ 1 & 0 & \frac{1}{3} & \frac{2}{3} & 4 \\ 0 & 0 & −\frac{2}{3} & −\frac{1}{3} & −5 \end{matrix}

In particular, look at the columns 3 and 4, corresponding to variables z_{1} and z_{2}. We can see that Final row (1) = \frac{1}{3} initial row (1) −\frac{1}{3} initial row (2) Final row (2) = \frac{1}{3} initial row (1) +\frac{2}{3} initial row (2) Final objective row = initial objective row −\frac{2}{3} initial row (1) −\frac{1}{3} initial row (2). In particular, suppose we want to make a small change in b, so we replace ( \begin{matrix}6 \\ 3 \end{matrix} ) by ( \begin{matrix}6 + ϵ_{1} \\ 3 + ϵ_{2} \end{matrix} ). Providing ϵ_{1}, ϵ_{2} are small enough they will not affect the sequence of simplex operations. Thus if the constraints move just a little the optimum will still occur with the same variables in the basis. The argument above indicates that the final tableau will be

\begin{matrix}0 & 1 & \frac{1}{3} & −\frac{1}{3} & 1 + \frac{1}{3}ϵ_{1} − \frac{1}{3}ϵ_{2} \\ 1 & 0 & \frac{1}{3} & \frac{2}{3} & 4 + \frac{1}{3}ϵ_{1} + \frac{2}{3}ϵ_{2} \\ 0 & 0 & −\frac{2}{3} & −\frac{1}{3} & −5 − \frac{2}{3}ϵ_{1} − \frac{1}{3}ϵ_{2} \end{matrix}

with corresponding solution x_{1} = 4 + \frac{1}{3}ϵ_{1} − \frac{1}{3}ϵ_{2} and x_{2} = 1 + \frac{1}{3}ϵ_{1} + \frac{2}{3}ϵ_{2} and objective function value 5 + \frac{2}{3}ϵ_{1} + \frac{1}{3}ϵ_{2}. If ϵ_{1}, ϵ_{2} are such that we have x_{1} < 0 or x_{2} < 0 then vertex C is no longer optimal. The objective function row of the final tableau shows the sensitivity of the optimal solution to changes in b and how the optimum value varies with small changes in b. For this reason the values in the objective row are sometimes known as shadow prices. Notice also that being able to see how the final tableau is related to the initial one without looking at the intermediate steps provides a useful way of checking your arithmetic if you suspect you have got something wrong!

2. Notice that for problem P the objective rows at the vertices A, B, C and D are:

\begin{matrix}A & 1 & 1 & 0 & 0 & 0 \\ B & 0 & 2 & 0 & −1 & −3 \\ C & 0 & 0 & −\frac{2}{3} & −\frac{1}{3} & −5 \\ D & \frac{1}{2} & 0 & −\frac{1}{2} & 0 & −3 \end{matrix}

Compare these values with the basic solutions of the dual problem (on page 7). You will see that the objective row of the simplex tableau corresponding to each b.f.s. of problem P is given by the dual variables for the corresponding basic solution to problem D (after a sign change). Shadow prices and dual variables are the same thing. Thus the simplex algorithm can (and should) be viewed as searching amongst primal basic feasible solutions looking for dual feasibility.

5 Lagrangian Methods

5.1 The Lagrangian

Shortly, we will begin the course afresh in a rather more theoretical style. But first, we begin with some preliminary remarks, using problems P and D as examples. Recall that so far we have been more-or-less treating the relationships between P and D as interesting coincidences. In later lectures we establish more exactly the theoretical relationships between P and D. If you feel comfortable with the simplex algorithm, the following remarks may give a clue to what is going on, as well as introducing some of the theory from the next section. Skip this bit on first reading if it does not make sense. Recall

P: maximize  x_{1} + x_{2}
subject to x_{1} + 2x_{2} + z_{1} = 6
x_{1} − x_{2} + z_{2} = 3
x_{1}, x_{2} ≥ 0

The simplex algorithm, in proceeding to the optimum, just added multiples of constraints to the objective function, and stopped if it managed to get all coefficients ≤ 0. So, it was trying to find scalar multiples λ_{1} and λ_{2} of the two constraints so that

max x_{1} + x_{2} − λ_{1}(x_{1} + 2x_{2} + z_{1} − 6) − λ_{2}(x_{1} − x_{2} + z_{2} − 3)

had a nice form. We write the Lagrangian

L(x, λ) = x_{1} + x_{2} − λ_{1}(x_{1} + 2x_{2} + z_{1} − 6) − λ_{2}(x_{1} − x_{2} + z_{2} − 3)

and rearranging terms

L(x, λ) = x_{1}(1 − λ_{1} − λ_{2}) + x_{2}(1 − 2λ_{1} + λ_{2}) − λ_{1}z_{1} − λ_{2}z_{2} + 6λ_{1} + 3λ_{2}.

The λs are called Lagrange multipliers.

More specifically, the simplex algorithm was looking for a b.f.s. x, z to P and the multipliers λ_{1}, λ_{2} such that the objective row was zero for basic variables and negative for non-basic variables.

\begin{aligned}x_{1} = 0, & ⟹ & 1 − λ_{1} − λ_{2} < 0 \\ x_{1} > 0, & ⟹ & 1 − λ_{1} − λ_{2} = 0 \end{aligned}  \begin{aligned}x_{2} = 0, & ⟹ & 1 − 2λ_{1} + λ_{2} < 0 \\ x_{2} > 0, & ⟹ & 1 − 2λ_{1} + λ_{2} = 0 \end{aligned}

and

\begin{aligned}z_{1} = 0, & ⟹ & λ_{1} > 0 \\ z_{1} > 0, & ⟹ & λ_{1} = 0 \end{aligned}  \begin{aligned}z_{2} = 0, & ⟹ & λ_{2} > 0 \\ z_{2} > 0, & ⟹ & λ_{2} = 0 \end{aligned}

But these conditions can be summarised as dual feasibility:

−v_{1} = 1 − λ_{1} − λ_{2} ≤ 0, −v_{2} = 1 − 2λ_{1} + λ_{2} ≤ 0, λ_{1}, λ_{2} ≥ 0,

and complementary slackness:

x_{1}(1 − λ_{1} − λ_{2}) = 0, x_{2}(1 − 2λ_{1} + λ_{2}) = 0, λ_{1}z_{1} = 0, λ_{2}z_{2} = 0.

If all these are satisfied, then objective = λ^{⊤}b = c^{⊤}x. See Theorem 2.4.

This is where the dual problem comes from. Lagrange multipliers are the same thing as dual variables (and shadow prices, and two more names we will encounter before the course ends.)

5.2 The Lagrangian sufficiency theorem

Given a general optimization under constraints problem,

P: min f(x) s.t. g(x) = b, x ∈ X,

with x ∈ ℝ^{n}, b ∈ ℝ^{m} (n variables and m constraints), the Lagrangian is

L(x, λ) = f(x) − λ^{⊤}(g(x) − b),

with λ ∈ ℝ^{m} (one for each constraint). The Lagrangian is a function L : ℝ^{n} × ℝ^{m} → ℝ.

(Some motivation for this definition can be found at the end of the last section, but note that we do not require anything to be linear for this definition, or the results that follow.) The following theorem will be useful, and is easy to prove.

Theorem 5.1 (Lagrangian sufficiency theorem) If x^{∗} and λ^{∗} exist such that x^{∗} is feasible for P and

L(x^{∗}, λ^{∗}) ≤ L(x, λ^{∗}) ∀x ∈ X,

then x^{∗} is optimal for P.

Proof. Define

X_{b} = {x : x ∈ X and g(x) = b}.

Note that X_{b} ⊆ X and that for any x ∈ X_{b}

L(x, λ) = f(x) − λ^{⊤}(g(x) − b) = f(x).

Now

f(x^{∗}) = L(x^{∗}, λ^{∗}) ≤ L(x, λ^{∗}) = f(x), ∀x ∈ X_{b}.

The inequality in the above holds for all x ∈ X by assumption, and so in particular for x ∈ X_{b}. Thus x^{∗} is optimal for P. Remarks.

1. Note the ‘If’ which starts the statement of the theorem. There is no guarantee that we can find a λ^{∗} satisfying the conditions of the theorem for general problems P. (However, we will see that there is a large class of problems for which λ^{∗} do exist.)

2. At first sight the theorem offers us a method for testing that a solution x^{∗} is optimal for P without helping us to find x^{∗} if we don’t already know it. Certainly we will sometimes use the theorem this way. But for some problems, there is a way we can use the theorem to find the optimal x^{∗}.

5.3 Strategy to solve problems with the Lagrangian sufficiency theorem

Attempt to find x^{∗}, λ^{∗} satisfying the conditions of the theorem as follows.

1. For each λ solve the problem

minimize L(x, λ) subject to x ∈ X.

Note that the only constraints involved in this problem are x ∈ X so this should be an easier problem to solve than P.

2. Define the set

Y = {λ : min_{x∈X} L(x, λ) > −∞}.

If we obtain −∞ for the minimum in step 1 then that λ is no good. We consider only those λ ∈ Y for which we obtain a finite minimum.

3. For λ ∈ Y, the minimum will be obtained at some x(λ) (that depends on λ in general). Typically, x(λ) will not be feasible for P.

4. Adjust λ ∈ Y so that x(λ) is feasible. If λ^{∗} ∈ Y exists such that x^{∗} = x(λ^{∗}) is feasible then x^{∗} is optimal for P by the theorem.

5.4 Example: use of the Lagrangian sufficiency theorem

minimize x_{1} − x_{2} − 2x_{3}
s.t. x_{1} + x_{2} + x_{3} = 5
x_{1}^{2} + x_{2}^{2} = 4
x ∈ X = ℝ^{3}.

Solution. Since we have two constraints we take λ ∈ ℝ^{2} and write the Lagrangian

L(x, λ) = x_{1} − x_{2} − 2x_{3} − λ_{1}(x_{1} + x_{2} + x_{3} − 5) − λ_{2}(x_{1}^{2} + x_{2}^{2} − 4)
= [x_{1}(1 − λ_{1}) − λ_{2}x_{1}^{2}] + [x_{2}(−1 − λ_{1}) − λ_{2}x_{2}^{2}]
+[x_{3}(−2 − λ_{1})] + 5λ_{1} + 4λ_{2}.

We first try to minimize L(x, λ) for fixed λ in x ∈ ℝ^{3}. Notice that we can minimize each square bracket separately.

First notice that −x_{3}(2 + λ_{1}) has minimum −∞ unless λ_{1} = −2. So we only want to consider λ_{1} = −2.

Now observe that the terms in x_{1}, x_{2} have a finite minimum only if λ_{2} < 0 in which case the minimum occurs at

\begin{aligned}1 − λ_{1} − 2λ_{2}x_{1} = 0 & ⟹ & x_{1} = 3/2λ_{2} \\ −1 − λ_{1} − 2λ_{2}x_{2} = 0 & ⟹ & x_{2} = 1/2λ_{2}. \end{aligned}

So Y is the set

Y = {λ : λ_{1} = −2, λ_{2} < 0},

and for λ ∈ Y the minimum of L(x, λ) occurs at x(λ) = (3/2λ_{2}, 1/2λ_{2}, x_{3})^{⊤}. Now to find a feasible x(λ) we need

x_{1}^{2} + x_{2}^{2} = 4 ⟹ \frac{9}{4λ_{2}^{2}} + \frac{1}{4λ_{2}^{2}} = 4 ⟹ λ_{2} = −\sqrt{5/8}.

So x_{1} = −3\sqrt{2/5}, x_{2} = −\sqrt{2/5} and x_{3} = 5 − x_{1} − x_{2} = 5 + 4\sqrt{2/5}.

The conditions of the Lagrangian sufficiency theorem are satisfied by

x^{∗} = (−3\sqrt{2/5}, −\sqrt{2/5}, 5 + 4\sqrt{2/5})^{⊤} and λ^{∗} = (−2, −\sqrt{5/8})^{⊤} .

So x^{∗} is optimal.

6 The Lagrangian Dual

6.1 Example: further use of the Lagrangian sufficiency theorem

Consider the problem

minimize \frac{1}{1 + x_{1}} + \frac{1}{2 + x_{2}}

s.t. x_{1} + x_{2} = b, x_{1}, x_{2} ≥ 0.

Solution. We define X = {x : x ≥ 0} and the Lagrangian

L(x, λ) = \frac{1}{1 + x_{1}} + \frac{1}{2 + x_{2}} − λ(x_{1} + x_{2} − b)
= (\frac{1}{1 + x_{1}} − λx_{1}) + (\frac{1}{2 + x_{2}} − λx_{2}) + λb.

Note that

(\frac{1}{1 + x_{1}} − λx_{1}) and (\frac{1}{2 + x_{2}} − λx_{2})

do not have a finite minimum in x ≥ 0 unless λ ≤ 0. So we take λ ≤ 0. Observe that in the range x ≥ 0 a function of the form (\frac{1}{c+x} − λx) will have its minimum either at x = 0, if this function is increasing at 0, or at the stationary point of the function, occurring where x > 0, if the function is decreasing at 0. So the minimum occurs at

x = \begin{matrix}0 \\ −c + \sqrt{−1/λ} \end{matrix} as \sqrt{−1/λ} \begin{matrix}≤ \\ ≥ \end{matrix} c.

Thus for \sqrt{−1/λ} ≤ 1 we have a minimum at x_{1} = 0 and otherwise the minimum is at x_{1} = −1 + \sqrt{−1/λ}. Similarly, if \sqrt{−1/λ} ≤ 2 we have a minimum at x_{2} = 0 and otherwise the minimum is at x_{2} = −2 + \sqrt{−1/λ}. So the set Y is \{λ : λ ≤ 0\} and for λ ∈ Y, the minimum of L(x, λ) in X is

x_{1} = 0,  x_{2} = 0  ≤ 1
x_{1} = −1 + \sqrt{−1/λ}, x_{2} = 0  as \sqrt{−1/λ} ∈ [1, 2]
x_{1} = −1 + \sqrt{−1/λ}, x_{2} = −2 + \sqrt{−1/λ}  ≥ 2

At first sight it appears that we don’t know which values of x_{1}, x_{2} to substitute into the constraint until we know λ, and we don’t know λ until we substitute x_{1}, x_{2} into the constraint. But notice that x_{1}(λ) + x_{2}(λ) satisfies

x_{1}(λ) + x_{2}(λ) = \{ \begin{matrix}0 & & ≤ −1 \\ −1 + 1/\sqrt{−λ} & as λ & ∈ [−1, −1/4] \\ −3 + 2/\sqrt{−λ} & & ∈ [−1/4, 0] \end{matrix}

PSfrag replacements

λ

x_{1}(λ) + x_{2}(λ)

−1 + \frac{1}{\sqrt{−λ}}

−3 + \frac{2}{\sqrt{−λ}}

−1 −\frac{1}{4}

1

So we can see that x_{1}(λ) + x_{2}(λ) is an increasing and continuous function (although it is not differentiable at λ = −1 and at λ = −1/4). Thus (by the Intermediate Value Theorem) for any b > 0 there will be a unique value of λ, say λ^{∗}, for which x_{1}(λ^{∗}) + x_{2}(λ^{∗}) = b. This λ^{∗} and corresponding x^{∗} will satisfy the conditions of the theorem and so x^{∗} is optimal.

Examples of this kind are fairly common.

6.2 The Lagrangian dual problem

We have already defined

Y = {λ : min_{x∈X} L(x, λ) > −∞}.

Now, for λ ∈ Y define

L(λ) = min_{x∈X} L(x, λ).

The following lemma is almost as easy to prove as the sufficiency theorem.

Lemma 6.1 (weak duality) For any feasible x ∈ X_{b} and any λ ∈ Y

L(λ) ≤ f(x).

Proof. For x ∈ X_{b}, λ ∈ Y,

f(x) = L(x, λ) ≥ min_{x∈X_{b}} L(x, λ) ≥ min_{x∈X} L(x, λ) = L(λ)

Thus, provided the set Y is non-empty, L(λ), λ ∈ Y provides a lower bound for the value of the objective function f(x) at the optimum for P.

It seems reasonable to try and get this lower bound as good as possible, i.e., to consider the problem

D: maximize L(λ) subject to λ ∈ Y, equivalently,

D: maximize_{λ∈Y} {min_{x∈X} L(x, λ)} .

This is known as the dual problem. (The original is the primal problem.)

In a moment we will return to linear problems and justify some of what we did in the first part of the course. For now, notice that the idea of a dual problem is quite general. For example, we can look again at the two examples we just studied.

In the example of Section 5.4 we had Y = {λ : λ_{1} = −2, λ_{2} < 0} and that min_{x∈X} L(x, λ) occurred for x(λ) = (3/2λ_{2}, 1/2λ_{2}, x_{3}). Thus

L(λ) = L(x(λ), λ) = \frac{10}{4λ_{2}} − 10 + 4λ_{2}.

The dual problem is thus

maximize_{λ_{2}<0} \{\frac{10}{4λ_{2}} − 10 + 4λ_{2}\} .

Note that the max occurs at λ_{2} = −\sqrt{5/8} and that the primal and dual have the same optimal value.

In the the example of Section 6.1, Y = {λ : λ ≤ 0}. By substituting the optimal value of x into L(x, λ) we obtain

L(λ) = \begin{aligned}3/2 + λb & & ≤ −1 \\ 1/2 + 2\sqrt{−λ} + (b + 1)λ & as λ & ∈ [−1, −1/4] \\ 4\sqrt{−λ} + (b + 3)λ & & ∈ [−1/4, 0] \end{aligned}

We can solve the dual problem max L(λ) s.t. λ ≤ 0. The solution lies in −1 ≤ λ ≤ −1/4 if 0 ≤ b ≤ 1 and in −1/4 ≤ λ ≤ 0 if 1 ≤ b. You can confirm for all b that primal and dual here have the same optimal values.

6.3 The dual problem for LP

Construction of the dual problem for LP is relatively straightforward. Consider the primal problem P:

maximize c^{⊤}x
subject to Ax ≤ b, x ≥ 0
equivalently Ax + z = b, x, z ≥ 0.

Write the Lagrangian

L(x, z, λ) = c^{⊤}x − λ^{⊤}(Ax + z − b) = (c^{⊤} − λ^{⊤}A)x − λ^{⊤}z + λ^{⊤}b.

As in the general case, we can find the set Y such that λ ∈ Y implies max_{x,z≥0} L(x, z, λ) is finite, and for λ ∈ Y we compute the minimum of L(λ).

Consider the linear term −λ^{⊤}z. If any coordinate λ_{i} < 0 we can make −λ_{i}z_{i} as large as we like, by taking z_{i} large. So there is only a finite maximum in z ≥ 0 if λ_{i} ≥ 0 for all i.

Similarly, considering the term (c^{⊤} − λ^{⊤}A)x, this can be made as large as we like unless (c^{⊤} − λ^{⊤}A)_{i} ≤ 0 for all i. Thus

Y = {λ : λ ≥ 0, λ^{⊤}A − c^{⊤} ≥ 0}.

If we pick a λ ∈ Y then max_{z≥0} −λ^{⊤}z = 0 (by choosing z_{i} = 0 if λ_{i} > 0 and any z_{i} if λ_{i} = 0) and also max_{x≥0}(c^{⊤} − λ^{⊤}A)x = 0 similarly. Thus for λ ∈ Y, L(λ) = λ^{⊤}b. So the dual problem for P is, as we have defined previously,

D: min λ^{⊤}b s.t. A^{⊤}λ ≥ c, λ ≥ 0.

Remark. Notice that in this section we have taken the primal to be a maximization problem and the dual to be a minimization problem. In Sections 5.2 and 6.2 these were the other way around. In general, the dual problem is always obtained by extremizing L(λ) with respect to λ ∈ Y in the opposite direction to that in which L(x, λ) is extremized with respect to x ∈ X.

6.4 The weak duality theorem in the case of LP

We can now apply Lemma 6.1 directly to P and D to obtain one of the theorems which was first stated in Lecture 2. Theorem 2.3 (weak duality) If x is feasible for P (so Ax ≤ b, x ≥ 0) and λ is feasible for D (so λ ≥ 0, A^{⊤}λ ≥ c) then c^{⊤}x ≤ λ^{⊤}b. Since this is an important (though easy) result it may be worth writing out a proof for this particular case which does not appeal to the general Lemma 6.1. Naturally enough, the proof is very similar. Proof. Write

L(x, z, λ) = c^{⊤}x − λ^{⊤}(Ax + z − b)

where Ax + z = b. Now for x and λ satisfying the conditions of the theorem,

c^{⊤}x = L(x, z, λ) = (c^{⊤} − λ^{⊤}A)x − λ^{⊤}z + λ^{⊤}b ≤ λ^{⊤}b.

7 Shadow Prices and Lagrangian Necessity

7.1 Sufficient conditions for optimality

Theorem 2.3 can be used to give a quick proof of the sufficient conditions for optimality. Theorem 2.4 (sufficient conditions for optimality) If x^{∗}, z^{∗} is feasible for P and λ^{∗} is feasible for D and (c^{⊤} − λ^{∗⊤}A)x^{∗} = λ^{∗⊤}z^{∗} = 0 (complementary slackness) then x^{∗} is optimal for P and λ^{∗} is optimal for D. Furthermore c^{⊤}x^{∗} = λ^{∗⊤}b. Proof. Write L(x^{∗}, z^{∗}, λ^{∗}) = c^{⊤}x^{∗} − λ^{∗⊤}(Ax^{∗} + z^{∗} − b). Now

c^{⊤}x^{∗} = L(x^{∗}, z^{∗}, λ^{∗})
= (c^{⊤} − λ^{∗⊤}A)x^{∗} − λ^{∗⊤}z^{∗} + λ^{∗⊤}b
= λ^{∗⊤}b

But for all x feasible for P we have c^{⊤}x ≤ λ^{∗⊤}b (by the weak duality theorem) and this implies that for all feasible x, c^{⊤}x ≤ c^{⊤}x^{∗}. So x^{∗} is optimal for P. Similarly, λ^{∗} is optimal for D (and the problems have the same optimums).

Remarks.

1. As an alternative to the above proof we could have just applied the Lagrangian sufficiency theorem directly to P using the λ^{∗} given in the theorem. The proof comes out very similar.

2. The theorems above explain most of what we observed in Lecture 2. There we were looking at two problems, P and D, that were dual to one another. We paired basic solutions to P and D by insisting on complementary slackness. Then since our solution to P satisfied Ax + z = b the same argument as in the proof above shows that we will obtain identical values of the objective function. (We don’t use dual feasibility in the first line of the proof.) The theorem on sufficient conditions for optimality explains why it was that we only obtained a pair of solutions feasible for P and D at the optimum of both.

7.2 Shadow prices

We have mentioned before that dual variables or Lagrangian multipliers are sometimes known as shadow prices. First, an explanation of shadow prices.

Suppose you require an amount b_{1}, …, b_{m} of m different vitamins. There are n foodstuffs available. Let a_{ij} =amounts of vitamin i in one unit of foodstuff j, and suppose foodstuff j costs c_{j} per unit. Your problem is therefore to choose the amount x_{j} of foodstuff j you buy to solve the LP

min ∑_{j} c_{j}x_{j}

subject to ∑_{j} a_{ij}x_{j} ≥ b_{i}, each i

and x_{j} ≥ 0 each j.

Now suppose that a vitamin company decides to market m different vitamin pills (one for each vitamin) and sell them at price p_{i} per unit for vitamin i. Assuming you are prepared to switch entirely to a diet of vitamin pills, but that you are not prepared to pay more for an artificial carrot (vitamin equivalent) than a real one, the company has to maximize profit by choosing prices p_{i} to

max ∑_{i} b_{i}p_{i}

subject to ∑_{i} a_{ij}p_{i} ≤ c_{j}, each j

and p_{i} ≥ 0 each i.

Note that this is the LP which is dual to your problem. The dual variable p_{i} is the price you are prepared to pay for a unit of vitamin i and is called a shadow price. By extension, dual variables are sometimes called shadow prices in problems where their interpretation as prices is very hard (or impossible) to see.

The dual variables tell us how the optimum value of our problem changes with changes in the right-hand side (b) of our functional constraints. This makes sense in the example given above. If you require an amount b_{i} + ϵ of vitamin i instead of an amount b_{i} you would expect the total cost of your foodstuff to change by an amount ϵp_{i}, where p_{i} is the value to you of a unit of vitamin i, even though in your problem you cannot buy vitamin i separately from the others.

Thus in the example in Section 5.4 we had constraints

x_{1} + x_{2} + x_{3} = 5
x_{1}^{2} + x_{2}^{2} = 4

and obtained values of Lagrange multipliers of λ_{1}^{∗} = −2, λ_{2}^{∗} = −\sqrt{5/8}. If we replace the constraints by

x_{1} + x_{2} + x_{3} = b_{1}
x_{1}^{2} + x_{2}^{2} = b_{2}

and write ϕ(b) = optimal value of the problem with b = (b_{1}, b_{2})^{⊤} then you can check that

\frac{∂ϕ}{∂b_{1}} |_{b=(5,4)} = −2 \frac{∂ϕ}{∂b_{2}} |_{b=(5,4)} = −\sqrt{5/8}.

Let us explain a little more why this is true.

In linear problems we can use what we know about the optimal solutions to see how this works. Let us assume the primal problem

P (b) : min c^{⊤}x
s.t. Ax + z = b, x ≥ 0.

has the optimal solution ϕ(b), depending on b. Consider two close together values of b, b′ and b′′ say, and suppose that optimal solutions have the same basic variables (so optimums occur with the same variables non-zero though the values of the variables will change slightly). The optimum still occurs at the same vertex of the feasible region though it moves slightly. Now consider the dual problem. This is

max λ^{⊤}b
s.t. A^{⊤}λ ≤ c, λ ≤ 0.

In the dual problem the feasible set does not depend on b, so the optimum of the dual will occur with the same basic variables and the same values of the dual variables λ. But the value of the optimum dual objective function is λ^{⊤}b′ in one case and λ^{⊤}b′′ in the other and we have seen that the primal and dual have the same solutions. Hence

ϕ(b′) = λ^{⊤}b′ and ϕ(b′′) = λ^{⊤}b′′

and the values of the dual variables λ give the rate of change of the objective value with b. The change is linear in this case. Remarks.

1. You can also test this argument on the non-linear problems we have studied.

2. Recall that we made an argument that the objective row of the final tableau contains constants which give the rate of change of the optimum with b (based on the fact that the final tableau is the initial tableau + multiples of the original constraints). But we also observed that the objective row contains values of the dual variables at the optimum. That argument and the argument above are two versions of the same thing.

7.3 Lagrangian necessity

In all the examples we have studied we have been able to find λs such that we can apply the Lagrangian Sufficiency Theorem. Furthermore, we have observed that the primal and dual problems have the same optimum values. It is outside the scope of the course to establish conditions under which we expect these results to hold, but it seems helpful to give a brief summary.

Theorem 7.1 If the functions f, g are convex, X is convex and x^{∗} is an optimal solution to P, then there exist Lagrange multipliers λ ∈ ℝ^{m} such that L(x^{∗}, λ) ≤ L(x, λ) ∀x ∈ X.

In particular, Lagrange multipliers will always exist for LP problems provided the problem is bounded (and primal and dual problems will have the same optimal value).

Define ϕ(b) as min f(x), subject to x ∈ X and g(x) = b. The proof of the theorem proceeds by first showing that ϕ(b) is convex. This implies that for each b^{∗} there is a tangent hyperplane to ϕ(b) at b^{∗} with the graph of ϕ(b) lying entirely above it. This uses the ‘supporting hyperplane theorem’ for convex sets, which is geometrically obvious, but some work to prove. Note the equation of the hyperplane will be y = ϕ(b^{∗}) + λ^{⊤}(b − b^{∗}) for some multipliers λ. This λ can be shown to be the required vector of Lagrangian multipliers and the picture below gives some intuition of a geometric kind into why the Lagrange multipliers λ exist and why these λs give the rate of change of the optimum ϕ(b) with b.

PSfrag replacements

λ min L

L = f − λ^{⊤}(g − b^{∗}) f g b^{∗} ϕ(b) ϕ(b^{∗})

feasible (g, f) ϕ(b^{∗}) + λ^{⊤}(g − b^{∗})

λ^{⊤}(g − b^{∗})

7.4 Convexity of the feasible set in the case of LP

In Lecture 1 we advertised three general theorems for LPs. We should now prove these. We write the LP as min c^{⊤}x s.t. Ax = b, x ≥ 0.

Theorem 1.1 The feasible set X_{b} is convex. Proof. Suppose x ∈ X_{b} and y ∈ X_{b}. So x, y ≥ 0 and Ax = Ay = b. Consider z = λx + (1 − λ)y for 0 ≤ λ ≤ 1. Then z_{i} = λx_{i} + (1 − λ)y_{i} ≥ 0 for all i. So z ≥ 0. Secondly,

Az = A(λx + (1 − λ)y) = λAx + (1 − λ)Ay = λb + (1 − λ)b = b.

So z ∈ X_{b} and hence X_{b} is convex.

8 Algebra of Linear Programming

8.1 Basic feasible solutions and extreme points

To prove the equivalence of basic feasible solutions and extreme points which was stated in Theorem 2.1 we need to recall some results on the solution of simultaneous linear equations.

1. If A is a m × m invertible matrix, b ∈ ℝ^{m}, then Ax = b has a unique solution, x = A^{−1}b.

2. If A is a m×n with n > m and no row of A is a linear combination of other rows (so there are no redundant equations or inconsistent equations, i.e., rank(A) = m) then there are (a) infinitely many solutions to Ax = b, and (b) infinitely many solutions to Ay = 0.

We will also find it convenient to divide the variables into two sets, and to split the matrix A accordingly. For example, given

( \begin{matrix}a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \end{matrix} ) ( \begin{matrix}x_{1} \\ x_{2} \\ x_{3} \end{matrix} ) = ( \begin{matrix}b_{1} \\ b_{2} \end{matrix} )

we can partition the variables into two disjoint subsets (basic and non-basic) B = {1, 2} and N = {3} and rewrite the equation

( \begin{matrix}a_{11} & a_{12} \\ a_{21} & a_{22} \end{matrix} ) ( \begin{matrix}x_{1} \\ x_{2} \end{matrix} ) + ( \begin{matrix}a_{13} \\ a_{23} \end{matrix} ) x_{3} = ( \begin{matrix}b_{1} \\ b_{2} \end{matrix} )

or

A_{B}x_{B} + A_{N}x_{N} = b,

where x_{B} = ( \begin{matrix}x_{1} \\ x_{2} \end{matrix} ) contains the variables in B and x_{N} = ( x_{3} ) contains the vari-

ables in N, and A_{B} has the columns of A corresponding to variables in B (columns 1 and 2) and A_{N} has columns corresponding to variables from N (column 3).

You should convince yourself that the two forms of the linear equations are equivalent and that the same trick would work for general m × n matrices A and partition of variables into two sets. If A = (a_{1}, …, a_{n}), where a_{i} is the ith column, then

Ax = ∑_{i∈B} x_{i}a_{i} + ∑_{i∈N} x_{i}a_{i} = A_{B}x_{B} + A_{N}x_{N} = b.

We usually choose B to have m elements (a basis) and N to have n − m elements. Then setting x_{N} = 0, we solve A_{B}x_{B} = b where A_{B} is an m × m matrix to find the basic solution x_{B} = A_{B}^{−1}b, x_{N} = 0. We also make the following assumptions

Non-degeneracy assumptions for A.

(a) The matrix A has linearly independent rows (i.e., rank(A) = m), and

(b) Any m × m matrix formed from m columns of A is non-singular (i.e., any set of m columns from A is linearly independent). Non-degeneracy of basic solutions.

(c) All basic solutions to x_{B} = A_{B}^{−1}b, x_{N} = 0 have exactly m non-zero coordinates, i.e., x_{i} ≠ 0 for i ∈ B.

Remarks.

(a) is a reasonable assumption for small problems such as we can do by hand, though it may be difficult to check if a large set of equations is inconsistent or contains redundancy. Assumptions (b) and (c) are just made to simplify the theory and there is no trouble in solving LPs when A does not satisfy these. For (b) there may be no solution corresponding to some choice of basis (e.g., the equations obtained by adding slack variables to x_{1} ≤ 2, x_{2} ≤ 3), but the simplex algorithm still works perfectly well. Or there may be many (e.g., x_{1} + x_{2} + z_{1} = 1, 2x_{1} + 2x_{2} + z_{2} = 2) and here there are solutions such as x_{1} = x_{2} = 1/2, z_{1} = z_{2} = 0 which are not extreme points. The answer here is to alter the definition of basic solution to exclude those obtained when A_{B} is not invertible. You will see on the examples sheet that the failure of (c) to hold causes no trouble. The proof which follows can be adapted to work if A does not satisfy the nondegeneracy assumptions as stated above.

Lemma 8.1 If x ∈ X_{b} and x = λy + (1 − λ)z for y, z ∈ X_{b}, y ≠ z, 0 < λ < 1, (i.e., x is not extreme) then x_{i} = 0 implies y_{i} = z_{i} = 0.

Proof. x_{i} = λy_{i} + (1 − λ)z_{i}. Suppose x_{i} = 0. Then since λ, 1 − λ, y_{i}, z_{i} are all ≥ 0, we must have y_{i} = z_{i} = 0. This result says that if x is on some face or edge of the simplex X_{b} then y and z must be in the same face or edge.

Lemma 8.2 If x is a basic feasible solution then x is an extreme point of X_{b}.

Proof. If x is a b.f.s. then ∃ a basis B and non-basis N such that the non-basic components of x satisfy x_{N} = 0. Suppose x is not extreme. Then ∃ y and z such that x lies on the line segment between y and z and hence y_{N} = z_{N} = 0. Furthermore, y and z are feasible, so A_{B}y_{B} = A_{B}z_{B} = b and this gives A_{B}(y_{B} − z_{B}) = 0. But as A_{B} is non-singular (by assumption (b)) the unique solution to A_{B}v_{B} = 0 is v_{B} = 0. But this implies y_{B} = z_{B}, and hence y = z. But this is a contradiction. So x is extreme.

Lemma 8.3 If x is an extreme point of X_{b} then x is a basic feasible solution.

Proof. Since x ∈ X_{b} we know it is feasible. We just need to show it is basic. Again, we argue by contradiction.

Suppose x is not a b.f.s. Then the number of non-zero coordinates, say p, is greater than m. (If the number is equal to m the solution is basic, and we can’t have less than m by the non-degeneracy assumption.) Intuitively, x is not at a vertex, but is in some face, or edge specified by setting n − p variables to zero.

Let P = {i : x_{i} > 0} and Q = {i : x_{i} = 0}. Since x is feasible, Ax = A_{P} x_{P} + A_{Q}x_{Q} = b ⟹ A_{P} x_{P} = b. But this is m equations in p > m variables, so ∃ non-zero

y_{P} s.t. A_{P} y_{P} = 0. We put y_{Q} = 0 and y = ( \begin{matrix}y_{P} \\ y_{Q} \end{matrix} ).

Consider x + ϵy and x − ϵy for small ϵ. For ϵ > 0 small enough these two points are both feasible since A(x ± ϵy) = Ax ± ϵAy = b, and x ± ϵy ≥ 0 (for small ϵ since y_{i} > 0 implies x_{i} > 0). Hence

x = \frac{1}{2}(x + ϵy) + \frac{1}{2}(x − ϵy),

and so x is not extreme.

Corollary 8.4 x is an extreme point of X_{b} if and only if x is a basic feasible solution of Ax = b. (This is a restatement of Theorem 2.1.)

We now need to prove that finite optimums occur at extreme points. The proof of this result is very like the proof of Lemma 8.3. The intuition is that if x is optimal and lies on a face (but not a vertex) then the whole face is optimal. So we construct a line segment through x which lies in the face (as in Lemma 8.3) and all solutions on this line segment are optimal. We extend the line segment till we come to an edge (or a face of one less dimension) and now have an optimal solution with one more zero coordinate. Continue until you end at an optimal vector.

Lemma 8.5 Suppose minimize c^{⊤}x s.t. x ∈ X_{b} has a finite optimal solution. Then there exists an extreme point which is optimal. (This gives Theorem 2.2.)

Proof. Suppose x is an optimal solution with the smallest number of non-zero components, but is not extreme. Then ∃ y, z such that x lies on the line segment between y and z. Since x is optimal, c^{⊤}x ≤ c^{⊤}y and c^{⊤}x ≤ c^{⊤}z. But c^{⊤}x = λc^{⊤}y + (1 − λ)c^{⊤}z and hence c^{⊤}x = c^{⊤}y = c^{⊤}z. Thus y and z are also optimal.

Now consider the line segment x(ϵ) = x + ϵ(y − z) for small ϵ. Note that c^{⊤}x(ϵ) = c^{⊤}x. So all points on the line are optimal. Furthermore, x_{i} = 0 ⟹ x(ϵ)_{i} = 0 for all ϵ and x_{i} > 0 ⟹ x(ϵ)_{i} > 0 for small ϵ. It is clear that by choosing ϵ large enough and of the right sign we can force at least one more component of x(ϵ) to become zero.

Now since x is not extreme, the number of non-zero components is p > m. But we have now constructed x(ϵ) which is optimal, feasible and has p − 1 non-zero components. This is a contradiction. So an optimal solution with the smallest possible number of non-zero components must be extreme.

8.2 Algebra of the simplex method

Let us take problem P in the form

P: max c^{⊤}x s.t. Ax = b, x ≥ 0.

Given a choice of basis B we can rewrite the problem as above

max {c_{B}^{⊤}x_{B} + c_{N}^{⊤}x_{N}} s.t. A_{B}x_{B} + A_{N}x_{N} = b, x_{B}, x_{N} ≥ 0.

At this stage we are just rewriting the equations in terms of the two sets of variables x_{B} and x_{N}. The equations hold for any feasible x. Now A_{B} is invertible by our non-degeneracy assumptions. Thus we can write

x_{B} = A_{B}^{−1}b − A_{B}^{−1}A_{N}x_{N},  (1)

and

f = c_{B}^{⊤}x_{B} + c_{N}^{⊤}x_{N}
= c_{B}^{⊤}(A_{B}^{−1}b − A_{B}^{−1}A_{N}x_{N}) + c_{N}^{⊤}x_{N}
= c_{B}^{⊤}A_{B}^{−1}b + (c_{N}^{⊤} − c_{B}^{⊤}A_{B}^{−1}A_{N})x_{N}  (2)

Equation (1) gives the constraints in a form where x_{B} appears with an identity matrix and (2) gives the objective function in terms of the non-basic variables. Thus the tableau corresponding to basis B will contain (after appropriate rearrangement of rows and columns)

basic  non-basic

I  A_{B}^{−1}A_{N}  A_{B}^{−1}b

0  c_{N}^{⊤} − c_{B}^{⊤}A_{B}^{−1}A_{N} −c_{B}^{⊤}A_{B}^{−1}b

Thus, given a basis (choice of variables) we can work out the appropriate tableau by inverting A_{B}. Note that for many choices of B we will find that A_{B}^{−1}b has negative components and so the basic solution x_{B} = A_{B}^{−1}b is not feasible; we need the simplex algorithm to tell which Bs to look at.

In programming a computer to implement the simplex algorithm you need only remember what your current basis is, since the whole tableau can then be computed from A_{B}^{−1}. The revised simplex algorithm works this way and employs various tricks to compute the inverse of the new A_{B} from the old A_{B}^{−1} (using the fact that only one variable enters and one leaves). This can be very much more efficient.

Now we know that the simplex algorithm terminates at an optimal feasible basis when the coefficients of the objective row are all ≤ 0. In other words, there is a basis B for which

c_{N}^{⊤} − c_{B}^{⊤}A_{B}^{−1}A_{N} ≤ 0.

Recall the dual problem is

D: min λ^{⊤}b s.t. A^{⊤}λ ≥ c.

Let us write λ = (A_{B}^{−1})^{⊤}c_{B}. Then we have A_{B}^{⊤}λ = c_{B} and A_{N}^{⊤}λ ≥ c_{N}, and hence λ is feasible. Furthermore, x_{B} = A_{B}^{−1}b, x_{N} = 0 is a basic solution for P and complementary slackness is satisfied since

c_{B} − A_{B}^{⊤}λ = 0 ⟹ (c_{B} − A_{B}^{⊤}λ)^{⊤}x_{B} = 0,
x_{N} = 0 ⟹ (c_{N} − A_{N}^{⊤}λ)^{⊤}x_{N} = 0.

Consequently, x_{B} = A_{B}^{−1}b, x_{N} = 0 and λ = (A^{−1})^{⊤}c_{B} are respectively optimal for the primal and dual. We also have that with these solutions

f = c_{B}^{⊤}x_{B} = c_{B}^{⊤}A_{B}^{−1}b = λ^{⊤}b.

So we have a proof of Theorem 2.5, that the primal and dual have the same objective value (if we accept that the simplex algorithm terminates at an optimum with ≤ 0 objective row) for the case of LP problems.

Remark We have shown that in general the objective row of the final (optimal) tableau will contain c_{N}^{⊤}−λ^{⊤}A_{N}, where λ are dual variables. This is consistent with our observation that the final tableau contains −λ when we start with a primal Ax ≤ b and add slack variables z. In this case the c_{i} corresponding to the slack variables are zero and columns of A_{N} corresponding to original variables are the columns of an identity matrix. So the part of the objective row beneath the original slack variables will contain −λ^{⊤}, and λ are the dual variables corresponding to the primal constraints. The rest of the objective row, beneath the original x, contains c_{B}^{⊤} − λ^{⊤}A_{B}, i.e., the slack variables for the dual problem. This is what we observed in our earlier example.

9 Two Person Zero-Sum Games

9.1 Games with a saddle-point

We consider games that are zero-sum, in the sense that one player wins what the other loses. The players make moves simultaneously. Each has a choice of moves (not necessarily the same). There is a pay-off matrix A = (a_{ij}).

II plays j
1 2 3 4
1 −5 3 1 20
I plays i 2 5 5 4 6 ← (a_{ij})
3 −4 6 0 −5

If player I makes move i and player II makes move j then player I wins (and player II loses) a_{ij}. What is the best strategy for the two players to adopt? We assume that they both know the matrix.

Let us ask what is the best that player I can do if player II plays move j.

II’s move: j =  1 2 3 4
I’s best response: i = 2 3 2 1
I wins  5 6 4 20 ← column maximums

Similarly, we ask what is the best that player II can do if I plays move i?

I’s move: i =  1 2 3
II’s best response: j = 1 3 4
I wins  −5 4 −5 ← row minimums

Here the minimal column maximum = maximal row minimum = 4, when player I plays 2 and player II plays 3. In this case we say that A has a saddle-point (2, 3) and the game is solved. Remarks. The game is solved in the sense that

1. Each player maximizes his minimum gain.

2. If either player announces any strategy (in advance) other than ‘I plays 2’ and ‘II plays 3’, he will do worse (including the announcement of a mixed strategy — when he plays 1, …, n with probabilities p_{1}, …, p_{n}).

3. Even if either player announces that he will play the saddle-point move is advance, the other player cannot improve on the saddle-point.

9.2 Example: Two-finger Morra, a game without a saddle-point

Each player displays either one or two fingers and simultaneously guesses how many fingers his opponent will show. If both players guess correctly or both guess incorrectly then the game is a tie. If only one player guesses correctly, then that player wins from the other player an amount equal to the total number of fingers shown. A strategy for a player is (a, b) =‘show a, guess b’. The pay-off matrix is

\begin{aligned}& (1,1) & (1,2) & (2,1) & (2,2) \\ (1,1) & 0 & 2 & −3 & 0 \\ (1,2) & −2 & 0 & 0 & 3 \\ (2,1) & 3 & 0 & 0 & −4 \\ (2,2) & 0 & −3 & 4 & 0 \end{aligned} = (a_{ij})

Column maximums are all positive and row minimums are all negative. So there is no saddle point (even though the game is symmetric and fair). If either player announces a fixed strategy (in advance), the other player will win. We assume each player will choose a mixed strategy where he plays move i with probability p_{i} and that each player attempts to maximize his minimum gain.

9.3 Determination of an optimal strategy

Each player must use a mixed strategy. We assume player I moves i with probability p_{i}, i = 1, …, n and player II moves j with probability q_{j}, j = 1, …, m. Player I’s expected payoff if player II plays move j is

∑_{i} a_{ij}p_{i}.

So player I attempts to

maximize {min_{j} ∑_{i} a_{ij}p_{i}} s.t. ∑_{i} p_{i} = 1, p_{i} ≥ 0.

Note that this is equivalent to

P: max v s.t. ∑_{i} a_{ij}p_{i} ≥ v, each j, and ∑_{i} p_{i} = 1, p_{i} ≥ 0,

since v on being maximized will increase until it equals the minimum of the ∑_{i} a_{ij}p_{i}. By similar arguments, player II’s problem is

D: min v s.t. ∑_{j} a_{ij}q_{j} ≤ v, each i, and ∑_{j} q_{j} = 1, q_{j} ≥ 0.

It is possible to show that P and D are duals to one another (by the standard technique of finding the dual of P). Consequently, the general theory gives sufficient conditions for strategies p and q to be optimal.

Theorem 9.1 If v, p and q, with ∑_{i} p_{i} = ∑_{j} q_{j} = 1, p_{i}, q_{j} ≥ 0 are such that p^{⊤}A ≥ v and Aq ≤ v and v = p^{⊤}Aq then p is optimal for P and q is optimal for D with common optimum (the value of the game) v.

Proof. p^{⊤}A ≥ v (which is understood to mean (p^{⊤}A)_{i} ≥ v in each component i) is just primal feasibility. Similarly, Aq ≤ v is dual feasibility. In addition,

v = p^{⊤}Aq ⟹ v = ∑_{i} p_{i} ∑_{j} a_{ij}q_{j}

⟹ 0 = ∑_{i} p_{i} (∑_{j} a_{ij}q_{j} − v) .

Similarly,

0 = ∑_{j} q_{j} (∑_{i} p_{i}a_{ij} − v) .

These are the complementary slackness conditions for this pair of problems. Thus p and q are optimal for P and D respectively. Clearly each has optimum v. Remarks.

1. These conditions allow us to check optimality. For small problems one can often use them to find the optimal strategies, but for larger problems it will be best to use some other technique to find the optimum (e.g., simplex algorithm). Note, however, that the problems P and D are not in a form where we can apply the simplex algorithm directly; the p_{i}s satisfy p_{i} ≥ 0, but v does not have a positivity constraint. Also the constraints are ∑_{i} a_{ij}p_{i} − v = 0 with r.h.s.= 0. It is possible, however, to transform the problem into a form amenable to the simplex algorithm.

(a) Add a constant k to each a_{ij} so that a_{ij} > 0 each i, j. This doesn’t change anything about the game except the value which is now guaranteed to be positive (v > 0).

(b) Change variables to x_{i} = p_{i}/v. We now have

max v s.t. ∑_{i} a_{ij}x_{i} ≥ 1, ∑_{i} x_{i} = 1/v, x_{i} ≥ 0,

which is equivalent to

min ∑_{i} x_{i} s.t. ∑_{i} a_{ij}x_{i} ≥ 1, x_{i} ≥ 0

and this is the type of problem we are used to.

2. Check that the general theorem above gives the right answer for a game with a saddle-point (a, b), (i.e., p_{a} = 1, q_{b} = 1 and other p_{i}, q_{j} = 0).

3. Two-finger Morra has an optimal solution p = q = (0, \frac{3}{5}, \frac{2}{5}, 0), v = 0, as can be easily checked. (But it is not uniquely optimal. Another optimal solution is p = q = (0, \frac{4}{7}, \frac{3}{7}, 0).) It is obvious that we expect to have p = q and v = 0 since the game is symmetric between the players (A = −A^{⊤}). A is called an anti-symmetric matrix.

9.4 Example: Colonel Blotto

Colonel Blotto has three regiments and his enemy has two regiments. Both commanders are to divide their regiments between two posts. At each post the commander with the greater number of regiments wins one for each conquered regiment, plus one for the post. If the commanders allocate equal numbers of regiments to a post there is a stand-off. This gives the pay-off matrix

Colonel Blotto Enemy commander

(2,0) (1,1) (0,2)
(3,0)  3  1  0
(2,1)  1  2 −1
(1,2) −1  2  1
(0,3)  0  1  3

Clearly it is optimal for Colonel Blotto to divide his regiments (i, j) and (j, i) with equal probability. So the game reduces to one with the payoff matrix

(2,0) (1,1) (0,2)
(3,0) or (0,3) \frac{3}{2}  1  \frac{3}{2}
(2,1) or (1,2) 0  2  0

You can try deriving the optimal solution by

(a) looking at player Colonel Blotto’s original problem: maximize {min_{j} ∑_{i} p_{i}a_{ij}}, i.e., maximize_{p} min\{\frac{3}{2}p, p + 2(1 − p)\},

(b) attempting to derive p, q, v from the conditions p^{⊤}A ≥ v, Aq ≤ v, p^{⊤}Aq = v, or

(c) converting the problem as explained above and using the simplex method.

For this game, p = (\frac{4}{5}, \frac{1}{5}), q = (\frac{1}{5}, \frac{3}{5}, \frac{1}{5}) and v = \frac{6}{5} is optimal. In the the original problem, this means that Colonel Blotto should distribute his regiments as (3,0), (2,1), (1,2), (0,3) with probabilities \frac{2}{5}, \frac{1}{10}, \frac{1}{10}, \frac{2}{5} respectively, and his enemy should distribute hers as (2,0), (1,1), (0,2) with probabilities \frac{1}{5}, \frac{3}{5}, \frac{1}{5} respectively.

10 Maximal Flow in a Network

10.1 Max-flow/min-cut theory

Consider a network consisting of n nodes labelled 1, …, n and directed edges between them with capacities c_{ij} on the arc from node i to node j. Let x_{ij} denote the flow in the arc i → j, where 0 ≤ x_{ij} ≤ c_{ij}.

PSfrag replacements v 1 n

Problem: Find maximal flow from node 1 to node n subject to the conservation of flow at nodes, i.e.,

maximize v s.t. 0 ≤ x_{ij} ≤ c_{ij},

and

flow out of ; node i − flow into ; node i = ∑_{j} x_{ij} − ∑_{j} x_{ji} = \begin{matrix}v & & i = 1 \\ 0 & as & i ≠ 1, n \\ −v & & i = n \end{matrix}

where the summations are understood to be over existing arcs only. v is known as the value of the flow.

This is obviously an LP problem, but with lots of variables and constraints. We can solve it more quickly (taking advantage of the special network structure) as follows.

Definition 10.1 A cut (S, S̄) is a partition of the nodes into two disjoint subsets S and S̄ with 1 ∈ S and n ∈ S̄.

Definition 10.2 The capacity of a cut

C(S, S̄) = ∑_{i∈S,j∈S̄} c_{ij}.

Thus given a cut (S, S̄) the capacity of the cut is the maximal flow from nodes in S to nodes in S̄. It is intuitively clear that any flow from node 1 to node n must cross the cut (S, S̄), since in getting from 1 to n at some stage it must cross from S to S̄. This holds for any flow and any cut. In fact, we have:

Theorem 10.1 (Max flow/min cut Theorem) The maximal flow value through the network is equal to the minimal cut capacity.

Example.

PSfrag replacements

2 1 3 4

Cut S = {1, 3}, S̄ = {2, 4}, C(S, S̄) = 3. Check that the maximal flow is 3. Proof. Let us define f(X, Y) = ∑_{i∈X,j∈Y} x_{ij}, the flow from X to Y. We will use N = {1, …, n} to represent the set of nodes in the network. Let (S, S̄) be a cut and {x_{ij}} a feasible flow with value v. Note that f(X, N) = f(X, S) + f(X, S̄). Since x_{ij} is feasible, we have

∑_{j} x_{ij} − ∑_{j} x_{ji} = \begin{matrix}v & & i = 1 \\ 0 & as & i ≠ 1, n \\ −v & & i = n \end{matrix} .

Summing this equality over i in S we obtain

v = ∑_{i∈S} ∑_{j}(x_{ij} − x_{ji})
= ∑_{i∈S} ∑_{j} x_{ij} − ∑_{i∈S} ∑_{j} x_{ji}
= f(S, N) − f(N, S)
= f(S, S) + f(S, S̄) − f(S̄, S) − f(S, S)
= f(S, S̄) − f(S̄, S)
≤ f(S, S̄)
≤ C(S, S̄)

So any flow ≤ any cut capacity, (and in particular max flow ≤ min cut). Now let f be a maximal flow, and define S ⊆ N recursively as follows:

(1) 1 ∈ S.

(2) If i ∈ S and x_{ij} < c_{ij}, then j ∈ S.

(3) If i ∈ S and x_{ji} > 0, then j ∈ S.

Keep applying (2) and (3) until no more can be added to S.

So S is the set of nodes to which we can increase flow. Now if n ∈ S we can increase flow along a path in S and f is not maximal. So n ∈ S̄ = N \ S and (S, S̄) is a cut. From the definition of S we know that for i ∈ S and j ∈ S̄, x_{ij} = c_{ij} and x_{ji} = 0, so in the formula above we get

v = f(S, S̄) − f(S̄, S) = f(S, S̄) = C(S, S̄).

So max flow = min cut capacity.

Corollary 10.2 If a flow value v = cut capacity C then v is maximal and C minimal.

The proof suggests an algorithm for finding the maximal flow.

10.2 Ford-Fulkerson algorithm

1. Start with a feasible flow (e.g., x_{ij} = 0).

2. Construct S recursively by the algorithm defined in the box above.

3. If n ∈ S then there is a path from 1 to n along which we can increase flow by

ϵ = min_{(ij)} max[x_{ji}, c_{ij} − x_{ij}] > 0.

where the minimum is taken with respect to all arcs i → j on the path.

Replace the flow by this increased flow. Return to 2.

If n ∉ S then the flow is optimal.

The algorithm is crude and simple; we just push flow through where we can, until we can’t do so anymore. There is no guarantee that it will be very efficient. With hand examples it is usually easy to ‘see’ the maximal flow. You just demonstrate that it is maximal by giving a cut with the same capacity as the flow and appeal to the min cut = max flow theorem.

The algorithm can be made not to converge if the capacities are not rational multiples of one another. However,

Theorem 10.3 If capacities and initial flows are rational then the algorithm terminates at a maximal flow in a finite number of steps. (Capacities are assumed to be finite.)

Proof. Multiply by a constant so that all capacities and initial flows are integers. The algorithm increases flow by at least 1 on each step.

Example: Failure to stop when capacities and initial flows are not rational

PSfrag replacements a b c d e f 0 1 w w^{2} w^{3} +w +w^{2}

1 − w − w^{2} = 0

The network consists of a square b c d e of directed arcs of capacity 1. The corners of the square are connected to a source at a and a sink at f by arcs of capacity 10. The initial flow of 1+w is shown in the first picture, where w = (\sqrt{5}−1)/2, so 1−w = w^{2}. The first iteration is to increase flow by w along a → c → b → e → d → f. The second increases it by w along a → d → e → b → f. The flow has increased by 2w and the resulting flow in the square is the same as at the start, but multiplied by w and rotated through 180^{o}. Hence the algorithm can continue in this manner forever without stopping and never reach the optimal flow of 40.

10.3 Minimal cost circulations

Definition 10.3 A network is closed if there is no flow into or out of the network.

Definition 10.4 A flow in a closed network is a circulation if ∑_{j} x_{ij} − ∑_{j} x_{ji} = 0 for each node i.

Most network problems can be formulated as the problem of finding a minimal cost circulation in a closed network where there are capacity constraints c_{ij}^{−} ≤ x_{ij} ≤ c_{ij}^{+} on arcs (i, j) and a cost per unit flow of d_{ij} in arcs (i, j). The full problem is

minimize ∑_{ij} d_{ij}x_{ij}

subject to ∑_{j} x_{ij} − ∑_{j} x_{ji} = 0, each i, and c_{ij}^{−} ≤ x_{ij} ≤ c_{ij}^{+}.

Definition 10.5 A circulation which satisfies the capacity constraints is called a feasible circulation.

There is a beautiful algorithm called the out-of-kilter algorithm which will solve general problems of this kind. It does not even require a feasible solution with which to start. (See Part IIA ‘Algorithms and Networks’.) In the next lecture we shall just derive conditions for a flow to be optimal.

We shall also see, although it should be obvious already, that the max flow problem that is studied in this lecture can be formulated as a minimal cost circulation problem.

11 Minimum Cost Circulation Problems

11.1 Sufficient conditions for a minimal cost circulation

Recall the minimum cost circulation problem:

minimize ∑_{ij} d_{ij}x_{ij}

subject to ∑_{j} x_{ij} − ∑_{j} x_{ji} = 0, each i, and c_{ij}^{−} ≤ x_{ij} ≤ c_{ij}^{+}.

Consider the Lagrangian for the problem of finding the minimum cost circulation. We shall treat the capacity constraints as the region constraints, so

X = {x_{ij} : c_{ij}^{−} ≤ x_{ij} ≤ c_{ij}^{+}}.

We introduce Lagrange multipliers λ_{i} (one for each node) and write

L(x, λ) = ∑_{ij} d_{ij}x_{ij} − ∑_{i} λ_{i} (∑_{j} x_{ij} − ∑_{j} x_{ji})

Rearranging we obtain

L(x, λ) = ∑_{ij}(d_{ij} − λ_{i} + λ_{j})x_{ij}.

We attempt to minimize L(x, λ) in X.

Provided c_{ij}^{−}, c_{ij}^{+} are finite we see that there is a finite minimum for all λ, achieved such that

x_{ij} = \{ \begin{matrix}c_{ij}^{−} if d_{ij} − λ_{i} + λ_{j} > 0 \\ c_{ij}^{+} if d_{ij} − λ_{i} + λ_{j} < 0 \end{matrix}  (1)

c_{ij}^{−} ≤ x_{ij} ≤ c_{ij}^{+} if d_{ij} − λ_{i} + λ_{j} = 0.  (2)

Theorem 11.1 If (x_{ij}) is a feasible circulation and there exists λ such that (x_{ij}), λ satisfy conditions (1) and (2) above, then (x_{ij}) is a minimal cost circulation.

Proof. Apply the Lagrangian sufficiency theorem.

Definition 11.1 The Lagrange multipliers λ_{i} are usually known as node numbers or potentials in network problems.

Definition 11.2 λ_{i} − λ_{j} is known as the tension in the arc (i, j).

11.2 Max-flow as a minimum cost circulation problem

The maximal flow problem we studied earlier can be set up as a minimal cost circulation problem. For each arc in the network we assign a capacity constraint 0 ≤ x_{ij} ≤ c_{ij}^{+} and all cost d_{ij} = 0. Add an arc from node n to 1 with no capacity constraint and cost −1.

PSfrag replacements

1 n

cost −1

v S S̄

The cost of the circulation is −v, so minimizing −v is the same as maximizing v. Let us seek node numbers λ_{i} which will satisfy optimality for this problem. Since arc (n, 1) has no capacity constraints, for a finite optimum we will require

d_{n1} − λ_{n} + λ_{1} = 0 ⟹ λ_{1} = λ_{n} + 1.

Let us set λ_{n} = 0, λ_{1} = 1. (Since it is only the differences in the λs that matter, we can pick one arbitrarily.) Let (S, S̄) be a minimal cut. Assign λ_{i} = 1 for i ∈ S and λ_{i} = 0 for i ∈ S̄. Now check that

(a) For i, j ∈ S or i, j ∈ S̄ ⟹ d_{ij} − λ_{i} + λ_{j} = 0 so x_{ij} can take any feasible value.

(b) For i ∈ S, j ∈ S̄ we have

d_{ij} − λ_{i} + λ_{j} = 0 − 1 + 0 = −1 ⟹ x_{ij} = c_{ij}^{+}.

(c) For i ∈ S̄, j ∈ S we have

d_{ij} − λ_{i} + λ_{j} = 0 − 0 + 1 = 1 ⟹ x_{ij} = 0.

But conditions (a)–(c) are precisely those satisfied by a maximal flow and minimal cut.

If we like, we can say that the Ford-Fulkerson algorithm in looking for a cut is trying to find node numbers and a flow to satisfy optimality conditions. Remark In many problems it is natural to take c_{ij}^{−} = 0, c_{ij}^{+} = ∞. In this case we will achieve a finite optimum only if d_{ij} − λ_{i} + λ_{j} ≥ 0 for each arc.

Theorem 11.2 For a minimal cost circulation problem with capacity constraints 0 ≤ x_{ij} < ∞ on each arc (i, j), if we have a feasible circulation (x_{ij}) and node numbers λ_{i} such that

d_{ij} − λ_{i} + λ_{j} ≥ 0, each (i, j), and
x_{ij} = 0 if d_{ij} − λ_{i} + λ_{j} > 0,

then (x_{ij}) is optimal.

Proof. Apply the Lagrangian sufficiency theorem.

Note: The optimality conditions imply (d_{ij} − λ_{i} + λ_{j})x_{ij} = 0 in this case (complementary slackness).

11.3 The transportation problem

Consider a network representing the problem of a supplier who has n supply depots from which goods must be shipped to m destinations. We assume there are quantities s_{1}, …, s_{n} of the goods at depots {S_{1}, …, S_{n}} and that the demands at destinations {D_{1}, …, D_{m}} are given by d_{1}, …, d_{m}. We also assume that ∑_{i} s_{i} = ∑_{j} d_{j} so that total supply = total demand. Any amount of goods may be taken directly from source i to destination j at a cost of d_{ij} (i = 1, …, n; j = 1, …, m) per unit. One formulation of the problem is

minimize ∑_{ij} d_{ij}x_{ij}

subject to ∑_{j} x_{ij} = s_{i} each i, ∑_{i} x_{ij} = d_{j} each j

with x_{ij} ≥ 0 each i, j

Here x_{ij} is the flow from S_{i} to D_{j}. The network looks like:

PSfrag replacements s_{1} s_{2} s_{n} d_{1} d_{2} d_{m}

n sources m destinations

with arcs (i, j), 0 ≤ x_{ij} < ∞, and cost d_{ij} per unit flow.

The Lagrangian for the problem can be written

L(x, λ, μ) = ∑_{ij} d_{ij}x_{ij} − ∑_{i} λ_{i} (∑_{j} x_{ij} − s_{i}) + ∑_{j} μ_{j} (∑_{i} x_{ij} − d_{j}) ,

where we label Lagrange multipliers (node numbers) λ_{i} for sources and μ_{j} for destinations. (We choose the sign of μ_{j} in this apparently unusual way to be consistent with our previous treatment. The demands d_{j} at the destination are really negative supplies.) Rearranging,

L(x, λ, μ) = ∑_{ij}(d_{ij} − λ_{i} + μ_{j})x_{ij} + ∑_{i} λ_{i}s_{i} − ∑_{j} μ_{j}d_{j}.

This will have a finite minimum in x_{ij} ≥ 0, and the minimum occurs with (d_{ij} − λ_{i} + μ_{j})x_{ij} = 0 on each arc. Thus the Lagrangian sufficiency theorem give the same optimality conditions as before.

Theorem 11.3 A flow x_{ij} is optimal for the transportation problem if ∃ λ_{i}, μ_{j} such that d_{ij} − λ_{i} + μ_{j} ≥ 0 each (i, j) and (d_{ij} − λ_{i} + μ_{j})x_{ij} = 0.

Proof. The Lagrangian sufficiency theorem applies. Remark It is no surprise that the same optimality conditions appear as in the minimal cost circulation problem. If we augment the transportation network by connecting all sources and all destinations to a common ‘artificial node’ by arcs where the flow is constrained to be exactly that which is required (and zero cost) we obtain the same problem as in minimal cost circulation form.

PSfrag replacements

cost d_{ij} 0 ≤ x_{ij} < ∞ cost d_{0i} = 0 cost d_{j0} = 0 arcs s_{i} ≤ x_{0i} ≤ s_{i} arcs d_{j} ≤ x_{j0} ≤ d_{j}

arcs 0 ≤ x_{ij} 0

The optimality conditions on the extra arcs are automatically satisfied by a feasible flow since x_{j0} = d_{j}, x_{0i} = s_{i} regardless of node numbers.

12 Transportation and Transshipment Problems

12.1 The transportation algorithm

1. Set out the supplies and costs in a table as below

D_{1} D_{2} D_{3} D_{4} s_{i}
S_{1}  5  3  4  6  8
S_{2}  2  7  4  1 10
S_{3}  5  6  2  4  9
d_{j} 6  5  8  8 27

2. Allocate an initial feasible flow (by North-West corner rule or any other sensible method). NW corner rule says start at top left corner and dispose of supplies and fulfill demands in order i, j increasing. In our case we get

6  2
3  7
1  8

In the absence of degeneracy (which we assume) we will obtain (m + n − 1) non-zero entries in a ‘stair-case’ arrangement.

Remark. In our network picture we have constructed a feasible flow on a spanning tree of m + n − 1 arcs connecting n sources and m destinations.

PSfrag replacements

1 2 3 6 7 8 0 4 −5 −3

A set of arcs is spanning if it connects all nodes. It is a tree if it contains no circuits. A spanning tree is the equivalent of a basic solution for this problem.

3. For optimality we require d_{ij}−λ_{i}+μ_{j} = 0 on any arc with non-zero flow. Set λ_{1} = 0 (arbitrarily) and then compute the remaining λ_{i}, μ_{j} by using d_{ij} − λ_{i} + μ_{j} = 0 on arcs for which x_{ij} > 0. On the table we have

\begin{matrix}λ_{i} \ μ_{j} & −5 & −3 & 0 & −2 \\ 0 & 6 5 & 2 3 & 4 & 6 \\ 4 & 2 & 3 7 & 7 4 & 1 \\ 2 & 5 & 6 & 1 2 & 8 4 \end{matrix}

The node numbers are also shown on the network version above. With non-zero flows forming a spanning tree we will always be able to compute uniquely all node numbers given one of them.

4. We now compute d_{ij} − λ_{i} + μ_{j} for all the remaining boxes (arcs). We obtain

\begin{matrix}0 & 0 & 4 & 4 \\ −7 & 0 & 0 & −5 \\ −2 & 1 & 0 & 0 \end{matrix}

which can be written elsewhere in the boxes if desired.

5. If all d_{ij} − λ_{i} + μ_{j} ≥ 0, then the flow is optimal. Stop.

6. If not, (e.g., i = 2, j = 1 gives −7) we attempt to increase the flow in arc (i, j) for some (i, j) such that d_{ij} − λ_{i} + μ_{j} < 0. We seek an adjustment of +ϵ to the flow in arc (i, j) which keeps the solution feasible (and therefore preserves total supplies and demands). In our case we do this by

\begin{matrix}6 − ϵ & 2 + ϵ & 0 & 0 \\ +ϵ & 3 − ϵ & 7 & 0 \\ 0 & 0 & 1 & 8 \end{matrix}

and pick ϵ as large as possible (without any flow going negative) to obtain the new flow (for ϵ = 3) of

3  5  0  0
3  0  7  0
0  0  1  8

There is only one way to do this in a non-degenerate problem. The operation is perhaps clearer in the network picture.

PSfrag replacements

1 2

2 + ϵ 3 − ϵ 6 − ϵ

7 8 0 4 −5 −3 +ϵ

We attempt to increase flow in the dotted arc. Adding an arc to a spanning tree creates a circuit. Increase flow around the circuit until one arc drops out, leaving a new spanning tree. The new solution is

PSfrag replacements

1 2 3 5 7 8

7. Now return to step 3 and recompute node numbers.

In our example we obtain λ_{i} = 0, −3, −5 and μ_{j} = −5, −3, −7, −9 at the next stage. The expression d_{ij} − λ_{i} + μ_{j} < 0 for (i, j) = (1, 3), (2, 4) and (1, 4). Increase the flow in (2, 4) by 7 to obtain the new flow below. This is now optimal. (Check!)

3  5  0  0
3  0  0  7
0  0  8  1

Remark. The route around which you may need to alter flows can be quite complicated though it is always clear how you should do it. For example, had we tried to increase the flow in arc (3, 1) instead of (2, 3) at step 5 we would have obtained

\begin{matrix}6 − ϵ & 2 + ϵ & 0 & 0 \\ 0 & 3 − ϵ & 7 + ϵ & 0 \\ +ϵ & 0 & 1 − ϵ & 8 \end{matrix}

and ϵ = 1. Try this out on a picture of the network and complete the algorithm from this position.

To summarise:

1. Pick initial feasible solution with m + n − 1 non-zero flows (NW corner rule).

2. Set λ_{1} = 0 and compute λ_{i}, μ_{j} using d_{ij} − λ_{i} + μ_{j} = 0 on arcs with non-zero flows.

3. If d_{ij} − λ_{i} + μ_{j} ≥ 0 for all (i, j) then flow is optimal.

4. If not, pick (i, j) for which d_{ij} − λ_{i} + μ_{j} < 0.

5. Increase flow in arc (i, j) by as much as possible without making the flow in any other arc negative. Return to 2.

12.2 *Simplex-on-a-graph*

The transportation algorithm can easily be generalised to a problem of minimizing costs in a general network with constraints 0 ≤ x_{ij} < ∞ on each arc and flows b_{i} into the network at each i (though it is hard to keep track of all the numbers by hand). Here we don’t label sources and destinations separately, but do allow b_{i} ≥ 0 and b_{i} ≤ 0. Clearly, ∑_{i} b_{i} = 0 for conservation of flow. The simplex-on-a-graph algorithm solves this problem in an identical fashion to the transportation algorithm. Once again a basic solution is a spanning tree of non-zero flow arcs.

1. Pick an initial basic feasible solution. Obtain n − 1 non-zero flow arcs.

2. Set λ_{1} = 0 and compute λ_{i} on other nodes using d_{ij} − λ_{i} + λ_{j} = 0 on arcs of the spanning tree.

3. Compute d_{ij} − λ_{i} + λ_{j} for other arcs. If all these are ≥ 0 then optimal. If not, pick (i, j) such that d_{ij} − λ_{i} + λ_{j} < 0.

4. Add arc (i, j) to the tree. This creates a circuit. Increase flow around the circuit (in direction of arc (i, j)) until one non-zero flow drops to zero and a new basic solution is created. Return to 2.

If it is hard to find an initial basic solution then there is a two-phase version of the algorithm (just as for the ordinary simplex algorithm).

Phase I. Add artificial node 0 with arcs to all nodes with b_{i} > 0 and from all nodes with b_{i} < 0. For Phase I objective put costs of 1 on all arcs to and from node 0 and costs 0 elsewhere in the network. The initial basic solution for Phase I is the spanning tree consisting of the node 0 and all arcs joining it to the original nodes.

At the end of Phase I (if the original problem was feasible) you will have reduced the Phase I cost to 0 (no flow in arcs to or from node 0), so have a basic solution (Spanning tree solution) for the original problem.

Phase II. Solve the original problem using the initial solution found by Phase I.

There are several applications of the network theory to problems of graph theory and operations research on the further examples sheet.

12.3 Example: optimal power generation and distribution

The following real-life problem can be solved by the simplex-on-a-graph algorithm.

PSfrag replacements

1 2 3 4 5 6 7 8 9 10 11 12 London Scotland Cumbria Northeast Northwest Yorkshire Wales Central Southwest Southcoast Thames Midlands

The demand for electricity at node i is d_{i}. Node i has k_{i} generators, that can generate electricity at costs of a_{i1}, …, a_{ik_{i}}, up to amounts b_{i1}, …, b_{ik_{i}}. There are n = 12 nodes and 351 generators in all. The capacity for transmission from node i to j is c_{ij} (= c_{ji}).

Let x_{ij} = amount of electricity carried i → j and let y_{ij} = amount of electricity generated by generator j at node i. The LP is

minimize ∑_{ij} a_{ij}y_{ij}

subject to ∑_{j} y_{ij} − ∑_{j} x_{ij} + ∑_{j} x_{ji} = d_{i}, i = 1, …, 12,

0 ≤ x_{ij} ≤ c_{ij}, 0 ≤ y_{ij} ≤ b_{ij}.

In addition, there are constraints on the maximum amount of power that may be shipped across the cuts shown by the dotted lines in the diagram.