Currently I'm trying to write a practically fast LP solver for a sparse instance, which is by simplex method with LU decomposition and eta-matrix update. In the development I realized that I'm not familiar with how to choose the initial basis for simplex method, and would like some tips.
Think about a standard form of LP: maximize $cx$ s.t. $Ax \le b, x \ge 0$. In this case we can choose a obvious basis, i.e., columns of slack variables. For example, if $A = \begin{bmatrix} 1 & 2 \\ 3 & 4\end{bmatrix}$, then we extend $A$ with slack variables ($A' = \begin{bmatrix} 1 & 2 & 1 & 0\\ 3 & 4 & 0 & 1\end{bmatrix}$) and can choose the columns 3 and 4 as an initial basis. (It doesn't necessarily provide a basic feasible solution but it doesn't matter here.)
A rough version of my question is, what is the best way to choose the initial basis for equality constraints? I know a transform from an equality constraint to inequality ones at least works. For example, we can replace an equality constraint $x_1 + 2x_2 + 3x_3 = 4$ with two inequality contraints $x_1 + 2x_2 + 3x_3 \le 4, x_1 + 2x_2 + 3x_3 \ge 4$ and add two columns for slack variables. Naively thinking, however, we can directly apply simplex method against an equality form of LP (i.e. maximize $cx$ s.t. $Ax = b, x \ge 0$) without any transform, if we can choose a set of linearly independent columns of matrix $A$. Of course $A$ can be rank-deficient, but we can at least reduce the number of extended columns if we know the maximal linearly independent columns of $A$. Is this `basis' problem as hard as e.g. LP?
The reasons I'd like to avoid a transform to two constraints are:
- It increases the size of instance,
- and it produces yet another problem related to numerical errors. (That is, we need to choose appropriate epsilon to avoid infeasibility caused by two inequality constraints with the same bound.)
A bit detailed version of my questions are as follows:
- How do well-known LP solvers deal with equality constraints? Do they just replace the equality with inequality? Are my concerns above of little practical importance?
- (If the answer to 1 requires solving `basis' problem,) is there any practically efficient way to solve the basis problem above for a sparse instance? For a dense instance, I know e.g. Gaussian elimination works.
