Linear programming - the auxiliary LP

 Yesterday, I wrote a series of posts explaining how to solve a linear programming problem of the form

• Minimize the objective function $f(x)=c\cdot x$
• subject to the constraints
• $Ax=b$ where $A$ is an $m\times n$ matrix
• $x\geq 0$

The simplex algorithm works great, but you need an initial basic feasible solution to get started. Luckily, one can construct a Phase I program to decide if any feasible solutions exist, and provide such a solution if so. The last remaining piece is to turn this into a basic feasible solution, which we can do with the help of an auxiliary LP.

Recall that the phase I LP is the following

Phase 1 LP (standard)

Introduce new variables $y=(y_1,\ldots,y_n)$.

• Minimize $\sum_{i=1}^n y_i=e\cdot y$ where $e=(1,\ldots,1)$
• subject to
• $Ax+Iy=b$
• $x,y\geq 0$

If we have an optimal solution $(x^*,y^*)$ with $y^*=0$, then $x^*$ is a solution to the original problem. If all of the $y^*$ variables are non-basic, then we can easily construct a bfs for the original LP. However it is quite possible that this is not the case. To remedy this, we introduce the auxiliary LP.

Auxiliary LP

Start with the phase 1 LP and add one new variable $z$.

• Minimize $c\cdot x$
• subject to
• $Ax+Iy=b$
• $e\cdot y+z = 0$ 
• $x,y,z\geq 0$

The auxiliary LP is useful for the following reasons:

• Fact 1: If $(x,y,z)$ is optimal for the auxiliary LP, then $x$ is optimal for the original LP.
• Fact 2: If the auxiliary LP is unbounded, then the original LP is unbounded.
• Fact 3: If $(x^*,0)$ is an optimal solution to the phase 1 LP, then $(x^*,0,0)$ is a bfs for the auxiliary LP. The basic variables are the same as those for $(x^*,0)$ with $z$ added.

Facts 1 and 2 basically come from the fact that any solution to the auxiliary LP must have $y=0$ and $z=0$, because the sum of nonnegative terms $e\cdot y+z$  is equal to $0$. 

The main thing about Fact 3 is to check that the new matrix, after appending the column $\left(\begin{matrix} 0\\ 1\end{matrix}\right)$, is invertible, which I leave as an exercise.

So now we have a complete algorithm to solve LPs of the form indicated at the start of the post. I've implemented this algorithm in C# and it works quite well. The one issue that I kept grappling with is the issue of thresholds, and I essentially just kept tinkering with them until the program stopped giving me error messages. The linear systems I have been dealing with are often highly degenerate, and one thing I am picking up from these posts is that I may want to row reduce the coefficient matrix and throw away extraneous rows before implementing the algorithm. 

So now the question is, what's next? I would like to understand how to do quadratic programming, but I would also like to run through at least one real life (but smallish) example of the LP algorithm we just described both for illustrative purposes, but also so I can see how I might be able to improve it.