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.
No comments:
Post a Comment