Finding Multiple Extreme Rays

I've done a few posts relating to extreme rays before: "Finding Extreme Rays in CPLEX"; "Farkas Certificates in CPLEX"; "Infeasible LPs and Farkas Certificates". Today the following question arose: given an unbounded linear program, can we get the solver to provide us with multiple extreme rays in a single operation. (The "dual question" to this is whether, given an infeasible primal problem, we can obtain multiple Farkas certificates in one gulp.)

I'll answer in the context of CPLEX, making no claims as to whether this applies to other solvers (though I suspect it does). I don't think there's any way to get multiple rays from one solution attempt without additional computations, but there is a way to get multiple extreme rays (though possibly not all of them) by modifying the objective function of the original problem, while exploiting the feasibility of the last known primal solution (the one that convinced the solver that the problem is unbounded). Hot-starting from previous solutions should reduce the amount of work involved in finding extra rays, but the additional work will not be zero.

Suppose that we solve the linear program \[ \begin{array}{lrcl} \mathrm{maximize} & & c^{\prime}x\\ \mathrm{s.t.} & Ax & \le & b\\ & x & \ge & 0 \end{array} \] and discover that it is unbounded. The CPLEX getRay method will then provide us with one extreme ray. To find additional rays, we need to modify the objective function so that (a) other recession directions continue to be favorable and (b) the recession directions we have already obtained cease to be favorable. The former suggests that we want a new objective vector acute to the original one; the latter means we need the new direction to be orthogonal or obtuse to the recession directions we have already identified.

Let $\mathcal{D}$ be the set of recession directions we know at any given stage. One way to find a new objective coefficient vector $\tilde{c}$ meeting our criteria is to solve the following auxiliary linear program:\[ \begin{array}{lrclr} \mathrm{maximize} & & c^{\prime}\tilde{c}\\ \mathrm{s.t.} & d^{\prime}\tilde{c} & \le & 0 & \forall d\in\mathcal{D}\\ & -1 & \le\tilde{c}\le & 1 \end{array}. \]The objective function encourages the new direction $\tilde{c}$ to be acute to the original direction $c$, while the constraints ensure that the new objective function ($\tilde{c}^\prime x$), which we maximize subject to the original constraints $Ax\le b,x\ge 0$, does not increase along any of the recession directions we have already identified.

As an example, consider the linear program\[ \begin{array}{lrcl} \mathrm{maximize} & 3x+2y+z\\ \mathrm{s.t.} & 2x-y-3z & \le & 4\\ & x+2y-5z & \le & 3\\ & -4x+2y+z & \le & 6\\ & -x-y+2z & \le & 2\\ & x,y,z & \ge & 0. \end{array} \]Using the approach I just described, I am able to find three extreme directions:

Objective direction	Recession direction
(3, 2, 1)	(5, 1, 3)
(0.4, 1, -1)	(2.2, 1.4, 1.0)
(-0.1818, 1.0, -1.0)	(1.0, 1.5833, 0.8333)
(0.2, 0.4, -1.0)	OPTIMAL

With the fourth objective vector, the modified LP becomes bounded, and we halt the search.

Did we find all the extreme directions of the original problem? I think not. Since our example lives in $\mathbb{R}^3$, each extreme direction will correspond to an edge, the intersection of two of the seven bounding hyperplanes. (The seven bounding hyperplanes are the hyperplanes corresponding to the four functional constraints, plus the coordinate hyperplanes $x=0$, $y=0$, $z=0$ corresponding to the sign restrictions.) The recession directions in the second column of the table are the intersections of the first and fourth constraint hyperplanes, the first and second, and the second and third respectively. The intersection of the first and third constraint hyperplanes gives a direction $(1, 2, 0)$ which is not a recession direction (movement eventually is blocked by the second constraint). Similarly, the intersection of the second and fourth constraint hyperplanes yields a direction $(-1, 3, 1)$ that is eventually blocked by either the third constraint or the sign restriction on $x$, depending on where in the feasible region you start moving in that direction. The intersection of the third and fourth constraint hyperplanes, however, yields the direction $(5, 7, 6)$, which is not blocked by anything and along which the original objective increases.

So this method is not guaranteed to find all extreme directions of an unbounded linear program, but it will find some. The computational cost is repeated solution of a relatively small auxiliary linear program (the same number of variables as the original LP, but likely considerably fewer constraints in practice) and repeated re-solving of the original problem with modified objective functions (where each previous solution gives a feasible starting solution).

If anyone cares, I have uploaded Java code for this demonstration.