0

Let's define some simple functions

Clear["Global`*"];

x0 = 0; y0 = 0; 
x1 = 1/Sqrt[3]; y1 = 0;
x2 = -(x1/2); y2 = 1/2;
x3 = x2; y3 = -y2;

r0 = Sqrt[(x - x0)^2 + (y - y0)^2];
r1 = Sqrt[(x - x1)^2 + (y - y1)^2];
r2 = Sqrt[(x - x2)^2 + (y - y2)^2]; 
r3 = Sqrt[(x - x3)^2 + (y - y3)^2];

Ω = 1/(3*(1 + β*Sqrt[3]))*(β/r0 + 1/r1 + 1/r2 + 1/r3) + 1/2*(x^2 + y^2);

μ = 0.5;
β = 1/μ - 1;

Ωx = D[Ω, x];
Ωy = D[Ω, y];

Then the module FindRoots2D for finding the intersections of two curves

Options[FindRoots2D] = {WorkingPrecision -> 30, MaxRecursion -> 30};

FindRoots2D[funcs_, {x_, a_, b_}, {y_, c_, d_}, opts___] := 
Module[{fZero, seeds, signs, fy}, 
  fy = Compile[{x, y}, Evaluate[funcs[[2]]]];
  fZero = 
  Cases[Normal[
  ContourPlot[
  funcs[[1]] == 0, {x, a - (b - a)/97, b + (b - a)/103}, {y, 
   c - (d - c)/98, d + (d - c)/102}, 
  Evaluate[FilterRules[{opts}, Options[ContourPlot]]]]], 
  Line[z_] :> z, Infinity];
  seeds = Flatten[((signs = Sign[Apply[fy, #1, {1}]];
    #1[[
     1 + Flatten[
       Position[Rest[signs*RotateRight[signs]], -1]]]]) &) /@ 
  fZero, 1];
  If[seeds == {}, {}, 
  Select[Union[({x, y} /. 
     Check[FindRoot[{funcs[[1]], 
        funcs[[2]]}, {x, #1[[1]]}, {y, #1[[2]]}, 
       Evaluate[FilterRules[{opts}, Options[FindRoot]]]], {x -> 
        b + 999, y -> d + 999}] &) /@ seeds, 
  SameTest -> (Norm[#1 - #2] < 10^(-8) &)], 
  a <= #1[[1]] <= b && c <= #1[[2]] <= d &]]]

Then we apply this module to the first-order derivatives

pts = Chop[FindRoots2D[{Ωx, Ωy}, {x, -2, 2}, {y, -2, 2}, PlotPoints -> 100]];
nps = Length[pts];
Print["N = ", nps]

The final output is

cont = ContourPlot[{Ωx == 0, Ωy == 0}, {x, -2, 2}, {y, -2, 2}, 
        ContourShading -> False, ContourStyle -> {{Thick, Darker[Green]}, {Thick, Blue}}, 
        PlotPoints -> 200, PerformanceGoal :> "Speed", 
        Epilog -> {AbsolutePointSize[8], Red, Point@pts}, ImageSize -> 500]

enter image description here

We see that there are 9 intersections, the coordinates of which are contained in the list pts. Personally, I have no idea how Mathematica orders the 9 points inside the list.

My question is the following: how can I automatically order (sort) the points, according to the labels 1, 2, ..., 9? I mean, how can I create a new list pts2 with {{x_1, y_1}, {x_2, y_2}, {x_3, y_3}, ..., {x_9, y_9}}?

Any ideas?

Many thanks in advance!

gwr
  • 13,452
  • 2
  • 47
  • 78
Vaggelis_Z
  • 8,740
  • 6
  • 34
  • 79

1 Answers1

2

In answering this I am assuming we have the picture and the labels and would simply like to sort using approximate coordinates for the labels:

(* approx coordinates for the labels in numerical order *)
approxCoordinates = {
    {   0 , 0   },
    {   1 , 0   },
    {-0.5 , 0   },
    {   0 , 0.5 },
    {   0 ,-0.5 },
    {-0.5 , 1   },
    {-0.5 ,-1   },
    {-0.25, 0.25},
    {-0.25,-0.25}
}

We can now use this information to define a NearestFunction using Nearest:

nf = First @* Nearest[ approxCoordinates -> Range[9] ]

(* Version 9 *)
nfv9 = Nearest[ approxCoordinates -> Range[9] ]

We then apply this to the list of points:

pos = Map[ nf, pts ]

(* Version 9 *)
(* pos = Map[ nfv9, pts ] // Flatten *)

{3, 6, 7, 9, 8, 1, 5, 4, 2}

Now we simply have to use Ordering and use it for rearranging the list of points:

posOrdering = Ordering @ pos

{6, 9, 1, 8, 7, 2, 3, 5, 4}

pts[[ posOrdering ]]

{{0.289946, 0}, {0.984258, 0}, {-0.696828, 0}, {0.348414, 0.603471}, {0.348414, -0.603471}, {-0.492129, 0.852392}, {-0.492129, -0.852392}, {-0.144973, 0.2511}, {-0.144973, -0.2511}}

gwr
  • 13,452
  • 2
  • 47
  • 78