Finding the distance to a Hexagon

596 Views Asked by At

Consider the regular hexagon

enter image description here

of side size equal to $1$, in the $xy$-plane and centered at the origin $(0,0)$.

For a given point $P=(x,y,z)$, how do I find the minimum distance to the boundary of the hexagon? It seems I might need a branching function, any ideas are appreciated.

2

There are 2 best solutions below

3
On BEST ANSWER

Iterative solution

Solution valid for each regular polygon with $n$-sides long $L$, of center $(0,\,0)$ and with a vertex lying on the positive $x$-axis, but still easily translatable and rotatable as desired.

In Wolfram Mathematica, defining:

{n, L, dx, xP, yP} = {6, 1, 0, 0, 0};
R = L / (2 Sin[Pi/n]);
poly = Table[{R Cos[2 k Pi/n] (1 + dx), R Sin[2 k Pi/n] / (1 + dx)}, {k, n}];

and executing the following simple algorithm:

{perimeter, area, flag} = {0, 0, 0};
For[k = 1, k <= n, k++,
    {xA, yA} = poly[[k]];
    {xB, yB} = If[k < n, poly[[k + 1]], poly[[1]]];
    dAB = Sqrt[(xA - xB)^2 + (yA - yB)^2];
    perimeter = perimeter + dAB;
    area = area + (xA - xB) (yA + yB) / 2;
    xQ = ((xA - xB)^2 xP + (yA - yB) (xA (yP - yB) - xB (yP - yA))) / dAB^2;
    yQ = ((yA - yB)^2 yP + (xA - xB) ((yA - yB) xP + xA yB - xB yA)) / dAB^2;
    dAQ = Sqrt[(xA - xQ)^2 + (yA - yQ)^2];
    dBQ = Sqrt[(xB - xQ)^2 + (yB - yQ)^2];
    If[Abs[N[dAB - dAQ - dBQ]] > 10^-10, 
       {xQ, yQ} = If[N[dAQ - dBQ] < 0, {xA, yA}, {xB, yB}]
      ];
    dPQ = Sqrt[(xP - xQ)^2 + (yP - yQ)^2];
    If[flag == 0,
       flag = 1;
       dmin = dPQ;
       points = {{xQ, yQ}},
       If[Abs[N[dPQ - dmin]] < 10^-10,
          points = Union[points, {{xQ, yQ}}],
          If[N[dPQ - dmin] < 0,
             dmin = dPQ;
             points = {{xQ, yQ}}
            ]
         ]
      ]
   ];

it's soon achieved as desired:

Graphics[{Green, RegionBoundary[Polygon[poly]],
          Red, PointSize[Large], Point[{xP, yP}],
          Blue, PointSize[Large], Point[points]}]
Print["points = ", FullSimplify[points]]
Print["dmin = ", FullSimplify[dmin]]
Print["2p = ", FullSimplify[n L]]
Print["2p' = ", FullSimplify[perimeter]]
Print["A = ", FullSimplify[n R^2 Sin[2 Pi/n] / 2]]
Print["A' = ", FullSimplify[area]]    

enter image description here

If then, again by way of example, we re-run the same code with:

{n, L, dx, xP, yP} = {6, 1, 0, 2, 1};

we obtain:

enter image description here

while, if we impose:

{n, L, dx, xP, yP} = {6, 1, 1/10, 2, 1};

we obtain:

enter image description here

Of course, if the numerical results are sufficient, in the output phase replace FullSimplify[] with N[]; the code will be executed more quickly.


Analytical solution

Solution valid for each regular hexagon with long sides $L$, of center $(0,\,0)$ and with a vertex lying on the positive $x$-axis, but still easily translatable and rotatable as desired.

In Wolfram Mathematica, defining:

{L, dx, xP, yP} = {1, 0, 0, 0};
poly = Table[{L Cos[2 k Pi/6] (1 + dx), L Sin[2 k Pi/6] / (1 + dx)}, {k, 6}];

and executing the following simple algorithm:

f[x_, y_] := (x - xP)^2 + (y - yP)^2
g[x_, y_] := RealAbs[2 y] + RealAbs[y + Sqrt[3] x] + RealAbs[y - Sqrt[3] x] - 2 Sqrt[3] L
h[x_, y_, z_] := f[x, y] + z g[x / (1 + dx), y (1 + dx)]
grad = D[h[x, y, z], {{x, y, z}, 1}];
sol = Solve[grad == {0, 0, 0}, Reals];

flag = 0;
For[k = 1, k <= Length[sol], k++,
    d = Sqrt[f[x, y] /. sol[[k]]];
    point = {sol[[k, 1, 2]], sol[[k, 2, 2]]};
    If[flag == 0,
       flag = 1;
       dmin = d;
       points = {point},
       If[Abs[N[d - dmin]] < 10^-10,
          points = Union[points, {point}],
          If[N[d - dmin] < 0,
             dmin = d;
             points = {point}
            ]
         ]
      ]
   ];

perimeter = Integrate[1, {x, y} ∈ ImplicitRegion[g[x/(1 + dx), y(1 + dx)] == 0, {x, y}]];
area = Integrate[1, {x, y} ∈ ImplicitRegion[g[x/(1 + dx), y(1 + dx)] <= 0, {x, y}]];

it's soon achieved as desired:

Graphics[{Green, RegionBoundary[Polygon[poly]], 
          Red, PointSize[Large], Point[{xP, yP}], 
          Blue, PointSize[Large], Point[points]}]
Print["points = ", FullSimplify[points]]
Print["dmin = ", FullSimplify[dmin]]
Print["2p = ", FullSimplify[6 L]]
Print["2p' = ", FullSimplify[perimeter]]
Print["A = ", FullSimplify[3 Sqrt[3] L^2 / 2]]
Print["A' = ", FullSimplify[area]] 

Of course, if the numerical results are sufficient, in the output phase replace FullSimplify[] with N[], while in the body of the algorithm replace Solve[] and Integrate[] respectively with NSolve[] and NIntegrate[]; the code will be executed more quickly.

0
On

Let $(x_i, y_i),\>i =1...6$ be the six vertexes. Assume that $(x,y)$ lies outside the hexagon. (The inside case is straightforward by comparing distances to the six sides.) below are the systematic steps to find the minimum distance.

1) Find the vertex $k$ with the minimum $d_k^2=(x-x_k)^2+(y-y_k)^2$.

2) Construct vectors $Q(x-x_k, y-y_k)$, $A(x_{k+1} -x_k,y_{k+1}- y_k)$, $B(x_{k-1}-x_k, y_{k-1}- y_k)$ and evaluate the dot products $a = Q\cdot A=d_k\cos\theta_a$ and $b = Q\cdot B= d_k\cos\theta_b$.

3) As seen in the diagram below, there are three cases for the minimum distance $d_m$ in the $xy$-plane, depending on the angle between $Q$ and the two sides .

Case 1: $a < 0$ and $b < 0$, the minimum distance is to the vertex itself, $d_m=d_k$.

Case 2: $a >0$, the minimum distance is to the side $A$, $d_m = d_k\sin\theta_a=\sqrt{d_k^2-a^2}$.

Case 3: $b> 0$, the minimum distance is to the side $B$, $d_m = d_k\sin\theta_b= \sqrt{d_k^2-b^2} $.

4) The 3D minimum distance is then $\sqrt{d_m^2+z^2}$.

enter image description here