# Find intersection points with V10

My question departs from this one: Find intersection of pairs of straight lines

but now I want to find the points with the new V10 Region-functions.

lines =
{{Line[{{243.8, 77.}, {467.4, 12.}}], Line[{{356.8, 32.}, {363.2, 120.}}]},
{Line[{{291.8, 130.}, {476., 210.5}}], Line[{{346., 245.}, {393.8, 158.}}]},
{Line[{{103.2, 327.}, {245.2, 110.5}}], Line[{{163.8, 211.5}, {230.2, 250.}}]},
{Line[{{47.4, 343.}, {87.4, 108.5}}], Line[{{54.6, 225.}, {139.6, 220.}}]},
{Line[{{371., 506.5}, {384.6, 277.}}], Line[{{366., 394.5}, {451.8, 372.}}]},
{Line[{{264.6, 525.5}, {353.8, 294.5}}], Line[{{241., 398.}, {321., 411.5}}]},
{Line[{{113.2, 484.5}, {296., 304.5}}], Line[{{163.2, 347.}, {213.2, 406.5}}]},
{Line[{{459.6, 604.5}, {320.2, 466.5}}], Line[{{332.4, 596.5}, {402.4, 528.5}}]},
{Line[{{288.2, 630.5}, {199.6, 446.5}}], Line[{{176., 585.5}, {256., 530.5}}]},
{Line[{{138.8, 615.5}, {81.8, 410.}}], Line[{{38.2, 553.}, {122.4, 507.}}]},
{Line[{{232.4, 795.}, {461.8, 727.}}], Line[{{345.2, 774.5}, {345.2, 688.}}]},
{Line[{{27.4, 671.5}, {206.8, 763.5}}], Line[{{104.6, 728.}, {161.8, 647.}}]}};


This function finds the points but is horribly slow:

(points = Point /@ RegionCentroid /@ DiscretizeRegion /@ RegionUnion @@@ lines); //
Timing // First


1.591210

Graphics[{lines, {Red, [email protected], points}}, Frame -> True]


The next function finds the same points in less time (0.28 seconds), but is ugly and probably not general enough.

points = Cases[Show[DiscretizeRegion /@ RegionIntersection @@@ lines],
{a_Real, b_Real} :> Point[{a, b}], Infinity];


I hope somebody can suggest a fast and terse V10-method to find intersection points in 2 dimensions

• Is something like Point[{x, y}] /. Solve[{x, y} \[Element] #, {x, y}] & /@ RegionIntersection @@@ lines the sort of thing you want? You have lots of requirements, fast, terse, V10, which might not be satisfied simultaneously. But this one is about four times faster than your last one. Jul 31, 2014 at 12:24
• @Michael E2 Yes, that's what I was hoping for. Why don't you put it as an answer?
– eldo
Jul 31, 2014 at 12:38

If you know all pairs of line segments intersect, then the following is about four times faster than using DiscretizeRegion:

points = (Point[{x, y}] /. Solve[{x, y} ∈ #, {x, y}] & /@ RegionIntersection @@@ lines); //
AbsoluteTiming // First
Graphics[{lines, {Red, [email protected], points}}, Frame -> True]
(*
0.040551
*)


Bug?

[Tested on beta V10 -- my institution is a little slow updating its software distribution.]

I feel this should work, but Solve emits error messages when RegionIntersection is used will Apply. Solve works perfectly on an explicit individual region intersection.

SeedRandom[1];
lines2 = Map[Line, RandomReal[{500, 800}, {5, 2, 2, 2}], {2}];
points = Point[p] /. Flatten[Solve[p ∈ #, p] & /@ RegionIntersection @@@ Join[lines2, lines], 1];
Graphics[{lines,
MapIndexed[{Hue[i/First[#2]], #1} &, lines2], {Red, [email protected], points}},
Frame -> True]


Solve::elemc: Unable to resolve the domain or region membership condition p ∈ RegionIntersection[Line[{{618.802,710.142},{563.548,724.597}}],Line[{{626.855,574.248},{<<17>>,<<18>>}}]]. >>

...rest omitted...

The first message corresponds to the second pair:

lines2[[2]]
(*
{Line[{{618.802, 710.142}, {563.548, 724.597}}],
Line[{{626.855, 574.248}, {793.152, 747.549}}]}
*)


This gives the same error message:

Solve[p ∈ RegionIntersection @@ lines2[[2]], p]


However, explicit code works fine:

Solve[p ∈
RegionIntersection[
Line[{{618.8018244645762, 710.1421345826734},
{563.5477937162382, 724.5970644448844}}],
Line[{{626.8551948019847, 574.2484342592024},
{793.151528522233, 747.5488818465542}}]],
p]
(*
{}
*)


A workaround is to use DeleteCases (and Quiet if desired):

points = Point[p] /.
Flatten[DeleteCases[
Solve[p ∈ #, p] & /@
RegionIntersection @@@ Join[lines2, lines], _Solve], 1];
Graphics[{lines,
MapIndexed[{Hue[First[#2]/5], #1} &, lines2], {Red, [email protected],
points}}, Frame -> True]


• am too slow...mine essentially a variant of yours...+1 Jul 31, 2014 at 13:29
• I'd be will to post the second part as a separate is-this-a-bug question, if someone can confirm it seems like misbehavior and that it can be reproduced. Jul 31, 2014 at 13:43
• @ubpdqn Thanks! I thought I was slow because the second part was annoying me no end. Jul 31, 2014 at 13:45
• very interesting update...will test my approach when get chance...with this full set of lines wrote my code to remove line segments which do not intersect but have not tested larger set...off to bed now...perhaps tomorrow Jul 31, 2014 at 13:48

With Mathematica 10 you can use:

ps[{l1_, l2_}] := Solve[p ∈ l1 ∧ p ∈ l2, p]


and then

points2 = Point[p] /. Map[ps[#] &, lines] // Flatten


For a large number of lines ParallelMap should give an additional speedup.

• You might want to use _?NumericQ for exact solutions (in case any exact input is present). Jul 31, 2014 at 12:37
• @Karsten Thanks for your answer. Even with these few lines ParallelMap more than doubles speed (0.038 seconds)
– eldo
Jul 31, 2014 at 12:59
ints[l1_, l2_] := Module[{rg, x, y, res},
rg = RegionIntersection[l1, l2];
res = {x, y} /. ToRules[Reduce[RegionMember[rg, {x, y}]]];
If[And @@ (NumericQ /@ res), res = res, res = Sequence[]]; res];
pair[dat_] := Subsets[dat, {2}]


As the lines have already been paired:

li01 = Graphics[{{Hue[RandomReal[]], Thick, #} & /@ lines, {Red,
PointSize[0.03], Point[ints @@@ lines]}}]


Or pretending you just had a bunch of lines:

li02 = Graphics[{{Thick, Blue, lines}, {Red, PointSize[0.04],
Point[ints @@@ pair[Flatten@lines]]}}]