# NMaximize selects duplicate inputs and a problem with EvaluationMonitor

The aim of the following is to maximize 'penalty' and monitor it for convergence during the process:

n = 4;
lines = n (n - 1)/2;
optimalelements = n - 1;
gridpoints = 100;
grid = Round[Sqrt[gridpoints]];
fitness[x2_, y2_, x3_, y3_, x4_, y4_] :=
(Clear[fitness, points, linepoints, d, c, penalty, threepoints,
threetest, fourpoints, fourtest, p];
points = {{0, 0}, {x2, y2}, {x3, y3}, {x4, y4}};
linepoints = Subsets[points, {2}];
d = EuclideanDistance @@@ linepoints;
c = Sort[Tally[d][[All, 2]]];
penalty = Total[Abs[Differences[c]]] - Length[c];
If [c == Range[optimalelements], penalty = penalty + 1];
p = penalty);
{sol, pts} = Reap[
NMaximize[{fitness[x2, y2, x3, y3, x4, y4], 0 <= x2 <= grid,
0 <= y2 <= grid, 0 <= x3 <= grid, 0 <= y3 <= grid,
0 <= x4 <= grid,
0 <= y4 <= grid}, {{x2, 0, grid}, {y2, 0, grid}, {x3, 0,
grid}, {y3, 0, grid}, {x4, 0, grid}, {y4, 0, grid}}, Integers,
Method -> {"SimulatedAnnealing", "SearchPoints" -> 1,
"PerturbationScale" -> 1, "RandomSeed" -> 1},
EvaluationMonitor :> Sow[{{x2, y2}, {x3, y3}, {x4, y4}, c, p}]]] //
AbsoluteTiming

Out[7] {0.865882, {{-6., {x2 -> 5, y2 -> 6, x3 -> 2, y3 -> 4, x4 -> 7,
y4 -> 4}}, {{{{0, 4}, {0, 4}, {2, 6}, {1, 1, 1, 1, 1,
1}, -6}, {{0, 6}, {0, 2}, {4, 5}, {1, 1, 1, 1, 1, 1}, -6}, {{2,
6}, {0, 4}, {4, 8}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 6}, {0,
2}, {4, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{3, 7}, {1, 4}, {3,
7}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 9}, {2, 5}, {2, 7}, {1, 1, 1,
1, 1, 1}, -6}, {{5, 10}, {4, 7}, {0, 4}, {1, 1, 1, 1, 1,
1}, -6}...etc


There are two problems:

1. How does one stop NMaximize selecting duplicate input points (e.g. {0, 4}, {0, 4}, {2, 6})?
2. Why are c, p (i.e. {1, 1, 1, 1, 1, 1}, -6) incorrectly identical in each iteration?
• Re #1: Try adding the desired constraints: x2 != x3 || y2 != y3, x2 != x4 || y2 != y4, x3 != x4 || y3 != y4 – Michael E2 May 30 '19 at 16:58
• These constraints do not appear to stop the duplicates. Also, it increases timing by almost an order of magnitude. I'll also eventually need to use more than 10 points. Is there another way to do this? – Friasco Jun 2 '19 at 10:30
• I don't have a suggestion at this point. Related: mathematica.stackexchange.com/questions/22359/nminimize-usage/… and mathematica.stackexchange.com/questions/59706/… --They suggest that an occasional violation of the constraints may be unavoidable, but the final answer should satisfy the constraints. Why do you want to avoid the duplications? Perhaps there's another way to approach the problem. – Michael E2 Jun 2 '19 at 21:59
• Thanks for digging up the references. I'll try allowing it to generate duplicates and then kill them off later. It's important because the duplicates turn the n-point problem into a <n-point case and the optimization then often chooses these cases as the (wrong) solution. – Friasco Jun 3 '19 at 14:36
• I added an answer that should prevent the solution from degenerating to lower dimensions. – Michael E2 Jun 3 '19 at 18:29

## 2 Answers

One can enforce solutions with distinct points (or other custom constraint) by adding one's own penalty to the objective function. Here we subtract 1000 times the number of duplicates. The factor 1000 should be chosen large enough so that the inputs to be excluded are not local maxima. Some knowledge of the objective function is necessary to do this. I did a little refactoring, but none of it is important. It makes changing the number of points a little easier.

n = 4;
lines = n (n - 1)/2;
optimalelements = n - 1;
gridpoints = 100;
grid = Round[Sqrt[gridpoints]];

ClearAll[fitness, iFitness];
fitness[xy__Integer] :=  (* sets up main objective function iFitness[] *)
iFitness[Join[{{0, 0}}, Partition[{xy}, 2]]];
iFitness[points_] := Module[{linepoints, d, c, penalty, p},
linepoints = Subsets[points, {2}];
d = EuclideanDistance @@@ linepoints;
c = Sort[Tally[d][[All, 2]]];
penalty = Total[Abs[Differences[c]]] - Length[c];
If[c == Range[optimalelements], penalty = penalty + 1];
p = penalty;
Sow[{Sequence @@ points, c, p}];  (* I moved the EvaluationMonitor here *)
p - 1000 (Length[points] - Length@DeleteDuplicates[points]) (* extra penalty *)
];
vars = {x2, y2, x3, y3, x4, y4};
{sol, pts} = Reap[
NMaximize[{fitness @@ vars, And @@ Thread[0 <= vars < grid]},
vars \[Element] Integers,
Method -> {"SimulatedAnnealing", "SearchPoints" -> 1,
"PerturbationScale" -> 1,
"RandomSeed" -> 1}]]; // AbsoluteTiming


The list of variables vars can be generated automatically like this:

vars = Flatten@Table[{"x" <> ToString[i], "y" <> ToString[i]}, {i, 2, n}];
Clear @@ Flatten@vars;   (* optional safety measure *)
vars = ToExpression@Evaluate@vars
(*  {x2, y2, x3, y3, x4, y4}  *)

• Great, thank you, and the vars code will be most useful in scaling up to more points. Initially it is only necessary to show that one solution exists for any n, one can find the other ones later. Is it possible to stop the optimization when the first solution is found? i.e. Set "SearchPoints" -> \[Infinity] and in the Sow, if p == scoreforasolution, stop looking for more maxima (perhaps using Catch/Throw?) – Friasco Jun 4 '19 at 10:58
• @Friasco You could put something like If[p == scoreforasolution, Return[{Sequence @@ points, c, p}, NMaximize]]; in iFitness. Using more "SearchPoints" won't necessarily make it faster. There is probably a balance point, but someone would have to understand your problem pretty well to have a heuristic to suggest. – Michael E2 Jun 4 '19 at 17:47

Your code for the fitness function is extremely unorthodox, and I don't think it works the way you intended. I would instead use a Module:

fitness[x2_?NumberQ, y2_, x3_, y3_, x4_, y4_] := Module[
{points,linepoints,d,penalty},

points = {{0,0},{x2,y2},{x3,y3},{x4,y4}};
linepoints = Subsets[points,{2}];
d = EuclideanDistance@@@linepoints;
c = Sort[Tally[d][[All,2]]];
penalty = Total[Abs[Differences[c]]]-Length[c];
If[c==Range[optimalelements],penalty=penalty+1];
p=penalty
]


I don't modularize the c and p variables because you want to Sow them during evaluation. Then, your NMaximize call uses {{x2, 0, 10}, ..} syntax in the second argument. I don't think this does what you intended, as the bounds only suggest a starting point, and don't restrict the range of the variables. Also, using a domain as a 3rd argument isn't supported (at least in M12). So, I would use the following NMaximize call instead:

{sol,pts} = Reap[
NMaximize[
{fitness[x2,y2,x3,y3,x4,y4],
0<=x2<=grid,0<=y2<=grid,0<=x3<=grid,0<=y3<=grid,0<=x4<=grid,0<=y4<=grid},
{x2,y2,x3,y3,x4,y4} ∈ Integers,
Method->{"SimulatedAnnealing","SearchPoints"->1,"PerturbationScale"->1,"RandomSeed"->1},
EvaluationMonitor:>Sow[{{x2,y2},{x3,y3},{x4,y4},c,p}]
]
]//AbsoluteTiming


{0.060854, {{-2., {x2 -> 0, y2 -> 10, x3 -> 8, y3 -> 4, x4 -> 6, y4 -> 8}}, {{{{0, 4}, {0, 4}, {2, 6}, {1, 1, 2, 2}, -3}, {{0, 6}, {0, 2}, {4, 5}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 8}, {0, 6}, {4, 2}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 7}, {1, 6}, {4, 3}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 8}, {0, 5}, {5, 0}, {1, 1, 1, 1, 2}, -4}, {{5, 10}, {2, 8}, {4, 1}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 9}, {0, 5}, {3, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 6}, {0, 3}, {2, 0}, {1, 1, 1, 1, 2}, -4}, {{0, 8}, {0, 3}, {0, 0}, {1, 1, 2, 2}, -3}, {{4, 10}, {0, 2}, {0, 0}, {1, 1, 2, 2}, -3}, {{6, 10}, {2, 4}, {0, 0}, {1, 1, 2, 2}, -3}, {{6, 10}, {1, 1}, {4, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 10}, {0, 4}, {3, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 7}, {2, 7}, {0, 0}, {1, 1, 2, 2}, -3}, {{6, 10}, {2, 9}, {0, 0}, {1, 1, 2, 2}, -3}, {{2, 9}, {0, 8}, {2, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{7, 10}, {0, 10}, {1, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{6, 8}, {1, 9}, {0, 0}, {1, 1, 2, 2}, -3}, {{4, 10}, {0, 7}, {2, 1}, {1, 1, 1, 1, 1, 1}, -6}, {{3, 7}, {1, 8}, {3, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 8}, {1, 9}, {3, 1}, {1, 1, 1, 1, 2}, -4}, {{4, 8}, {0, 9}, {1, 1}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 10}, {2, 10}, {6, 1}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 10}, {4, 8}, {8, 0}, {1, 1, 1, 1, 2}, -4}, {{0, 10}, {6, 9}, {7, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 9}, {5, 9}, {7, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 10}, {3, 10}, {4, 0}, {1, 1, 2, 2}, -3}, {{3, 10}, {4, 10}, {5, 0}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 10}, {4, 10}, {6, 2}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 10}, {5, 7}, {6, 4}, {1, 1, 1, 1, 1, 1}, -6}, {{3, 8}, {7, 6}, {9, 6}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 10}, {8, 4}, {6, 8}, {1, 1, 1, 3}, -2}, {{0, 8}, {8, 3}, {7, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 7}, {8, 6}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 7}, {9, 8}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 5}, {10, 8}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 7}, {8, 9}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 8}, {6, 9}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 6}, {4, 8}, {10, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 7}, {5, 9}, {9, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 5}, {5, 8}, {10, 10}, {1, 1, 2, 2}, -3}, {{3, 4}, {6, 8}, {10, 10}, {1, 1, 1, 1, 2}, -4}, {{3, 5}, {5, 9}, {9, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 4}, {6, 8}, {8, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 4}, {7, 7}, {9, 10}, {1, 1, 1, 1, 2}, -4}, {{4, 5}, {6, 8}, {10, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{4, 2}, {5, 6}, {9, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{7, 4}, {4, 5}, {9, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{7, 3}, {5, 4}, {9, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{6, 4}, {4, 3}, {8, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{8, 3}, {2, 4}, {8, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{7, 4}, {2, 4}, {7, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{6, 3}, {1, 4}, {5, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 3}, {1, 5}, {6, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 3}, {0, 4}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{6, 3}, {0, 5}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{6, 3}, {0, 6}, {4, 10}, {1, 1, 1, 1, 2}, -4}, {{6, 4}, {0, 6}, {4, 10}, {1, 1, 1, 1, 2}, -4}, {{8, 5}, {0, 5}, {6, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{5, 4}, {0, 7}, {5, 10}, {1, 1, 1, 1, 2}, -4}, {{5, 3}, {0, 7}, {5, 10}, {1, 1, 2, 2}, -3}, {{6, 4}, {0, 7}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{3, 4}, {0, 6}, {5, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{3, 2}, {1, 6}, {5, 9}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 2}, {1, 6}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 3}, {1, 6}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 2}, {0, 6}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{2, 1}, {0, 7}, {2, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 0}, {1, 6}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 1}, {0, 6}, {3, 10}, {1, 1, 1, 1, 2}, -4}, {{1, 1}, {1, 7}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 1}, {0, 6}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 0}, {0, 6}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 1}, {0, 6}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 9}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 8}, {3, 9}, {1, 1, 2, 2}, -3}, {{0, 0}, {1, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {1, 8}, {4, 10}, {1, 1, 2, 2}, -3}, {{1, 1}, {1, 7}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {1, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {1, 8}, {2, 10}, {1, 1, 2, 2}, -3}, {{0, 1}, {1, 9}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {1, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {1, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{1, 0}, {1, 8}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{1, 0}, {1, 8}, {2, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{1, 0}, {0, 7}, {4, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 1}, {0, 7}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 8}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 1}, {0, 8}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 7}, {4, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 1}, {0, 7}, {3, 10}, {1, 1, 1, 1, 1, 1}, -6}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}, {{0, 0}, {0, 7}, {3, 10}, {1, 1, 2, 2}, -3}}}}}