[Help-glpk] Re: convex objective function, linear integer constraints

[Peter T. Breuer]

"Also sprach ptb:"
> I'm interested in maximising a certain convex objective function over a
> set of linear constraints. I want to find the x which maximises (or
> rather, the maximum value of)
> 
>    MIN_j(MAX(0, a_j1 - x_1, x_1 - b_j1), ..., MAX(0, a_jN - x_n, x_n - b_jN))
> 
> as vector x varies over a convex linearly bounded integer n-dim
> subspace.  The a_jk < b_jk are integer constants.
> 
> In fact, x is usually constrained to a simple cube [a_j0,b_j0], and the
> expression above is a measure of the minimum distance to a set of n
> other cubes.

I'll follow up on this. As I said ...

> The objective is piecewise linear.

But I was wrong about the convexity.

Andrew Makhorin tried to help, mistakenly believing my idiocity:

  > Since the objective is piece-wise linear and convex, you can
  > reformulate the problem as follows:

(incidentally, how DOES one canonically do minimax problems using GLPK?
For minimax this problem certainly is, being NP-complete; travelling
salesman is easy to encode in it).

I've since managed to make some progress by reading up on what appear
to be called SOS2 problems (a recent thread here drew my attentioon to
that), which have given me some ideas.

As I wrote to Andrew, the distance yij of an arbitrary point x from cube i
in the jth dimension looks like this:

 yij
  |
  | \       /            yij = 0                         (aij <= xj <= bij)
  |  \     /             yij = aij - xj for j = 1 ... n  (aij >= xj)
  |   \___/              yij = xj - bij for j = 1 ... n  (bij <= xj)
  .__________xj
  L  aij bij R

So one can see that I got convex/concave mixed up, probably.

The diagram shows how the distance to the cube with sides at xj=aij
and xj=bij varies as one moves further or closer in that (jth) dimension
normal to those sides.

L and R are just some fixed bounds on the problem. But the idea I got
from SOS2 formulations was to invent variables wij1 wij2 wij3 which
are 1 respectively when xj is in the range L to aij, aij to bij, bij to
R. And to reparametrise the problem on extra variables zij1 zij2 zij3
zij4 which are weights for where xj falls in terms of an average of the
4 points L, aij, bij, R.

That is

    let  xj = zij1 L + zij2 aij + zij3 bij + zij4 R
    with      zij1, zij2, zij3, zij4 >= 0, zij1 + zij2 + zij3 + zij4 = 1

and force only two consecutive of zij1 zij2 zij3 zij4 to be nonzero
(SOS2?). One can force that by writing

     wij1 = zij1 + zij2
     wij2 = zij2 + zij3
     wij3 = zij3 + zij4

     wij1, wij2, wij3 >= 0
     wij1 + wij2 + wij3 = 1

and specifying that wij1, wij2, wij3 are integer.

Then the idea is that wijk will solve as being 1 or 0, and exactly one
of wij1 wij2 wij3 is 1. That means that xj is in the corresponding
arc L to aij, aij to bij, bij to R, and yij can then be written

    yij = zij1 (aij - L) + zij4 (R - bij)

exactly.

Now we're pumping, because the distance of an arbitrary point x from
the ith cube (of N) (in n dimensions) is then just

    yi = yi1 + ... + yin

and if one wants to find the point x which is furthest from ANY of the
cubes, then one just maximises v subject to

    v <= y1, ... v <= yN

while restricting x to the zone of interest. If v comes out > 0, then
there is part of the zone of interest not covered by the given cubes.
If v comes out as 0, then all of the zone of interest is covered by the
given cubes.


My difficulty seems to be that I can set up this problem for glpk,
albeit a bit awkwardly, and solve it using the general simplex method,
but the answer comes out with wijk not integer. Usually they're some
value like 0.9/0.1. If I then push the problem into lpx_intopt, it
complains about there being no "primal simplex solution", or some such.

And  I can't really debug 300-odd variable problems very well at this
stage in my expertise (30, I can probably cope with).

The way I have set it up is as follows, if anybody cares to know:

     maximise v

var        v
constraint 0 <= v        // can leave that unconstrained

var        zijk          for i = 1 to N,  j = 1 to n, k = 1 to 4
constraint 0 <= zijk <= 1
                         for i = 1 to N,  j = 1 to n, k = 1 to 4

intvar     wijk          for i = 1 to N,  j = 1 to n, k = 1 to 3
constraint 0 <= wijk <= 1

var        xj            for j = 1 to n
constraint a0j <= xj <= b0j
                         for j = 1 to n

var        yij           for i = 1 to N,  j = 1 to n
constraint 0 <= yij      // can leave this unconstrained



Then I have to awkwardly express the parameterisation of the xij in terms of
more fundamental zij1 .. zij4 in this roundabout way:


row        eij           for i = 1 to N, for j = 1 to n
constraint 0 = eij
                         for i = 1 to N, for j = 1 to n
equation   eij = zij1 (L - 0.5) + zij2 (aij - 0.5) + zij3 (bij + 0.5) + zij4 (R 
+ 0.5) - xj
                         for i = 1 to N, for j = 1 to n


The extra 0.5's are in there because this is really an integer problem and
the cubes have to be deemed to extend out to the next half-integer mark
when embedded in a real-valued domain.  If two cubes are as close as 1
apart, then notionally they touch!

Saying foo = a - b where foo = 0 is a pretty roundabout way of getting
a = b, but I don't seem to have much choice in glpk.  If I directly had
the xj on the left (as 'a'), then they'd be rows.  But that's notionally
wrong since they're the free variable that oen has to move around to
find the optimum with.

OK, they've been parameterised in terms of the zijk, so maybe it would
make sense here, but it doesn't make sense in the next set of equations,
where I want to write wij1 = zij1 + zij2 and then apply more equations
to the wijk.

So what I'm grousing about is that I can't write a = bar(b) where b =
gum(c).  I.e.  I can't use row names on the right hand side of
equations.  Instead I have to expand out the substitutions till I get
down to column variables.

Here's the clumsy method of writing wij1 = zij1 + zij2, via
foo = wij1 - zij1 - zij2, where foo = 0.

row        vijk          for i = 1 to N,  j = 1 to n, k = 1 to 3
equation   vij1 = wij1 - zij1 - zij2
           vij2 = wij2 - zij2 - zij3
           vij3 = wij3 - zij3 - zij4
constraint vijk = 0      for i = 1 to N,  j = 1 to n, k = 1 to 3

Since the wijk are miraculously preserved asvariables (columns), not
rows, this way I can apply the restriction I wanted, that wij1 + wij2 +
wij3 = 1, as an equations.

row        rij           for i = 1 to N,  j = 1 to n
constraint 1 = rij       for i = 1 to N,  j = 1 to n
equation   rij = wij1 + wij2 + wij3
                         for i = 1 to N,  j = 1 to n


Now I get to express the yij in the same dorky way, just because I
can't write yij = and then apply further equational constraints to the
yij.


row        dij           for i = 1 to N,  j = 1 to n
constraint 0 = dij       for i = 1 to N,  j = 1 to n
equation   dij = zij1 (aij - L) + zij4 (R - bij) - yij
                         for i = 1 to N,  j = 1 to n

and finally I have the fun of setting v to lie below the measure of
distance to each of the cubes:


row        ui            for i = 1 to N
constraint 0 <= ui       for i = 1 to N
equation   ui = yi1 + .. + yin - v
                         for i = 1 to N


and that works. With teh exception that the integer constraint on the wijk
seems not to go down a treat.




PTB