Drawing pairs (planar graph, circle packing)

Question

Circle packing theorem states:

For every connected simple planar graph G there is a circle packing in the plane whose intersection graph is (isomorphic to) G.

How to use Mathematica to draw examples of such pairs (planar graph, circle packing)?

More specifically, the requirements of my question could be formulated as follows:

Generate a random connected simple planar graph with given number of nodes $N$ (this can be done lets say using adjacency matrix; in a way another parameter affecting resulting graph could be "density" of links). Draw its usual graphical representation containing nodes and links. Then draw its circle packing representation.

---

In other words, a solution could contain a slider for $N$, and a slider for the "density" (those two variables are parameters). Based on values of these sliders, get regular drawing of a random graph with given $N$ and density, and also drawing of its circle pack.

---

NOTE ON GENERATING RANDOM PLANAR GRAPHS:

To my knowledge, this is not a trivial problem, and there is no built-in Mathematica solution. Though, there are several questions and answers here. One of best is Create triangular mesh from random list of points :

SeedRandom[1];
pts = RandomReal[{0, 12}, {100, 2}];
Needs["ComputationalGeometry`"];
dt = DelaunayTriangulation[pts];
dt // Column
toPairs[{m_, ns_List}] := Map[{m, #} &, ns];
edges = Flatten[Map[toPairs, dt], 1];
Graphics[GraphicsComplex[pts, {Line[edges], 
  Red, PointSize[Large], Point[pts]}]]

---

NOTE ON GENERATING CIRCLE PACK CORRESPONDENT TO A PLANAR GRAPH:

See: Collins, Stephenson: A circle packing algorithm.

---

Visual examples:

---

---

---

---

This question is not an assignment. It is an intellectual experiment, related to aestetics and structure. I am not requesting code. I am seeking ideas, approaches, and solutions.

Answer 1

I bought 'Introduction to Circle Packing' by Kenneth Stephenson a few years ago, and have been trying to understand it since. I found theFindMaximumanswer by @SimonWoods very good, but larger graphs slow that algorithm down, as he said.

@VividD said "I am not requesting code" and "I am seeking ideas, approaches, and solutions". I finally have the code, but there is a lot of it. Therefore, here is a description of one approach based on my interpretation of Stephenson.

The key is to break the problem into two halves.

First: Find radii without placing any circles.

Find the radii of all circles such that the combinatorics of the graph are satisfied. That is, find radii such that the distance between adjacent nodes i and j is radii[[i]]+radii[[j]], for all i and j, i $\ne$ j.

Stephenson (page 244) describes the uniform neighbour algorithm to iteratively find these radii. For any node $V$ in the graph, adjust its current radius by multiplying it by $b (1-d) /(d (1-b))$, where $b=\sin(t/(2k))$, $d=\sin(\pi/k)$, $t$ is the angle sum, and $k$ is the number of neighbouring nodes surrounding $V$. The angle sum (pages 15, 57) is the total angle subtended at the central node $V$ by adjacent pairs of surrounding circles labelled with their current, as yet approximate, radii. An angle sum greater than $2\pi$ means the central circle radius is too small, causing the surrounding circles to overlap. An angle sum less than $2\pi$ means the central circle radius is too large, producing gaps between the surrounding circles. An angle sum of $2\pi$ is just right, the $k$ surrounding circles are pair-wise tangent.

The uniform neighbour algorithm of Stephenson cycles through each interior node, adjusting the radius corresponding to that node. Multiple cycles through the interior nodes enables convergence to a stable set of radii.

Second: Find the circle centres corresponding to their radii.

Once the radii are found, place any two adjacent circles, a and b, such that their centres are separated by the sum of their radii, that is, they are tangent. Then place a (positively oriented) triangular face by finding the position of the centre of the third circle c which makes a mutually tangent triple with these two placed circles. This process generates two new pairs of placed circles, (c,b) and (a,c), leading to placements of new third circles, triangular face by face.

Following this outline, plus some compiling and memoization, I was able to construct circle packings of 1500 circles in about 125 s. Finding the radii is by far the most time-consuming part.

Rendering is the easy part:

Answer 2

Two versions of a solution to this question will be presented. Both are applications of built-in Mathematica optimization tools. The key point of both version is the choice of optimization goal and constraints (or, in another words, model formulation).

Interested reader may find article Configuration Analysis and Design by Using Optimization Tools in Mathematica by FRANK J. KAMPAS and JÁNOS D. PINTÉR useful. It is not about the problem from this question, but deals with solving geometric problems using Mathematica's optimization tools. There are also several related questions on thus site.

First version

This will surely be hopeless for large graphs, but here's a crude approach using NMinimize. The aim is to minimize the difference between the inter-vertex distance and the sum of their radii, and also the offset distance between the vertex positions in the original and final graph.

Let's say heptahedral-5 graph is chosen as the original graph:

g = GraphData[{"Heptahedral", 5}]

This code will produce its circle packing, using NMinimize in appropriate setting:

n = VertexCount[g];    
pts = Table[PropertyValue[{g, i}, VertexCoordinates], {i, n}];    
rads = Array[r, n];
xs = Array[x, n];
ys = Array[y, n];

f1 = Total[(r[#1] + r[#2] - EuclideanDistance[{x[#1], y[#1]}, {x[#2], y[#2]}])^2 & @@@ EdgeList[g]];
f2 = Total[MapThread[({#1, #2} - #3)^2 &, {xs, ys, pts}], -1];

rules = Last@NMinimize[f1 + 0.1 f2, Join[rads, xs, ys]];

circles = MapThread[Circle[{#1, #2}, #3] &, {xs, ys, rads} /. rules];
g2 = Fold[SetProperty[{##}, VertexCoordinates -> {x[#2], y[#2]} /. rules] &, g, Range[n]];

Show[Graphics[circles], g2]

Improved version

This is an improved version which seems to cope with larger graphs. This is still based on black box optimization but with some changes:

• Use FindMaximum instead of NMaximize to find a local maximum.

• Circles which should be touching (i.e. graph edges) are put in as constraints, i.e. sum of radii == distance between centres.

• All other pairs of vertices (non-edges) are placed as far apart as possible - the sum of these distances is what we maximize.

Also, in this example, I've used a Delaunay triangulation (as suggested in the question) to create a simple random planar graph (however, using an undocumented DelaunayMesh function.

Other steps are also marked and described by comments in the code:

Graphics`Region`RegionInit[];
dist[{a_, b_}] := Sqrt[(a - b).(a - b)];
n = 30;

(* get a random mesh *)
pts = RandomReal[1, {n, 2}];
m = DelaunayMesh[pts]["MeshObject"];

(* get edges and non-edges *)
edges = Sort /@ m["Edges"];
nonedges = Complement[Subsets[Range[n], {2}], edges];

(* variable lists *)
rads = Array[r, n];
pos = Array[{x[#], y[#]} &, n];

(* set up constraints for edges *)
cons = (Total@rads[[#]] == dist@pos[[#]]) & /@ edges;
cons = Join[cons, Thread[rads > 0]];

(* define function to maximize *)
f = Total[dist@pos[[#]] & /@ nonedges];

(* variable initializations *)
vars = Join[rads, Flatten[{pos, pts}, {2, 3}]];

(* do the maximization *)
{posv, radsv} = {pos, rads} /. Last@FindMaximum[{f, cons}, vars];

(* result *)
Graphics[{MapThread[Circle, {posv, radsv}],
  Opacity[0.5, Green], Line[posv[[#]] & /@ edges]}]