# Generating Doyle spiral painting

I recently came across an interesting paining by Nicola Sutcliffe:

This painting is actually related to Doyle spirals.

From author's website:

The central part of the picture shows the Doyle spiral cicle packing P = 2 Q = 12. Three spirals pass through each circle, dividing it into six segments, each of which is a different colour, organised so that no touching segments have the same colour.

There are Wolfram demos devoted to Doyle spirals:

Can Mathematica somehow generate the above painting (just the central part with circles)? Furthermore, can the code be written to generate a parametrized version of the painting so that one could produce and enjoy an unlimited number of pieces of art?

The solution does not even have to be an exact copy of the painting, but it should preserve its spirit.

The Wolfram demo contains the interesting code with many options that may serve as a good starting point.

Related Wolfram demonstrations:

By the way, the paining was part of the art exhibition at the Bridges (maths and arts) 2008 conference in Leeuwarden.

Interesting math notes related to Doyle spirals:

• It would be good if you could add sources. A link to where you got the image, a link to support the claim that "the painting is actually related to Doyle spirals", and a link to the "Wolfram demo" that you mention. Maybe also a link to a page about Doyle spirals if you have one you think is good, because why not. And the usual question, to which you are not new: What are your thoughts on this problem? Have you tried anything? – C. E. May 18 '17 at 23:01
• A small discovery: the case in the painting corresponds to parameter P=4, Q=9 in the second Wolfram demo. – Wjx May 27 '17 at 15:17
• @Wjx I really should pay attention to comments. After seeing that my spirals were off, I printed out the original and counted the circles only to discover what you already have found out. Did you see that she explicitly says that she used P=2 and Q=12 on her webpage? Very nice finding :) – halirutan Jun 2 '17 at 23:12

TL;DR Yes, this can be done!

If you read the article "Hexagonal circle packings and Doyle spirals" by Leys, you will see that for a choice p and q, we need to find the complex values A, B and r. For that purpose, we can steal this part from the demonstration you linked:

doyle[pi_, qi_] := Module[{p = pi, q = qi, s, t, r},
r[s_, t_, p_,
q_] := (s^2 + s^(2 p/q) -
2 s^((p + q)/q) Cos[(2 \[Pi] + p t - q t)/q])/(s + s^(p/q))^2;
{s, t} = {s, t} /. FindRoot[
{r[s, t, 0, 1] - r[s, t, p, q] == 0,
r[s, t, 0, 1] - r[s^(p/q), (p*t + 2. Pi)/q, 0, 1] == 0},
{{s, 1.5}, {t, 0}}];
{s*Exp[I*t], s^(p/q)*Exp[I (p*t + 2*Pi)/q], Sqrt[r[s, t, 0, 1]]}
]

{a, b, r} = doyle[2, 12]


Now we have the centers for both additional complex circles. Knowing a bit of complex analysis, one understands that for creating all packed circles, we only need to iteratively multiply by e.g. a. So we could write down a function that does this multiplication. I'm keeping complex values until the final visualization, where we use the real part for x and the imaginary part for y:

iterate[a_, b_, j_, n_] := Module[{start = b^j},
Table[a^i*start, {i, Range[-n, n]}]
]
Graphics[Circle[ReIm[#], r*Abs[#]] & /@ iterate[a, b, 0, 3]]


This shows the $0th$ spiral of packed circles, $3$ circles inward to our base circle and $3$ circles outward. To create the complete plane, we have to create $12$ columns, since q was $12$:

toCirle[z_, r_] := Disk[{Re[z], Im[z]}, Abs[z]*r];
pack = Table[iterate[a, b, j, 5], {j, 12}];
gr = Graphics[{EdgeForm[Black],
Map[{RandomColor[], toCirle[#, r] & /@ #} &, pack]},
PlotRange -> {{-10, 10}, {-10, 10}}]


Unfortunately, that is not enough, because the artist chose to use the logarithmic spirals through the circle centers as guide to divide each circle into different parts. In order to do this, we need to go further.

Let us make a cut here and divide the following into small sections where we look at the details. These details will be important for the overall approach

## Connection between a, b, r and the circles

As pointed out in the article, the complex numbers a and b are the generators of the circles. This means, all circle centers can be obtained by repeated multiplication. The base circle with the center {1,0} is given by $a^0\cdot b^0$ which is 1 (meaning re=1 and im=0).

Now each multiplication by a or by b gives the center of the circle that is next to the base circle. So 1*a*a=a^ gives the second circle in the direction of a. a^2*b shifts this last circle in the direct of b.

Note, that even a^(-3) is perfectly OK and gives the 3rd circle in the oposite direction. These are the small circles the fill the center. OK, one Manipulate says more than a thousand words. Let us create a dynamic table of all circles in a range for a^i*b^j. Note that, as pointed out in the article, the correct radius for each circle is Abs[a^i*b^j]*r where r is the radius we got from the solution of doyle.

Manipulate[
Graphics[{EdgeForm[Black], FaceForm[Gray], Table[
toCirle[a^i*b^j, r],
{i, i1, i2}, {j, j1, j2}]}, PlotRange -> {{-10, 10}, {-10, 10}}],
{{i1, 0}, -5, 3, 1},
{{i2, 0}, i1, 5, 1},
{{j1, 0}, -5, 3, 1},
{{j2, 0}, j1, 5, 1}
]


## Circles and spirals

We have seen that we can go from one circle to any neighbour by increasing (or decreasing) either i or j by 1. But what if we don't make jump from say a^3 to a^4 but a smooth transition? Well, the function for such a thing is easy because a^0=1 and a^1=a, so we can make a function a^3*a^t and let t run from 0 to 1.

Show[
Graphics[{
FaceForm[Gray],
toCirle[#, r] & /@ {a^3, a^4}}
],
ParametricPlot[ReIm[a^3*a^t], {t, 0, 1}, PlotStyle -> White]
]


This looks very much like the spirals that were used to divide the circles in the original art. So it seems if we pick out the center any circle next to our base circle, we can create spiral functions that go through the circles. Note that the approach of shifting a spiral to its neighbouring spiral is similar to shifting circles. Here is an example:

Show[
gr,
ParametricPlot[Table[ReIm[b^i a^t], {i, 12}], {t, -10, 10},
PlotRange -> {{-10, 10}, {-10, 10}}, PlotStyle -> White],
ParametricPlot[Table[ReIm[a^j b^t], {j, 5}], {t, -10, 10},
PlotRange -> {{-10, 10}, {-10, 10}}, PlotStyle -> White]
]


## Spiral functions inside circles

For our later approach, I want to be able to draw the spiral only inside a circle. As we have seen, going from t=0 to t=1 will connect the centers of the circles. This is not what we want. We want values for t that start and end with the circle. Let's make the plot we did earlier again, but use values for t between -1/3 and 1/3

OK, that looks promising. Remember, we know the center of this circle with a^3 and we know its radius with Abs[a^3]*r. What are the bounds where our spiral is exactly on the radius? Let us ask FindRoot:

tb = t /. FindRoot[Abs[1 - a^t] - r, {t, #}] & /@ {-1/3, 1/3}
(* {-0.565183, 0.433533} *)


But wait! I haven't used a^3 at all! Correct. The good thing is that the bounds for the circles apply to each circle of the same spiral. Therefore I'm using the next neighbour of the base circle which is a for FindRoot. Look here:

Show[
Graphics[{
FaceForm[Gray],
toCirle[#, r] & /@ {a^3, a^4}}
],
ParametricPlot[ReIm[a^3*a^t], {t, tb[[1]], tb[[2]]},
PlotStyle -> White]
]


## What spirals did the artist use?

As it turns out she used the spirals of the following direct neighboring circles of the base circle:

spoints = {a*b^-1, a, b}
(* {1.46301 - 0.54185 I, 1.67036 + 0.343254 I, 0.927594 + 0.578172 I} *)


Let's make a small function that calculates their bounds and returns them with a spiral function. The spiral function will directly incorporate the i and j so that we can easily draw it on every circle we like

spiral[pt_] := Module[{t1, t2},
{t1, t2} =
Block[{t},
t /. FindRoot[Abs[1 - pt^(t)] - r, {t, #}] & /@ {-1/3, 1/3}];
{t1, t2, Function[{i, j, t}, a^i*b^j*pt^t]}
]


Now let's plot these 3 spirals inside our base circle {1,0}

Show[{
Graphics[Circle[{1, 0}, r]],
ParametricPlot[ReIm@#3[0, 0, t], {t, #1, #2}] & @@@ spiral /@ spoints
}]


Now, we can calculate the points of the spirals inside each circle, we have the radius of each circle and through the spiral's start and end points, we have 6 points on each circle.

## Creating polygons points for the parts of a circle

For each cake-part of a circle, we can now create a polygon by

• starting in the center
• creating points along a spiral outwards to the circle boundary
• going counterclockwise along the circle to the endpoint of the next spiral
• create points along this next spiral from the outer point to the center

However, one tiny point is missing. How do we create points along the circle from when we go from one spiral point to the next. That is not as hard as it sounds.

Assume you have two (complex) points that lie on a circle around a center. Then you can subdivide them and create arbitrarily many points between them that all lie on the circle.

circle[z1_, z2_, cent_] :=
Module[{zz1 = z1 - cent, zz2 = z2 - cent, r},
r = Abs[Mean[{zz1, zz2}]];
# + cent & /@
Nest[Riffle[#,
Function[zz, With[{m = Mean[zz]}, m/Abs[m]*Abs[zz1]]] /@
Partition[#, 2, 1]] &, {zz1, zz2}, 5]
]


Having this, we can create the points for all cake-parts of circle i, j defined by the provided spirals that divide the circle:

createCircleParts[spirals_, i_, j_] :=
Module[{center, outward, inward},
outward = Table[#3[i, j, t], {t, 0, #2, #2/10.}] & @@@ spirals;
inward = Table[#3[i, j, t], {t, 0, #1, #1/10.}] & @@@ spirals;
center = a^i*b^j;
{i, j,
Join[#1[[;; -2]], circle[#1[[-1]], #2[[-1]], center],
Reverse[#2[[;; -2]]]] & @@@
Partition[Join[outward, inward, {First[outward]}], 2, 1]}
]


The function returns {i,j, {part1, part2, ...}} and we will use i and j later for the coloring as it gives us information about the position of the circle.

To test this function, let us see what happens with the circle i=1 and j=2:

Graphics[{RandomColor[], Polygon[ReIm[#]]} & /@
Last@createCircleParts[spiral /@ spoints, 1, 2]
]


## Coloring of circles

For one circle we have the information i, j which encodes the global position and of course we have n cake-parts. An easy way would be to provide a list of colors and select a color depending on the information we have.

I could not really find a pattern in the coloring of the artists image, so lets keep it simple but let us use equivalent colors:

cols = {Black, RGBColor[0.078, 0.71, 0.964], Orange, Red, Darker@Green, Purple};
colorCircleParts[{i_, j_, parts_}, col_List] :=
Table[{col[[Mod[i + j + n, Length[col]] + 1]],
Polygon[ReIm@parts[[n]]]}, {n, Length[parts]}]


## Putting everything together

The last thing we need to do is to create a table containing the circles and their parts for a range of i and j values. Then we color the circle parts and display them:

all = Table[colorCircleParts[createCircleParts[spiral /@ spoints, i, j],
cols], {i, -5, 6}, {j, 0, 12}];
Graphics[all, PlotRange -> {{-20, 20}, {-20, 20}}]


## Aftermath: Getting something close to the artist's work

The webpage of the artist suggests that

The central part of the picture shows the Doyle spiral circle packing P = 2 Q = 12.

That is not true. The values of P and Q define how many circles you need to close one loop. Additionally, the rotation of the circles in the artist's work is clockwise while in mathematics, we usually prefer to go counter-clockwise.

Lucky for us, this is no big deal because to go clockwise we just need to conjugate our complex values a and b. After printing the painting and counting the circles (and paying absolutely no attention to Wjx's comment who already found out that the values are off), I discovered that the painting uses P=3 and Q=8.

Let me show you what that means:

pqPlot[p_, q_] := Module[{a, b, r, c1, c2},
{a, b, r} = doyle[p, q];
{a, b} = Conjugate /@ {a, b};
c1 = toCirle[#, r] & /@ NestList[a*# &, 1, p - 1];
c2 = toCirle[#, r] & /@ NestList[b*# &, 1, q - 1];
Graphics[{EdgeForm[Black], FaceForm[LightYellow], c2,
FaceForm[LightBlue], c1, FaceForm[LightGreen], EdgeForm[Thick],
toCirle[1, r], toCirle[a^p, r]}]
]

pqPlot[3, 8]


If you include the inner base circle in your counting, you have 3 circles in the first and 8 circles in the other direction until you reach the outer end circle.

Taking this into account and including some of the colors in the original painting, we can come up with a very close optical copy of what the artist did. I played around with the plot-range to make it fit.

{a, b, r} = doyle[3, 8];
{a, b} = Conjugate /@ {a, b};
spoints = {a*b^-1, a, b};
cols = {GrayLevel[0.1], RGBColor[0.078, 0.71, 0.964],
RGBColor[0.95, 0.36, 0.09], RGBColor[0.77, 0.17, 0.12],
RGBColor[0.07, 0.6, 0.25], RGBColor[.32, .24, .55]};
range = 5.585;

Graphics[
Table[colorCircleParts[createCircleParts[spiral /@ spoints, i, j],
cols], {i, -5, 2}, {j, 0, 7}],
PlotRange -> {{-range, range}, {-range, range}}
]


• " I do realize that this answer is not as verbose...". I LOL'd. And +1 for a really neat answer! – ciao May 29 '17 at 8:28
• Great solution so far, looking forward to seeing even better version! Note that colored areas of the same color never touch on the original paining. (If I don't ask too much, in the final version maybe you can attach filtered pictures of areas of the same color only, lets say only orange areas from original; I believe there will be a neat visual effect in such picture) – VividD May 29 '17 at 8:32
• @VividD Look what we've got here. – halirutan May 29 '17 at 14:51
• @VividD Unfortunately, I cannot find any useful pattern in the drawing of the artist. I have the following information for coloring inside my plot: The circle part number, the number of the "polar column" if you go around the main spiral, the number of the "radial row" (because going around on the main spiral will overjump one set of spiral circles) – halirutan May 29 '17 at 14:55
• @ciao Now it's appropriately packed with information :) – halirutan May 29 '17 at 17:05

(Not an answer but an extended comment.)

## Simple mush-ups

In view of question's main part:

Furthermore, can the code be written to generate a parametrized version of the painting so that one could produce and enjoy an unlimited number of pieces of art?

using halirutan's great answer we can do some mush-ups.

Here is an example.

AbsoluteTiming[
doyleSpiralImages = Flatten@
Table[(lcols = RandomSample[cols];
all =
Table[colorCircleParts[
createCircleParts[spiral /@ spoints, i, j], lcols], {i, -5,
6}, {j, 0, 12}];
Image[Graphics[all,
PlotRange -> {{-20, 20}, {-20, 20}}]]), {20}];
]
doyleSpiralImagesBW =
ColorConvert[#, "Grayscale"] & /@ doyleSpiralImages;
doyleSpiralImagesBWBin =
Flatten@Table[
Binarize[#, b] & /@ doyleSpiralImagesBW, {b, {0.15, 0.45}}];

(* {31.5697, Null} *)

AbsoluteTiming[
nc = 16;
directBlendingImages =
Blend[Colorize[#,
ColorFunction ->
RandomChoice[{"BrightBands", "FruitPunchColors",
RandomChoice[doyleSpiralImagesBWBin, nc],
RandomReal[1, nc]]], {25}];
]

(* {42.7025, Null} *)

Multicolumn[ColorNegate /@ directBlendingImages, 6]


## Mush-ups with larger collections of images

Ideally, we can use procedures in the spirit of the one above over a large collection of Doyle spiral images. Such images can be made from scratch or obtained from the Web.

For example, if such an image collection is fed to the neural net functions in Mathematica / WL, in principle we will be able to obtain new spiral images by examining the layers.

• Yes, thanks, you address an important aspect of the question. – VividD Jun 2 '17 at 21:28
• Also, I would love to understand neural network application to this problem, but I will leave it for the future, when I learn more about such Mathematica capabilities. I liked the idea at first sight though! – VividD Jun 2 '17 at 21:34
• I am reading your blog! – VividD Jun 2 '17 at 21:43
• @VividD Thanks! I wanted to try the Doyle spiral images generation with machine learning algorithms, but currently it is too time consuming for me to do it... – Anton Antonov Jun 2 '17 at 22:09