# Kobon Triangles: number of nonoverlapping Δs from $n$ lines

Posted 3 years ago
3725 Views
|
|
8 Total Likes
|

This comes from an old post over on the SE about Kobon Triangles, but I've added some new functionality and thought I might share it here.

I enjoyed working on this, I learned how to construct a MeshRegion from a set of points, and to find polygons by using FindCycle. First I will give the code and then explain,

kobonTriangle[k_] :=
Module[{r0, r1, r2, pts, ilns, lines, edges, vertices, triangles, mesh},
r0 := RandomReal[{-1, 1}];
r1 := RandomReal[{-1, 0}];
r2 := RandomReal[{0, 1}];
pts = Transpose[
{Array[{r0, r1} &, k - 1],
Array[{r0, r2} &, k - 1]
}];
ilns = InfiniteLine /@ pts~Join~{{{0, 0}, {1, 0}}};
lines = Flatten[
Partition[Sort@#, 2, 1] & /@ Table[
Flatten[List @@@ (RegionIntersection[
ilns[[n]], #] & /@ Delete[ilns, n]), 1],
{n, Length@ilns}], 1];

vertices = Flatten[lines, 1] // DeleteDuplicates;
edges = lines /. MapIndexed[#1 -> First@#2 &, vertices];
triangles = FindCycle[Graph[#1 <-> #2 & @@@ edges], {3}, All];
mesh = MeshRegion[
vertices, {Line /@ edges,
triangles /. {a_ <-> b_, b_ <-> c_, c_ <-> a_} :>
Polygon[{a, b, c, a}]}];
(* Uncomment this section to have the number of triangles reported on the image *)
(* Labeled[
mesh,
Row[{"Number of lines = ", k, ", Number of Triangles = ",
Length@triangles}]] *)
mesh
];


Here are a couple of examples,

{kobonTriangle[5], kobonTriangle[8]}


<img src="http://i.stack.imgur.com/MsaOA.png" alt="enter image description here"/>

In any iteration, chances are you won't find the optimal solution. For example, for 5 and 8 lines, there are solutions with 5 and 15 triangles, respectively, rather than the 3 and 9. But if you run the code enough times, you can often find a near-optimal solution. There are most definitely better algorithms to search for these, but I like this one. I let it run for an hour and got these results:

<img src="http://i.stack.imgur.com/TJmmc.png" alt="enter image description here"/>

## How it works

I was inspired by Trevor Simonton's javascript code here. The idea is to generate k random lines that intersect so as to get a decent number of triangles. To that end, we start with one line that is oriented horizontally, and then generate k-1 lines that cross this line.

Here is the code to do this,

r0 := RandomReal[{-1, 1}];
r1 := RandomReal[{-1, 0}];
r2 := RandomReal[{0, 1}];
pts = Transpose[
{Array[{r0, r1} &, k - 1],
Array[{r0, r2} &, k - 1]
}];
ilns = InfiniteLine /@ pts~Join~{{{0, 0}, {1, 0}}};


You can see the lines via,

Graphics[ilns]


<img src="http://i.stack.imgur.com/IbQHV.png" alt="enter image description here"/>

We need to zoom out to see all the intersections

Graphics[ilns, PlotRange -> {{-2, 2.0}, {-2, 2.0}}]


<img src="http://i.stack.imgur.com/tietE.png" alt="enter image description here"/>

Now I would like to cut off the lines after the intersection points to create a closed shape. First I will use RegionIntersection to find all the intersection points. Then I create line segments between each intersection point, but first I sort the intersection points to make sure that we don't have any overlapping line segments.

lines = Flatten[
Partition[Sort@#, 2, 1] & /@ Table[
Flatten[List @@@ (RegionIntersection[
ilns[[n]], #] & /@ Delete[ilns, n]), 1],
{n, Length@ilns}], 1];
vertices = Flatten[lines, 1] // DeleteDuplicates;
Graphics[{Line /@ lines, {Red, PointSize[Medium], Point /@ vertices}}]


<img src="http://i.stack.imgur.com/iwGXI.png" alt="enter image description here"/>

So we have our basic shape, but how to find the triangles, and only the non-overlapping triangles? By making a Graph that is isomorphic to the shape above, we can take advantage of the Graph functions in Mathematica

edges = lines /. MapIndexed[#1 -> First@#2 &, ipts]
Graph[edges, VertexLabels -> "Name"]
(* {{1, 2}, {2, 3}, {3, 4}, {5, 3}, {3, 6}, {6, 7}, {8, 9}, {9,
5}, {5, 4}, {9, 10}, {10, 1}, {1, 7}, {8, 10}, {10, 2}, {2, 6}} *)


<img src="http://i.stack.imgur.com/eP38F.png" alt="enter image description here"/>

Now we can find the triangles easily enough, and only non-overlapping triangles will be found because we've cut the lines into non-overlapping segments already.

triangles = FindCycle[Graph[#1 <-> #2 & @@@ edges], {3}, All]
Length@triangles
(* {{8 <-> 9, 9 <-> 10, 10 <-> 8}, {2 <-> 3, 3 <-> 6,
6 <-> 2}, {1 <-> 10, 10 <-> 2, 2 <-> 1}, {3 <-> 4, 4 <-> 5, 5 <-> 3}} *)
(* 4 *)


Now we just wrap it all up into a MeshRegion for display purposes,

Labeled[
MeshRegion[
vertices, {Line /@ edges,
triangles /. {a_ <-> b_, b_ <-> c_, c_ <-> a_} :>
Polygon[{a, b, c, a}]}],
Row[{"Number of lines = ", k, ", Number of Triangles = ",
Length@triangles}]]


<img src="http://i.stack.imgur.com/Jd4sb.png" alt="enter image description here"/>

So this code is perhaps not efficient - I imagine that FindCycles and the routine to find the intersection both scale at or worse than $\mathcal{O}(n^2)$ but $n$ is small so that is no worry.

## One slight problem

For every time we get a decent shape like

kobonTriangle[14]


<img src="http://i.stack.imgur.com/eLTS8.png" alt="enter image description here"/>

we will get 5 that look like this,

<img src="http://i.stack.imgur.com/thIfx.png" alt="enter image description here"/>

where some of the lines are so long as to make the shape hard to view. The best plan I could come up with for discriminating against the latter shapes is by comparing the area of the triangles to the total area of a square encasing all the points. The shape with the best relative triangle area wins. Here is a routine that goes through 2000 9-line shapes and displays the best result alongside the current result,

evaluate[mesh_] := {(Area /@ MeshPrimitives[mesh, 2] //
Total)/(Times @@
Abs@*Subtract @@@ (MinMax /@
Transpose[MeshPrimitives[mesh, 0][[All, 1]]])),
Length@MeshPrimitives[mesh, 2]}
bestN = 0;
bestA = 0;
bestN = 0;
bestA = 0;

Dynamic[{bestN, bestTriangle, currentTriangle}]
With[{nt = 9},
Do[
currentTriangle = kobonTriangle[nt];
{area, ntriangles} = evaluate@currentTriangle;
If[ntriangles >= bestN,
bestN = ntriangles;
If[area > bestA,
bestA = area;
bestTriangle = currentTriangle];
];
, {2000}]
]