Welcome, guest | Sign In | My Account | Store | Cart

It draws Apollonian Gasket Fractal for any n using Descartes theorem.

This is not the standard way though. It simply randomly finds 3 tangent circles at each iteration and tries to add new circles. The good thing is it can start w/ any arbitrary configuration of main circles, unlike the standard way. The bad thing is some circles will stay missing because of random selections they never get a chance. You can increase the maxIt to get better result but it would slow it down a lot.

Python, 159 lines
 ``` 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159``` ```# Generalized Apollonian Gasket Fractal # http://en.wikipedia.org/wiki/Apollonian_Gasket # http://en.wikipedia.org/wiki/Descartes'_theorem # http://en.wikibooks.org/wiki/Fractals/Apollonian_fractals # FB - 20120131 import math import cmath import random from PIL import Image, ImageDraw imgx = 512 imgy = 512 image = Image.new("RGB", (imgx, imgy)) draw = ImageDraw.Draw(image) maxIt = 2000 # of iterations n = random.randint(2, 7) # of main circles > 1 print 'Number of main circles: ' + str(n) eps = 0.0001 def appEq(u, v): # approximately equal excluding signs return abs(abs(u) - abs(v)) <= eps # test if 2 circles are approximately equal def appEqCir(a, b): return (abs(a - b) <= eps and appEq(a, b)) # test if 2 circles are tangent def isTanCir(a, b): return (appEq(a - b, abs(a) + abs(b)) or appEq(a - b, abs(a) - abs(b))) # test if 2 circles are intersecting def isIntCir(a, b): dist = abs(a - b) # distance between 2 centers rmin = min(abs(a), abs(b)) # min radius rmax = max(abs(a), abs(b)) # max radius return (dist + eps < rmin + rmax and dist + rmin > rmax + eps) # test if 3 circles are tangent to each other def isTangent(a, b, c): flag = True if appEqCir(a, b) or appEqCir(a, c) or appEqCir(b, c): # all circles must be different flag = False if not isTanCir(a, b): flag = False if not isTanCir(a, c): flag = False if not isTanCir(b, c): flag = False return flag # test if the 4th circle is tangent to previous 3 def isTan4(a, b, c, d): return isTangent(a, b, d) and isTangent(a, c, d) and isTangent(b, c, d) # input: 3 circles w/ each as a tuple: (complex(x,y),r) # output: 4 tangent circles as a tuple def cif(a, b, c): # curvatures of the circles k1 = 1.0 / a k2 = 1.0 / b k3 = 1.0 / c # curvatures of tangent circles temp0 = 2.0 * math.sqrt(k1 * k2 + k1 * k3 + k2 * k3) temp1 = k1 + k2 + k3 k40 = temp1 + temp0 k41 = temp1 - temp0 # centers of tangent circles temp0 = 2.0 * cmath.sqrt(k1*k2*a*b+k1*k3*a*c+k2*k3*b*c) temp1 = k1 * a + k2 * b + k3 * c c40 = (temp1 + temp0) / k40 c41 = (temp1 - temp0) / k41 c42 = (temp1 - temp0) / k40 c43 = (temp1 + temp0) / k41 # radius of tangent circles r0 = 1.0 / k40 r1 = 1.0 / k41 # there are 4 solutions and only 2 are correct (tangent) circles # those are c0 and c1 or c2 and c3 c0 = (c40, r0) c1 = (c41, r1) c2 = (c42, r0) c3 = (c43, r1) if isTan4(a, b, c, c0) and isTan4(a, b, c, c1): return (c0, c1) else: return (c2, c3) # pixel coordinate calculation def pxy(x, y): kx = int((x - xa) / (xb - xa) * (imgx - 1)) ky = int((y - ya) / (yb - ya) * (imgy - 1)) return (kx, ky) # generate the main circles: n circles in a circle circles = [] a = 2.0 * math.pi / n # radius of the inner circles r = math.hypot(math.cos(a) - math.cos(0), math.sin(a) - math.sin(0)) / 2.0 r0 = 1.0 + r # radius of the outer circle r1 = r0 - 2.0 * r # radius of the center circle # origin ox = r0 oy = r0 rmin = 2.0 * r0 / math.hypot(imgx, imgy) # min radius that can be drawn # drawing area xa = ox - r0 xb = ox + r0 ya = oy - r0 yb = oy + r0 circles.append((complex(ox, oy), -r0)) # add the outer circle w/ - curvature for j in range(n): # add the main circles cx = math.cos(a * j) + ox cy = math.sin(a * j) + oy circles.append((complex(cx, cy), r)) if r1 > rmin: # if possible add a center circle circles.append((complex(ox, oy), r1)) p0 = 0 for i in range(maxIt): p = int(100 * i / (maxIt - 1)) if p > p0: p0 = p print "%" + str(p) # randomly select 3 circles while True: a = random.choice(circles) b = random.choice(circles) c = random.choice(circles) # test if all 3 circles are tangent to each other if isTangent(a, b, c): break # find inversion circles d = cif(a, b, c) for e in d: dx = e.real dy = e.imag dr = e # ignore the outer circle (- curvature) solutions and too small circles if dx > 0.0 and dy > 0.0 and dr > rmin: # add the new circle if it was not added before and not intersecting flag = True for cir in circles: # new circle must not be equal/intersecting to any existing circle if appEqCir(cir, e) or isIntCir(cir, e): flag = False break if flag: circles.append(e) # draw all circles for cir in circles: dx = cir.real dy = cir.imag dr = abs(cir) try: draw.ellipse((pxy(dx - dr, dy - dr), pxy(dx + dr, dy + dr))) except: pass image.save("ApollonianGasket_" + str(n) + ".png", "PNG") ``` yishu 5 years, 8 months ago

thank you so much for your code! However, there is one equation I have not understand. Could you explain the equation? How to understand this equation? r = math.hypot(math.cos(a) - math.cos(0), math.sin(a) - math.sin(0)) / 2.0

Thank you so much!! FB36 (author) 5 years, 8 months ago

If you put N circles (as big as possible) of same size inside a bigger circle, that equation gives you the radius of each smaller circles. It is the distance between 2 centers of circles inside (next to each other) divided by 2. yishu 5 years, 7 months ago

Thank you for your answer. But I still cannot figure out the geometry realtion here. Do you mind to explain it more clearly?

Besides, in your code, the generation of each circle is random. Do you have any idea on generating circles by generations, which means generating circles from big to small?

Thanks again! FB36 (author) 5 years, 7 months ago

I don't think I can explain better without drawings etc also it was years ago for me. But I think my idea was instead of starting with a unit circle and trying to figure out how to pack N circles inside of it, I realized the problem would be simpler if the N circles all have centers on the circumference of unit circle (not inside of it). That equation comes from that. You need to a drawing to see it. And once you have N circles, it it not hard to calculate the radius of a circle that surrounds those N circles and the radius of the small circle in the center of N circles.

Generating the circles in generations/levels is something I couldn't figure out myself. That is why I made it random. But generating them randomly of course is not a good efficient way. That is why I gave up on this method and switched to a much better way later. I don't know if you tried my other programs in this website titled "Circle Inversion Fractals". They also draw Apollonian Gasket fractals but in a much better way. yishu 5 years, 7 months ago

I found a paper that introduced the packing method generation by generation. http://www.math.columbia.edu/~staff/RTGPapers/ApollonianSums.pdf

It is not easy for me to write my own packing code because I am a beginner just startin learning the packing method. I noticed you assign colors (or images) to different circles, is there any special tips here? FB36 (author) 5 years, 7 months ago

I think that paper you found has the solutions for both creating generations and coloring. (But it requires spending sometime to understand :-)

I usually avoid using any advanced coloring methods for any fractals to keep the code short and easy to understand. So I use either randomly assigned colors or just create colors using mod operation:

if k is fractal iteration number/level etc then:

red=(k%a) * (256/a)

green=(k%b) * (256/b)

blue=(k%c) * (256/c)

a,b,c values selected from 4,8,16,32,64,128 and a<>b,a<>c,b<>c yishu 5 years, 7 months ago

Thank you so much for your patience! I will take time to read the papers and try to understrand the coloring methods.

In 2D packing ,the ceter is straightforward and can be thought of as a point in the 2D-plane and be described by a complex number. Do you have any idea on how to represent the center of sphere in 3D packing? I am trying to find a way, but cannot figure it out now. FB36 (author) 5 years, 7 months ago

Since there are no 3d extension of complex numbers, I am guessing 3d solution must use 4d numbers: https://en.wikipedia.org/wiki/Quaternion John Park 5 years, 6 months ago

There are some missing circles after generation.. Created by FB36 on Wed, 1 Feb 2012 (MIT)