Inconsistent behavior from RegionMeasure on the RegionIntersection between a Polygon and a Line

Question

Consider a Polygon with lines drawn between two of its vertices,

pgon = Get[
    "https://gist.githubusercontent.com/jasondbiggs/\
59cb7d4aa802bde68c9f6a5203ef698f/raw/\
8b3768422fdcc5c28d03cd1d7a5361c0a1d2d15a/gistfile1.txt"][[2]];
pts = pgon[[1]];
pair1 = {2, 4};
pair2 = {2, 7};
Graphics[{Red, 
  pgon, {Blue, Line[pts[[#]]], Green, PointSize[Large], 
     Point[pts[[#]]]} & /@ {pair1, pair2}}]

Mathematica graphics

Clearly one of the lines lies entirely within the polygon, while the other touches the polygon at its endpoints. We can ask for the intersections between these regions

intersections = 
 RegionIntersection[Line[pts[[#]]], pgon] & /@ {pair1, pair2}
(* {RegionIntersection[
  Line[{{7.52229, 9.31563}, {7.93484, 4.23025}}], 
  Polygon[{{9.23481, 8.12637}, {7.52229, 9.31563}, {8.39443, 
     4.79338}, {7.93484, 4.23025}, {6.00149, 1.75529}, {7.75656, 
     2.24465}, {9.438, 3.48553}}]], 
 RegionIntersection[Line[{{7.52229, 9.31563}, {9.438, 3.48553}}], 
  Polygon[{{9.23481, 8.12637}, {7.52229, 9.31563}, {8.39443, 
     4.79338}, {7.93484, 4.23025}, {6.00149, 1.75529}, {7.75656, 
     2.24465}, {9.438, 3.48553}}]]} *)

Ideally, the system would return either Line or Point objects, but it doesn't. We could then ask for the RegionMeasure of these intersections, to determine whether either of the lines is contained by the polygons,

RegionMeasure /@ intersections
(* {2., 6.13677} *)

The first result is the number of 0-dimensional points at which an intersection occurs, and the second result is the length of the line.

This seems inconsistent. In trying to define a function to determine whether the line lies within the polygon, I would like to get a zero for the case where the line has zero length inside.

I get no help from ArcLength, RegionDimension, or RegionEmbeddingDimension

ArcLength /@ intersections
RegionDimension /@ intersections
RegionEmbeddingDimension /@ intersections
(* {2., 6.13677} *)
(* {1, 1} *)
(* {2, 2} *)

Is this expected behavior?

In the case of RegionMeasure it is expected, it is mentioned under "possible issues" in the documentation. — C. E., May 27 '16 at 12:49
Couldn't you say with a very high probability that the region measure will not be exactly 2 if the line intersects the polygon? — C. E., May 27 '16 at 12:58
That was what I was doing before, but then I decided to just say that if the ArcLength of the Intersection was the same as the length of the line, then it is in the polygon. It still seems like a non-foolproof method, as there can be cases where the line has a length of exactly 2. — Jason B., May 27 '16 at 13:00
@C.E. - I could delete the question (since I didn't read the Possible Issues). I still think the system could respond better to this request, and it seems that ArcLength and RegionDimension are just wrong in this case. I had seen this method but it fails for the case above. I could rephrase the question to be about RegionDimension — Jason B., May 27 '16 at 13:07

score 1 · Accepted Answer · edited Apr 13 '17 at 12:55

This is very similar to this question. In that case you are looking at the 1-dimensional overlap between 2-dimensional regions and taking its measure. Here you are looking at the intersection between a 1D and a 2D region and taking its measure. Both of these cases can be considered "near singular", to borrow a phrase from one of the developers, and it is even mentioned as one of the "possible issues" of RegionMeasure.

But the result of ArcLength is definitely wrong, and the fact that it doesn't even give an error message means this is a bug. I can guess that the function that decides how the arc length should be calculated and the function that is called to calculate the arc length react differently to finite-precision representations of near-singular regions.

But we can use the same method as suggested by illian here to get around this issue. Essentially we need to convert the numbers in the Polygon and the Line to exact numbers via exact[p_]:=SetPrecision[p,∞] or exact[p_]:= Rationalize[p, 0] and get the ArcLength via

ArcLength /@ exact[intersections] // N
(* {0., 6.13677} *)

score 1 · Answer 2 · answered Feb 01 '18 at 06:27

Looks like there has been some progress on this issue. This is what I get in M11.2:

intersections = RegionIntersection[Line[pts[[#]]],pgon]&/@{pair1,pair2}

{Point[{{7.52229, 9.31563}, {7.93484, 4.23025}}], Line[{{7.52229, 9.31563}, {9.438, 3.48553}}]}

Inconsistent behavior from RegionMeasure on the RegionIntersection between a Polygon and a Line

2 Answers2

Linked