1
$\begingroup$

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?
$\endgroup$
  • $\begingroup$ Re #1: Try adding the desired constraints: x2 != x3 || y2 != y3, x2 != x4 || y2 != y4, x3 != x4 || y3 != y4 $\endgroup$ – Michael E2 May 30 at 16:58
  • $\begingroup$ 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? $\endgroup$ – Friasco Jun 2 at 10:30
  • 1
    $\begingroup$ 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. $\endgroup$ – Michael E2 Jun 2 at 21:59
  • 1
    $\begingroup$ 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. $\endgroup$ – Friasco Jun 3 at 14:36
  • 1
    $\begingroup$ I added an answer that should prevent the solution from degenerating to lower dimensions. $\endgroup$ – Michael E2 Jun 3 at 18:29
1
$\begingroup$

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}  *)
$\endgroup$
  • $\begingroup$ 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?) $\endgroup$ – Friasco Jun 4 at 10:58
  • 1
    $\begingroup$ @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. $\endgroup$ – Michael E2 Jun 4 at 17:47
2
$\begingroup$

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}}}}}

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.