# Finding a set of maximally different solutions using linear programming or other optimization technique

Traditionally, linear programming is used to find the one optimal solution to a set of constraints, variables and a goal (all described as linear relationships). Sometimes, when the objective is parallel to a constraint, there are infinite or many equally good optimal solutions. I am not asking about this latter case.

I am more interested in finding many solutions that are in the feasible region generated by my set of constraints. But I would like the solutions I find to be 'scattered' around the feasible region in the sense that they are maximally far from one another. Is there a known way to, without running a solver multiple times, generate multiple solutions and use the objective function to enforce that the solutions should be separated?

For example, any linear program with decisions a and b and constraints w <= a <= x and y <= b <= z can be 'duplicated' to find two solutions. Our new linear program has variables a1, a2, b1, and b2 and the constraints w <= a1 <= x and w <= a2 <= x and similar for b1, b2. However, when it comes to forming an objective function we run into trouble as we can not use norms other than the L1-norm without discarding linearity and we can not truly even use the L1 norm because it is not possible (so far as I know) to encode absolute values.

Perhaps I should look into convex optimization or semidefinite programming or something?

Is there a known way to generate a set of solutions to a linear program, and using an objective that enforces "distance" between the solutions?

• Compute the smallest cube surrounding your feasible region (if it's unbounded, choose some bounded part), lay a hypergrid with the desired resolution over it and discard all points that don't fulfill the restrictions. Would that work for you? – Raphael Apr 11 '13 at 22:47
• That might work for me, although it is not clear to me how I would go about computing the hypercube, and I think the feasible region I am investigating is highly non-trivial - I expect that many points would have to be discarded. My particular application has tens of thousands of variables/decisions and hundreds of constraints. – Ross Apr 11 '13 at 23:00
• A basic feasible solution is in a vertex of the polytope. Can you not look at the normals of the incident faces to compute a direction "across" the polytope and follow it to the boundary of the feasible region? That should give you reasonably different solutions, but probably not the most different ones. – adrianN Aug 10 '13 at 16:22
• Don't use a cube. Use an ellipsoid (specifically, a small ellipsoid containing the polytope, which can be found by the ellipsoid method). That way, you're guaranteed to find a reasonable number of points in the region. – Peter Shor Sep 6 '13 at 23:32

## A heuristic, using linear programming

One approach might be to pick a random objective function, and maximize it. Then repeat, with a different set of objective functions each time.

In other words, suppose the unknowns are $x_1,x_2,\dots,x_n$, and you have some constraints $\mathcal{C}$. In each iteration you pick $c_1,c_2,\dots,c_n \in \mathbb{R}$ randomly, then search for a solution that maximizes $c_1 x_1 + \dots + c_n x_n$ subject to the constraints $\mathcal{C}$.

I would expect this heuristic might often find a somewhat scattered set of solutions -- not necessarily maximally scattered (maximally far from each other), but probably not too close to each other, either.

## Maximizing the average pairwise L2 distance, using quadratic programming

Alternatively, use quadratic programming. For simplicity, let's look at the problem of finding two solutions. Suppose you want two solutions $x,y$ that are as far apart from each other as possible, under the $L_2$ norm (Euclidean distance). Then this can be formulated as a quadratic programming problem.

Basically, you want to maximize the squared distance $d(x,y)^2 = (x_1-y_1)^2+\cdots+(x_n-y_n)^2$ between $x$ and $y$, subject to the requirement that both $x$ and $y$ must satisfy the constraints. This is the problem of maximizing a quadratic objective function, with linear constraints -- i.e., quadratic programming.

If you want $k$ points that are maximally scattered, this is also possible. Say the points are $x^1,\dots,x^k \in \mathbb{R}^n$. Then you could maximize the objective function

$$\sum_{i<j} d(x^i,x^j)^2,$$

i.e., the function

$$\sum_{i<j} \sum_\ell (x^i_\ell - x^j_\ell)^2.$$

This is a quadratic function, and you have linear constraints $\mathcal{C}$ on each of the points $x^i$, so this is a quadratic programming instance. It finds you points that are maximally scattered in the sense that the average pairwise distance is maximized.

You can also formulate a greedy variant of this algorithm, where you already have $k$ solutions, and you want to find a $k+1$th solution that satisfies all the linear inequalities and also maximizes the average L2 distance from it to the other $k$ solutions. That too can be formulated as a quadratic programming problem.

Quadratic programming is harder than linear programming, but there are off-the-self solvers that will solve quadratic programming problems for you.

## Maximizing the minimal pairwise L2 distance, using QCQP

Finally, let's say you want your $k$ points to be scattered in the sense that you want to maximize the minimum pairwise distance. In other words, let's say you want to find the largest possible threshold $t$ such that it is possible to find $k$ points $x^1,\dots,x^k \in \mathbb{R}^n$ that each satisfy the linear constraints, and such that each pair of points is at distance $t$ away from each other: $d(x^i,x^j) \ge t$ for all $i<j$. Then this can be formulated as a quadratic optimization program with quadratic constraints, i.e., QCQP. QCQP is even harder, but there are off-the-shelf solvers for QCQP you could try out, too.

I found an approach to generating absolute values.

Suppose we have the variables $a_1$, $a_2$, $b_1$ and $b_2$ and a bunch of constraints. Our objective functions looks something like: maximize $|a_1 - a_2| + |b_1 - b_2|$; the idea being we want to maximize the L1 norm of these two solutions (as per the original question).

We can introduce "slack variables" abs_a and abs_b and the constraints:

$$\mathrm{abs_a} + a_1 - a_2 \leq 0$$

$$\mathrm{abs_a} - a_1 + a_2 \leq 0$$

and similarly for $b_1$ and $b_2$. These constraints force $\mathrm{abs_a}$ to be at most the difference between $a_1$ and $a_2$, and possibly less. In other words $\mathrm{abs_a}$ can not be larger than the maximum difference between $a_1$ and $a_2$.

Then what is left is replace the objective function: maximize $\mathrm{abs_a} + \mathrm{abs_b}$.

• Actually, these encodings only work for the minimization of absolute values. Therein it doesn't solve my problem. More information here: lpsolve.sourceforge.net/5.5/absolute.htm – Ross Apr 12 '13 at 18:39