Yet Another Math Programming Consultant: The Hostile Brothers Problem (3)

Follow up on http://yetanothermathprogrammingconsultant.blogspot.com/2015/01/the-hostile-brothers-problem-1.html and http://yetanothermathprogrammingconsultant.blogspot.com/2015/01/the-hostile-brothers-problem-2.html.

Now some experiments that did not work.

Find global solutions with Baron

This is a very difficult job. One idea that may help is to add some symmetry breaking constraints:

This can reduce the search space considerably. We run the equation over all i–1: in GAMS this is notation that means we drop the first possible equation (the equation related to ord(i)=1).

In addition we can pass on a good solution we found with techniques described earlier. Unfortunately this helps but not enough to bring us close to proving global optimality:

Iteration Open nodes Total time Lower bound Upper bound 1 1 000:00:01 0.100000E-08 2.00000 1 1 000:00:07 0.200000E-08 2.00000 * 2 2 000:00:10 0.163965E-01 2.00000 * 2 2 000:00:18 0.173057E-01 2.00000 6 4 000:00:51 0.173057E-01 1.67193 12 8 000:01:29 0.173057E-01 1.56665 17 12 000:02:31 0.173057E-01 1.53265 21 15 000:03:09 0.173057E-01 1.50781 25 18 000:03:58 0.173057E-01 1.49840 29 20 000:04:29 0.173057E-01 1.49070 33 22 000:05:23 0.173057E-01 1.49070	Default. Note that Baron actually finds good solutions quickly. But after a little while the lower bound will not move again and the upper bound will move slower and slower.
Iteration Open nodes Total time Lower bound Upper bound 1 1 000:00:02 0.363949E-03 2.00000 * 1 1 000:00:24 0.111980E-02 1.01351 * 1 1 000:00:29 0.139747E-01 1.01351 * 1 1 000:00:43 0.151056E-01 1.01351 1 1 000:00:49 0.151056E-01 1.01351 3 2 000:01:23 0.151056E-01 1.01351 4 3 000:02:00 0.151056E-01 1.01351 6 4 000:02:42 0.151056E-01 1.00880 8 5 000:03:29 0.151056E-01 1.00802 11 7 000:04:17 0.151056E-01 1.00802 12 8 000:04:50 0.151056E-01 1.00802 15 8 000:05:35 0.151056E-01 1.00191	Here we added symmetry breaking constraints. This helps with the upper bound.
Iteration Open nodes Total time Lower bound Upper bound 1 1 000:00:04 0.175061E-01 2.00000 1 1 000:00:19 0.175061E-01 1.01351 4 3 000:00:52 0.175061E-01 1.01351 6 4 000:01:30 0.175061E-01 1.00930 8 6 000:02:13 0.175061E-01 1.00930 10 7 000:02:55 0.175061E-01 1.00783 11 7 000:03:26 0.175061E-01 1.00783 13 7 000:04:01 0.175061E-01 1.00783 15 8 000:04:39 0.175061E-01 1.00244 17 10 000:05:13 0.175061E-01 1.00000	Here we use a starting point and symmetry breakers. This gives best performance but still not sufficient to get close to a global solution.

Understandably global solvers can not handle just anything we throw at them. Heaving said this I think Baron does a good job finding good solutions if we don’t provide our own initial points.

Other Global solvers

None of the global solvers can solve this problem quickly. It is interesting to see which objective can be used with different global solvers.

Solver	Discontinuous Formulated allowed	Smooth Formulation Accepted
Baron	Can not handle function ‘min’	OK
LocalSolver	OK	OK (n=75 problem ran out of demo license)
OQNLP	OK	OK
LindoGlobal	OK automatic reformulated into bigM constraints	OK (n=75 problem: Model size exceeds LindoGlobal limits of (2000,3000)).
LGO	OK	OK (n=75 problem:Reduced Model Size exceeds LGO limits of (2000,3000)).

Can we solve as a MIP?

The quadratic expressions form a problem. Although Cplex can solve non-convex QPs it does not handle non-convex quadratic constraints (I believe SCIP does). We can use some piecewise linear approximation, but it seems to make sense to first try a linear version with a different distance function: the Manhattan distance. Here the distance between two points is the sum of the distance in the x direction and in the y direction. The formulation can look like: