Recursively create apollonian gaskets [With solution]

Question

Apollonian gaskets = They are planar fractals generated from triples of circles, where each circle is tangent to the other two. In his drawing of the gasket, we start with two externally tangent circles which diameter is D1 and D2. Then we add a third circle which diameter is D1+D2 and to which the two original circles are internally tangent. This is the first generation of circles. Each subsequent generation of circles is constructed by applying the following scheme: For any three circles A, B C of any previous generations which are tangent to each other a new circle is constructed which is tangent to A,B,C. The new circle must differ from all circles constructed so far. When a generation is complete, i.e no other circle can be added, then the next generation of circles can start being constructed.

There is an additional stopping rule which prevents from generating infinitesimally small circles. A circle can be added to the gasket if and only if the lenght of its diameter is least minD which is a fixed positive value.

Input consists of one line with three decimal numbers D1, D2 and minD. The number are separated by spaces. The format is usual decimal format (see also the examples bellow) with no exponent part. It holds that 1.0 ≤ D1, D2 ≤ 1000.0, 0.001 ≤ minD ≤ D1+D2.

Ouput consists of one text line containing two decimal numbers L1 and L2. L1 represents the sum of areas of all circles in the gasket except for the bigggest circle. L2 represents the sum of perimeters of all circles in tin the gasket except for the bigggest circle. Both output values are rounded to 6 decimal digits. Decimal digits must be always present in the output even if some of them are zeros. Maximim output value is less than 107.

Input

17.000000 40.000000 1.000000

Output

2439.258588 835.263228

2

For given D1 and D2, I create this two circles like this (first iteration):

   double D1 = 17.00;
   double D2 = 40.00;
   double minD = 1.00;
   int i = 250, j = 350;
   comp.addCircle(i, j, (int) D2, randomColor);
   comp.addCircle(i + (int) D2 / 2 + (int) D1 / 2, j, (int) D1, randomColor);
   comp.addCircle(i + (int) D1 / 2, j, (int) (D1 + D2), randomColor);

UPDATE:

So, solution is based on Descartes' theorem. We well work with radius, not diameter, and Curvature, with is 1/r. We will use double for all calculation, but if you work with significantly small numbers, I would prefer BigDecimal. It will slow algorithm, and you should use external method for finding square root, because BigDecimal doesn't have any.

For given D1, D2, minD we modify code above for efficiency:

Some preparation:

 double D1 = sc.nextDouble() / 2;
 double D2 = sc.nextDouble() / 2;
 minD = sc.nextDouble() / 2;
 double D3 = D1 + D2;

So, first step looks like this:

Next step looks a little bit more complicated.

Assume we want to write a recursion to solve this problem, and according to Descartes' theorem, for given curvatures of three circles, tangent to each other, (pic. below) , we could find curvatures of two circles, but for our purposes, we need only small one, so, we can simplify formula to

this.curve = a.curve + b.curve + c.curve + 2 * Math.sqrt(Math.abs(a.curve * b.curve + a.curve * c.curve + b.curve * c.curve));

Lets take a look at Apollonian gaskets again: try to play with it. See? It is same gaskets, but with different start condition. And whats more important for us, is that it is symmetrical! So, we will calculate just a half, and then multiply result by two! Lets write a recursion! Inputs will be curvatures of three circles. No output, we will use change our global variables.

double radius_sum = 0.0;
double square_radius_sum = 0.0;

void createAG(double a, double b, double c){
 double n = a + b + c + Math.sqrt(a*b + a*c + b*c + 4.0); 
 if ((minD * n) < 1){
     radius_sum += 2. / n; //Remember about symmetry? 
     square_radius_sum += 2. * (1. / n) * (1. / n); //Remember about symmetry? 
     createAG(a, b, n);
     createAG(a, c, n);
     createAG(b, c, n);
   }
}

To find the result, we will use formulas to calculate area and perimeter of circle. Perimeter is length of circumference and equal to . Area is equal to , as you already know, because we already calculated it in previous step, otherwise we had to store every radius and do more calculations.

radius_sum  = 2 * Math.Pi * radius_sum;
square_radius_sum = Math.Pi * square_radius_sum;

But we forget about our first two circles! Let's fix it!

   radius_sum  += D1*2 + D2*2;
   square_radius_sum += D1*D1 + D2*D2;
   radius_sum  = 2 * Math.Pi * radius_sum;
   square_radius_sum = Math.Pi * square_radius_sum;

And there is always a room for improvement. For example, to use IEEE 754 in better way, I assume you will use 1. / x instead of 1 / x.

Thank you!

P.S. Copyright! This task (text and first picture of Apollonian gasket) is created by teachers at CTU, for course ALG. Picture of formulas is from Wikipedia. Everything else is public domain, if not patented, registered e.t.c.

Answer 1

So, solution is based on Descartes' theorem. We well work with radius, not diameter, and Curvature, with is 1/r. We will use double for all calculation, but if you work with significantly small numbers, I would prefer BigDecimal. It will slow algorithm, and you should use external method for finding square root, because BigDecimal doesn't have any.

For given D1, D2, minD we modify code above for efficiency:

Some preparation:

 double D1 = sc.nextDouble() / 2;
 double D2 = sc.nextDouble() / 2;
 minD = sc.nextDouble() / 2;
 double D3 = D1 + D2;

So, first step looks like this:

Next step looks a little bit more complicated.

Assume we want to write a recursion to solve this problem, and according to Descartes' theorem, for given curvatures of three circles, tangent to each other, (pic. below) , we could find curvatures of two circles, but for our purposes, we need only small one, so, we can simplify formula to

this.curve = a.curve + b.curve + c.curve + 2 * Math.sqrt(Math.abs(a.curve * b.curve + a.curve * c.curve + b.curve * c.curve));

Lets take a look at Apollonian gaskets again: try to play with it. See? It is same gaskets, but with different start condition. And whats more important for us, is that it is symmetrical! So, we will calculate just a half, and then multiply result by two! Lets write a recursion! Inputs will be curvatures of three circles. No output, we will use change our global variables.

double radius_sum = 0.0;
double square_radius_sum = 0.0;

void createAG(double a, double b, double c){
 double n = a + b + c + Math.sqrt(a*b + a*c + b*c + 4.0); 
 if ((minD * n) < 1){
     radius_sum += 2. / n; //Remember about symmetry? 
     square_radius_sum += 2. * (1. / n) * (1. / n); //Remember about symmetry? 
     createAG(a, b, n);
     createAG(a, c, n);
     createAG(b, c, n);
   }
}

To find the result, we will use formulas to calculate area and perimeter of circle. Perimeter is length of circumference and equal to . Area is equal to , as you already know, because we already calculated it in previous step, otherwise we had to store every radius and do more calculations.

radius_sum  = 2 * Math.Pi * radius_sum;
square_radius_sum = Math.Pi * square_radius_sum;

But we forget about our first two circles! Let's fix it!

   radius_sum  += D1*2 + D2*2;
   square_radius_sum += D1*D1 + D2*D2;
   radius_sum  = 2 * Math.Pi * radius_sum;
   square_radius_sum = Math.Pi * square_radius_sum;

And there is always a room for improvement. For example, to use IEEE 754 in better way, I assume you will use 1. / x instead of 1 / x.

Thank you!

P.S. Copyright! This task (text and first picture of Apollonian gasket) is created by teachers at CTU, for course ALG. Picture of formulas is from Wikipedia. Everything else is public domain, if not patented, registered e.t.c.