Back a few years I challenged people to
give an actual algorithmic IFS for the Mandelbrot book cover: http://mathweb.scranton.edu/monks/pubs/pharaoh/pharaoh.html

Back in 2002 and 2004 I failed even with Dr. Ken Monks' circle drawing definition in hand.
The main reason was that you actually need the Klein group
to get a good algorithm.

My efforts have yielded some results. https://lh5.googleusercontent.com/-tjrQm_ICW-I/VgRjNkwIhdI/AAAAAAAAnso/NFisSPGgC9U/w610-h534-no/IFS_4part_Apollonian_strip_xequalssqrt5_inversion_3000000.png

Still not perfect,
but much closer:

The first important thing that this IFS version seems to show is that
the circles drawing approach is not giving the actually Klein group limit set ( the spiraling sections are left out).
The Ken Monks picture is appears to be a bad approximation of the actual Klein group?!
Now, why would that be important in a Mathematical sense?

(* mathematica*)
(* Klein group from Indra's Pearl's page 265*)
(* Gaussian sphere projections*)
Clear[f, dlst, pt, cr, ptlst, ptlst2, cr2]
dlst = Table[ Random[Integer, {1, 4}], {n, 3000000}];
allColors = ColorData["Legacy"][[3, 1]];
firstCols = {"Tomato", "DeepNaplesYellow", "DodgerBlue", "Cyan", "Red",
"Blue", "DeepNaplesYellow", "Green", "Magenta", "DarkOrchid",
"LightSalmon", "LightPink", "Sienna", "Green", "Mint", "DarkSlateGray",
"ManganeseBlue", "SlateGray", "DarkOrange", "MistyRose",
"DeepNaplesYellow", "GoldOchre", "SapGreen", "Yellow", "LimeGreen"};
cols = ColorData["Legacy", #] & /@
Join[firstCols, Complement[allColors, firstCols]];
s = Sqrt[3];
z = x + I*y;
(* Klein group {a,b} generators*)
f1a[z_] = ComplexExpand[ z/(-2*I*z + 1)];
f2a[z_] = ComplexExpand[ ((1 - I)*z + 1)/(z + (1 + I))];

(* projection {a,b,a^(-1),b^(-1)}*)

f[1, {x_, y_}] = N[{-Re[f2a[z]], -Im[f2a[z]]}];
f[2, {x_,
y_}] = (N[{Re[f1a[z]],
Im[f1a[z]]}/({Re[f1a[z]], Im[f1a[z]]}.{Re[f1a[z]], Im[f1a[z]]})]);
f[3, {x_, y_}] =
N[{Re[f2a[z]],
Im[f2a[z]]}/({Re[f2a[z]], Im[f2a[z]]}.{Re[f2a[z]], Im[f2a[z]]})];
f[4, {x_, y_}] = N[{-Re[f1a[z]], -Im[f1a[z]]}];


pt = {0.5, 0.5};
cr[n_] := cr[n] = cols[[12 + n]]
cr2[n_] := cr2[n] = cols[[4 + n]]
cr3[n_] := cr3[n] = cols[[8 + n]]
cr4[n_] := cr4[n] = cols[[n]]
aa = Table[pt = f[dlst[[j]], pt], {j, Length[dlst]}];
dd = Table[{aa[[i,
2]]/(1 + {aa[[i, 1]] + Sqrt[5], aa[[i, 2]]}.{aa[[i, 1]] + Sqrt[5],
aa[[i, 2]]}), (Sqrt[5] +
aa[[i, 1]])/(1 + {aa[[i, 1]] + Sqrt[5],
aa[[i, 2]]}.{aa[[i, 1]] + Sqrt[5], aa[[i, 2]]})}, {i, Length[aa]}]; dd1 = Table[{-aa[[i,
2]]/(1 + {aa[[i, 1]] + Sqrt[5], aa[[i, 2]]}.{aa[[i, 1]] + Sqrt[5],
aa[[i, 2]]}), (Sqrt[5] +
aa[[i, 1]])/(1 + {aa[[i, 1]] + Sqrt[5],
aa[[i, 2]]}.{aa[[i, 1]] + Sqrt[5], aa[[i, 2]]})}, {i, Length[aa]}]; ListPlot[{dd, dd1}, PlotStyle -> {PointSize -> 0.001, {Red, Blue}},
PlotRange -> All]
ptlst4 :=
Point[Developer`ToPackedArray[dd],
VertexColors -> Developer`ToPackedArray[cr4 /@ dlst]];
ptlst5 :=
Point[Developer`ToPackedArray[dd1],
VertexColors -> Developer`ToPackedArray[cr2 /@ dlst]];
g4 = Graphics[{PointSize[.001], ptlst4, ptlst5}, AspectRatio -> Automatic,
PlotRange -> {{-0.4, 0.4}, {-0.1, 0.6}}, ImageSize -> 1000]

(* end*)

--- SoupGate-Win32 v1.05
* Origin: fsxNet Usenet Gateway (21:1/5)