1 | /* -*- Mode: C; tab-width: 4 -*- */ |
---|
2 | /* apollonian --- Apollonian Circles */ |
---|
3 | |
---|
4 | #if 0 |
---|
5 | static const char sccsid[] = "@(#)apollonian.c 5.02 2001/07/01 xlockmore"; |
---|
6 | #endif |
---|
7 | |
---|
8 | /*- |
---|
9 | * Copyright (c) 2000, 2001 by Allan R. Wilks <allan@research.att.com>. |
---|
10 | * |
---|
11 | * Permission to use, copy, modify, and distribute this software and its |
---|
12 | * documentation for any purpose and without fee is hereby granted, |
---|
13 | * provided that the above copyright notice appear in all copies and that |
---|
14 | * both that copyright notice and this permission notice appear in |
---|
15 | * supporting documentation. |
---|
16 | * |
---|
17 | * This file is provided AS IS with no warranties of any kind. The author |
---|
18 | * shall have no liability with respect to the infringement of copyrights, |
---|
19 | * trade secrets or any patents by this file or any part thereof. In no |
---|
20 | * event will the author be liable for any lost revenue or profits or |
---|
21 | * other special, indirect and consequential damages. |
---|
22 | * |
---|
23 | * radius r = 1 / c (curvature) |
---|
24 | * |
---|
25 | * Descartes Circle Theorem: (a, b, c, d are curvatures of tangential circles) |
---|
26 | * Let a, b, c, d be the curvatures of for mutually (externally) tangent |
---|
27 | * circles in the plane. Then |
---|
28 | * a^2 + b^2 + c^2 + d^2 = (a + b + c + d)^2 / 2 |
---|
29 | * |
---|
30 | * Complex Descartes Theorem: If the oriented curvatues and (complex) centers |
---|
31 | * of an oriented Descrates configuration in the plane are a, b, c, d and |
---|
32 | * w, x, y, z respectively, then |
---|
33 | * a^2*w^2 + b^2*x^2 + c^2*y^2 + d^2*z^2 = (aw + bx + cy + dz)^2 / 2 |
---|
34 | * In addition these quantities satisfy |
---|
35 | * a^2*w + b^2*x + c^2*y + d^2*z = (aw + bx + cy + dz)(a + b + c + d) / 2 |
---|
36 | * |
---|
37 | * Enumerate root integer Descartes quadruples (a,b,c,d) satisfying the |
---|
38 | * Descartes condition: |
---|
39 | * 2(a^2+b^2+c^2+d^2) = (a+b+c+d)^2 |
---|
40 | * i.e., quadruples for which no application of the "pollinate" operator |
---|
41 | * z <- 2(a+b+c+d) - 3*z, |
---|
42 | * where z is in {a,b,c,d}, gives a quad of strictly smaller sum. This |
---|
43 | * is equivalent to the condition: |
---|
44 | * sum(a,b,c,d) >= 2*max(a,b,c,d) |
---|
45 | * which, because of the Descartes condition, is equivalent to |
---|
46 | * sum(a^2,b^2,c^2,d^2) >= 2*max(a,b,c,d)^2 |
---|
47 | * |
---|
48 | * |
---|
49 | * Todo: |
---|
50 | * Add a small font |
---|
51 | * |
---|
52 | * Revision History: |
---|
53 | * 25-Jun-2001: Converted from C and Postscript code by David Bagley |
---|
54 | * Original code by Allan R. Wilks <allan@research.att.com>. |
---|
55 | * |
---|
56 | * From Circle Math Science News April 21, 2001 VOL. 254-255 |
---|
57 | * http://www.sciencenews.org/20010421/toc.asp |
---|
58 | * Apollonian Circle Packings Assorted papers from Ronald L Graham, |
---|
59 | * Jeffrey Lagarias, Colin Mallows, Allan Wilks, Catherine Yan |
---|
60 | * http://front.math.ucdavis.edu/math.NT/0009113 |
---|
61 | * http://front.math.ucdavis.edu/math.MG/0101066 |
---|
62 | * http://front.math.ucdavis.edu/math.MG/0010298 |
---|
63 | * http://front.math.ucdavis.edu/math.MG/0010302 |
---|
64 | * http://front.math.ucdavis.edu/math.MG/0010324 |
---|
65 | */ |
---|
66 | |
---|
67 | #ifdef STANDALONE |
---|
68 | #define MODE_apollonian |
---|
69 | #define PROGCLASS "Apollonian" |
---|
70 | #define HACK_INIT init_apollonian |
---|
71 | #define HACK_DRAW draw_apollonian |
---|
72 | #define apollonian_opts xlockmore_opts |
---|
73 | #define DEFAULTS "*delay: 1000000 \n" \ |
---|
74 | "*count: 64 \n" \ |
---|
75 | "*cycles: 20 \n" \ |
---|
76 | "*ncolors: 64 \n" |
---|
77 | #include "xlockmore.h" /* in xscreensaver distribution */ |
---|
78 | #include "erase.h" |
---|
79 | #else /* STANDALONE */ |
---|
80 | #include "xlock.h" /* in xlockmore distribution */ |
---|
81 | #endif /* STANDALONE */ |
---|
82 | |
---|
83 | #ifdef MODE_apollonian |
---|
84 | |
---|
85 | #define DEF_ALTGEOM "True" |
---|
86 | #define DEF_LABEL "True" |
---|
87 | |
---|
88 | static Bool altgeom; |
---|
89 | static Bool label; |
---|
90 | |
---|
91 | static XrmOptionDescRec opts[] = |
---|
92 | { |
---|
93 | {(char *) "-altgeom", (char *) ".apollonian.altgeom", XrmoptionNoArg, (caddr_t) "on"}, |
---|
94 | {(char *) "+altgeom", (char *) ".apollonian.altgeom", XrmoptionNoArg, (caddr_t) "off"}, |
---|
95 | {(char *) "-label", (char *) ".apollonian.label", XrmoptionNoArg, (caddr_t) "on"}, |
---|
96 | {(char *) "+label", (char *) ".apollonian.label", XrmoptionNoArg, (caddr_t) "off"}, |
---|
97 | }; |
---|
98 | static argtype vars[] = |
---|
99 | { |
---|
100 | {(caddr_t *) & altgeom, (char *) "altgeom", (char *) "AltGeom", (char *) DEF_ALTGEOM, t_Bool}, |
---|
101 | {(caddr_t *) & label, (char *) "label", (char *) "Label", (char *) DEF_LABEL, t_Bool}, |
---|
102 | }; |
---|
103 | static OptionStruct desc[] = |
---|
104 | { |
---|
105 | {(char *) "-/+altgeom", (char *) "turn on/off alternate geometries (off euclidean space, on includes spherical and hyperbolic)"}, |
---|
106 | {(char *) "-/+label", (char *) "turn on/off alternate space and number labeling"}, |
---|
107 | }; |
---|
108 | |
---|
109 | ModeSpecOpt apollonian_opts = |
---|
110 | {sizeof opts / sizeof opts[0], opts, sizeof vars / sizeof vars[0], vars, desc}; |
---|
111 | |
---|
112 | #ifdef DOFONT |
---|
113 | extern XFontStruct *getFont(Display * display); |
---|
114 | #endif |
---|
115 | |
---|
116 | #ifdef USE_MODULES |
---|
117 | ModStruct apollonian_description = |
---|
118 | {"apollonian", "init_apollonian", "draw_apollonian", "release_apollonian", |
---|
119 | "init_apollonian", "init_apollonian", (char *) NULL, &apollonian_opts, |
---|
120 | 1000000, 64, 20, 1, 64, 1.0, "", |
---|
121 | "Shows Apollonian Circles", 0, NULL}; |
---|
122 | |
---|
123 | #endif |
---|
124 | |
---|
125 | typedef struct { |
---|
126 | int a, b, c, d; |
---|
127 | } apollonian_quadruple; |
---|
128 | |
---|
129 | typedef struct { |
---|
130 | double e; /* euclidean bend */ |
---|
131 | double s; /* spherical bend */ |
---|
132 | double h; /* hyperbolic bend */ |
---|
133 | double x, y; /* euclidean bend times euclidean position */ |
---|
134 | } circle; |
---|
135 | enum space { |
---|
136 | euclidean = 0, spherical, hyperbolic |
---|
137 | }; |
---|
138 | |
---|
139 | const char * space_string[] = { |
---|
140 | "euclidean", |
---|
141 | "spherical", |
---|
142 | "hyperbolic" |
---|
143 | }; |
---|
144 | |
---|
145 | /* |
---|
146 | Generate Apollonian packing starting with a quadruple of circles. |
---|
147 | The four input lines each contain the 5-tuple (e,s,h,x,y) representing |
---|
148 | the circle with radius 1/e and center (x/e,y/e). The s and h is propagated |
---|
149 | like e, x and y, but can differ from e so as to represent different |
---|
150 | geometries, spherical and hyperbolic, respectively. The "standard" picture, |
---|
151 | for example (-1, 2, 2, 3), can be labeled for the three geometries. |
---|
152 | Origins of circles z1, z2, z3, z4 |
---|
153 | a * z1 = 0 |
---|
154 | b * z2 = (a+b)/a |
---|
155 | c * z3 = (q123 + a * i)^2/(a*(a+b)) where q123 = sqrt(a*b+a*c+b*c) |
---|
156 | d * z4 = (q124 + a * i)^2/(a*(a+b)) where q124 = q123 - a - b |
---|
157 | If (e,x,y) represents the Euclidean circle (1/e,x/e,y/e) (so that e is |
---|
158 | the label in the standard picture) then the "spherical label" is |
---|
159 | (e^2+x^2+y^2-1)/(2*e) (an integer!) and the "hyperbolic label", is |
---|
160 | calulated by h + s = e. |
---|
161 | */ |
---|
162 | circle examples[][4] = { |
---|
163 | { /* double semi-bounded */ |
---|
164 | { 0, 0, 0, 0, 1}, |
---|
165 | { 0, 0, 0, 0, -1}, |
---|
166 | { 1, 1, 1, -1, 0}, |
---|
167 | { 1, 1, 1, 1, 0} |
---|
168 | }, |
---|
169 | #if 0 |
---|
170 | { /* standard */ |
---|
171 | {-1, 0, -1, 0, 0}, |
---|
172 | { 2, 1, 1, 1, 0}, |
---|
173 | { 2, 1, 1, -1, 0}, |
---|
174 | { 3, 2, 1, 0, 2} |
---|
175 | }, |
---|
176 | { /* next simplest */ |
---|
177 | {-2, -1, -1, 0.0, 0}, |
---|
178 | { 3, 2, 1, 0.5, 0}, |
---|
179 | { 6, 3, 3, -2.0, 0}, |
---|
180 | { 7, 4, 3, -1.5, 2} |
---|
181 | }, |
---|
182 | { /* */ |
---|
183 | {-3, -2, -1, 0.0, 0}, |
---|
184 | { 4, 3, 1, 1.0 / 3.0, 0}, |
---|
185 | {12, 7, 5, -3.0, 0}, |
---|
186 | {13, 8, 5, -8.0 / 3.0, 2} |
---|
187 | }, |
---|
188 | { /* Mickey */ |
---|
189 | {-3, -2, -1, 0.0, 0}, |
---|
190 | { 5, 4, 1, 2.0 / 3.0, 0}, |
---|
191 | { 8, 5, 3, -4.0 / 3.0, -1}, |
---|
192 | { 8, 5, 3, -4.0 / 3.0, 1} |
---|
193 | }, |
---|
194 | { /* */ |
---|
195 | {-4, -3, -1, 0.00, 0}, |
---|
196 | { 5, 4, 1, 0.25, 0}, |
---|
197 | {20, 13, 7, -4.00, 0}, |
---|
198 | {21, 14, 7, -3.75, 2} |
---|
199 | }, |
---|
200 | { /* Mickey2 */ |
---|
201 | {-4, -2, -2, 0.0, 0}, |
---|
202 | { 8, 4, 4, 1.0, 0}, |
---|
203 | { 9, 5, 4, -0.75, -1}, |
---|
204 | { 9, 5, 4, -0.75, 1} |
---|
205 | }, |
---|
206 | { /* Mickey3 */ |
---|
207 | {-5, -4, -1, 0.0, 0}, |
---|
208 | { 7, 6, 1, 0.4, 0}, |
---|
209 | {18, 13, 5, -2.4, -1}, |
---|
210 | {18, 13, 5, -2.4, 1} |
---|
211 | }, |
---|
212 | { /* */ |
---|
213 | {-6, -5, -1, 0.0, 0}, |
---|
214 | { 7, 6, 1, 1.0 / 6.0, 0}, |
---|
215 | {42, 31, 11, -6.0, 0}, |
---|
216 | {43, 32, 11, -35.0 / 6.0, 2} |
---|
217 | }, |
---|
218 | { /* */ |
---|
219 | {-6, -3, -3, 0.0, 0}, |
---|
220 | {10, 5, 5, 2.0 / 3.0, 0}, |
---|
221 | {15, 8, 7, -1.5, 0}, |
---|
222 | {19, 10, 9, -5.0 / 6.0, 2} |
---|
223 | }, |
---|
224 | { /* asymmetric */ |
---|
225 | {-6, -5, -1, 0.0, 0.0}, |
---|
226 | {11, 10, 1, 5.0 / 6.0, 0.0}, |
---|
227 | {14, 11, 3, -16.0 / 15.0, -0.8}, |
---|
228 | {15, 12, 3, -0.9, 1.2} |
---|
229 | }, |
---|
230 | #endif |
---|
231 | /* Non integer stuff */ |
---|
232 | #define DELTA 2.154700538 /* ((3+2*sqrt(3))/3) */ |
---|
233 | { /* 3 fold symmetric bounded (x, y calculated later) */ |
---|
234 | { -1, -1, -1, 0.0, 0.0}, |
---|
235 | {DELTA, DELTA, DELTA, 1.0, 0.0}, |
---|
236 | {DELTA, DELTA, DELTA, 1.0, -1.0}, |
---|
237 | {DELTA, DELTA, DELTA, -1.0, 1.0} |
---|
238 | }, |
---|
239 | { /* semi-bounded (x, y calculated later) */ |
---|
240 | #define ALPHA 2.618033989 /* ((3+sqrt(5))/2) */ |
---|
241 | { 1.0, 1.0, 1.0, 0, 0}, |
---|
242 | { 0.0, 0.0, 0.0, 0, -1}, |
---|
243 | {1.0/(ALPHA*ALPHA), 1.0/(ALPHA*ALPHA), 1.0/(ALPHA*ALPHA), -1, 0}, |
---|
244 | { 1.0/ALPHA, 1.0/ALPHA, 1.0/ALPHA, -1, 0} |
---|
245 | }, |
---|
246 | { /* unbounded (x, y calculated later) */ |
---|
247 | /* #define PHI 1.618033989 *//* ((1+sqrt(5))/2) */ |
---|
248 | #define BETA 2.890053638 /* (PHI+sqrt(PHI)) */ |
---|
249 | { 1.0, 1.0, 1.0, 0, 0}, |
---|
250 | {1.0/(BETA*BETA*BETA), 1.0/(BETA*BETA*BETA), 1.0/(BETA*BETA*BETA), 1, 0}, |
---|
251 | { 1.0/(BETA*BETA), 1.0/(BETA*BETA), 1.0/(BETA*BETA), 1, 0}, |
---|
252 | { 1.0/BETA, 1.0/BETA, 1.0/BETA, 1, 0} |
---|
253 | } |
---|
254 | }; |
---|
255 | |
---|
256 | #define PREDEF_CIRCLE_GAMES (sizeof (examples) / (4 * sizeof (circle))) |
---|
257 | |
---|
258 | #if 0 |
---|
259 | Euclidean |
---|
260 | 0, 0, 1, 1 |
---|
261 | -1, 2, 2, 3 |
---|
262 | -2, 3, 6, 7 |
---|
263 | -3, 5, 8, 8 |
---|
264 | -4, 8, 9, 9 |
---|
265 | -3, 4, 12, 13 |
---|
266 | -6, 11, 14, 15 |
---|
267 | -5, 7, 18, 18 |
---|
268 | -6, 10, 15, 19 |
---|
269 | -7, 12, 17, 20 |
---|
270 | -4, 5, 20, 21 |
---|
271 | -9, 18, 19, 22 |
---|
272 | -8, 13, 21, 24 |
---|
273 | Spherical |
---|
274 | 0, 1, 1, 2 |
---|
275 | -1, 2, 3, 4 |
---|
276 | -2, 4, 5, 5 |
---|
277 | -2, 3, 7, 8 |
---|
278 | Hyperbolic |
---|
279 | -1, 1, 1, 1 |
---|
280 | 0, 0, 1, 3 |
---|
281 | -2, 3, 5, 6 |
---|
282 | -3, 6, 6, 7 |
---|
283 | #endif |
---|
284 | |
---|
285 | typedef struct { |
---|
286 | int size; |
---|
287 | XPoint offset; |
---|
288 | int geometry; |
---|
289 | circle c1, c2, c3, c4; |
---|
290 | int color_offset; |
---|
291 | int count; |
---|
292 | Bool label, altgeom; |
---|
293 | apollonian_quadruple *quad; |
---|
294 | #ifdef DOFONT |
---|
295 | XFontStruct *font; |
---|
296 | #endif |
---|
297 | int time; |
---|
298 | int game; |
---|
299 | } apollonianstruct; |
---|
300 | |
---|
301 | static apollonianstruct *apollonians = (apollonianstruct *) NULL; |
---|
302 | |
---|
303 | #define FONT_HEIGHT 19 |
---|
304 | #define FONT_WIDTH 15 |
---|
305 | #define FONT_LENGTH 20 |
---|
306 | #define MAX_CHAR 10 |
---|
307 | #define K 2.15470053837925152902 /* 1+2/sqrt(3) */ |
---|
308 | #define MAXBEND 100 /* Do not want configurable by user since it will take too |
---|
309 | much time if increased. */ |
---|
310 | |
---|
311 | static int |
---|
312 | gcd(int a, int b) |
---|
313 | { |
---|
314 | int r; |
---|
315 | |
---|
316 | while (b) { |
---|
317 | r = a % b; |
---|
318 | a = b; |
---|
319 | b = r; |
---|
320 | } |
---|
321 | return a; |
---|
322 | } |
---|
323 | |
---|
324 | static int |
---|
325 | isqrt(int n) |
---|
326 | { |
---|
327 | int y; |
---|
328 | |
---|
329 | if (n < 0) |
---|
330 | return -1; |
---|
331 | y = (int) (sqrt((double) n) + 0.5); |
---|
332 | return ((n == y*y) ? y : -1); |
---|
333 | } |
---|
334 | |
---|
335 | static void |
---|
336 | dquad(int n, apollonian_quadruple *quad) |
---|
337 | { |
---|
338 | int a, b, c, d; |
---|
339 | int counter = 0, B, C; |
---|
340 | |
---|
341 | for (a = 0; a < MAXBEND; a++) { |
---|
342 | B = (int) (K * a); |
---|
343 | for (b = a + 1; b <= B; b++) { |
---|
344 | C = (int) (((a + b) * (a + b)) / (4.0 * (b - a))); |
---|
345 | for (c = b; c <= C; c++) { |
---|
346 | d = isqrt(b*c-a*(b+c)); |
---|
347 | if (d >= 0 && (gcd(a,gcd(b,c)) <= 1)) { |
---|
348 | quad[counter].a = -a; |
---|
349 | quad[counter].b = b; |
---|
350 | quad[counter].c = c; |
---|
351 | quad[counter].d = -a+b+c-2*d; |
---|
352 | if (++counter >= n) { |
---|
353 | return; |
---|
354 | } |
---|
355 | } |
---|
356 | } |
---|
357 | } |
---|
358 | } |
---|
359 | (void) printf("found only %d below maximum bend of %d\n", |
---|
360 | counter, MAXBEND); |
---|
361 | for (; counter < n; counter++) { |
---|
362 | quad[counter].a = -1; |
---|
363 | quad[counter].b = 2; |
---|
364 | quad[counter].c = 2; |
---|
365 | quad[counter].d = 3; |
---|
366 | } |
---|
367 | return; |
---|
368 | } |
---|
369 | |
---|
370 | /* |
---|
371 | * Given a Descartes quadruple of bends (a,b,c,d), with a<0, find a |
---|
372 | * quadruple of circles, represented by (bend,bend*x,bend*y), such |
---|
373 | * that the circles have the given bends and the bends times the |
---|
374 | * centers are integers. |
---|
375 | * |
---|
376 | * This just performs an exaustive search, assuming that the outer |
---|
377 | * circle has center in the unit square. |
---|
378 | * |
---|
379 | * It is always sufficient to look in {(x,y):0<=y<=x<=1/2} for the |
---|
380 | * center of the outer circle, but this may not lead to a packing |
---|
381 | * that can be labelled with integer spherical and hyperbolic labels. |
---|
382 | * To effect the smaller search, replace FOR(a) with |
---|
383 | * |
---|
384 | * for (pa = ea/2; pa <= 0; pa++) for (qa = pa; qa <= 0; qa++) |
---|
385 | */ |
---|
386 | |
---|
387 | #define For(v,l,h) for (v = l; v <= h; v++) |
---|
388 | #define FOR(z) For(p##z,lop##z,hip##z) For(q##z,loq##z,hiq##z) |
---|
389 | #define H(z) ((e##z*e##z+p##z*p##z+q##z*q##z)%2) |
---|
390 | #define UNIT(z) ((abs(e##z)-1)*(abs(e##z)-1) >= p##z*p##z+q##z*q##z) |
---|
391 | #define T(z,w) is_tangent(e##z,p##z,q##z,e##w,p##w,q##w) |
---|
392 | #define LO(r,z) lo##r##z = iceil(e##z*(r##a+1),ea)-1 |
---|
393 | #define HI(r,z) hi##r##z = iflor(e##z*(r##a-1),ea)-1 |
---|
394 | #define B(z) LO(p,z); HI(p,z); LO(q,z); HI(q,z) |
---|
395 | |
---|
396 | static int |
---|
397 | is_quad(int a, int b, int c, int d) |
---|
398 | { |
---|
399 | int s; |
---|
400 | |
---|
401 | s = a+b+c+d; |
---|
402 | return 2*(a*a+b*b+c*c+d*d) == s*s; |
---|
403 | } |
---|
404 | |
---|
405 | static Bool |
---|
406 | is_tangent(int e1, int p1, int q1, int e2, int p2, int q2) |
---|
407 | { |
---|
408 | int dx, dy, s; |
---|
409 | |
---|
410 | dx = p1*e2 - p2*e1; |
---|
411 | dy = q1*e2 - q2*e1; |
---|
412 | s = e1 + e2; |
---|
413 | return dx*dx + dy*dy == s*s; |
---|
414 | } |
---|
415 | |
---|
416 | static int |
---|
417 | iflor(int a, int b) |
---|
418 | { |
---|
419 | int q; |
---|
420 | |
---|
421 | if (b == 0) { |
---|
422 | (void) printf("iflor: b = 0\n"); |
---|
423 | return 0; |
---|
424 | } |
---|
425 | if (a%b == 0) |
---|
426 | return a/b; |
---|
427 | q = abs(a)/abs(b); |
---|
428 | return ((a<0)^(b<0)) ? -q-1 : q; |
---|
429 | } |
---|
430 | |
---|
431 | static int |
---|
432 | iceil(int a, int b) |
---|
433 | { |
---|
434 | int q; |
---|
435 | |
---|
436 | if (b == 0) { |
---|
437 | (void) printf("iceil: b = 0\n"); |
---|
438 | return 0; |
---|
439 | } |
---|
440 | if (a%b == 0) |
---|
441 | return a/b; |
---|
442 | q = abs(a)/abs(b); |
---|
443 | return ((a<0)^(b<0)) ? -q : 1+q; |
---|
444 | } |
---|
445 | |
---|
446 | static double |
---|
447 | geom(int geometry, int e, int p, int q) |
---|
448 | { |
---|
449 | int g = (geometry == spherical) ? -1 : |
---|
450 | (geometry == hyperbolic) ? 1 : 0; |
---|
451 | |
---|
452 | if (g) |
---|
453 | return (e*e + (1.0 - p*p - q*q) * g) / (2.0*e); |
---|
454 | (void) printf("geom: g = 0\n"); |
---|
455 | return e; |
---|
456 | } |
---|
457 | |
---|
458 | static void |
---|
459 | cquad(circle *c1, circle *c2, circle *c3, circle *c4) |
---|
460 | { |
---|
461 | int ea, eb, ec, ed; |
---|
462 | int pa, pb, pc, pd; |
---|
463 | int qa, qb, qc, qd; |
---|
464 | int lopa, lopb, lopc, lopd; |
---|
465 | int hipa, hipb, hipc, hipd; |
---|
466 | int loqa, loqb, loqc, loqd; |
---|
467 | int hiqa, hiqb, hiqc, hiqd; |
---|
468 | |
---|
469 | ea = (int) c1->e; |
---|
470 | eb = (int) c2->e; |
---|
471 | ec = (int) c3->e; |
---|
472 | ed = (int) c4->e; |
---|
473 | if (ea >= 0) |
---|
474 | (void) printf("ea = %d\n", ea); |
---|
475 | if (!is_quad(ea,eb,ec,ed)) |
---|
476 | (void) printf("Error not quad %d %d %d %d\n", ea, eb, ec, ed); |
---|
477 | lopa = loqa = ea; |
---|
478 | hipa = hiqa = 0; |
---|
479 | FOR(a) { |
---|
480 | B(b); B(c); B(d); |
---|
481 | if (H(a) && UNIT(a)) FOR(b) { |
---|
482 | if (H(b) && T(a,b)) FOR(c) { |
---|
483 | if (H(c) && T(a,c) && T(b,c)) FOR(d) { |
---|
484 | if (H(d) && T(a,d) && T(b,d) && T(c,d)) { |
---|
485 | c1->s = geom(spherical, ea, pa, qa); |
---|
486 | c1->h = geom(hyperbolic, ea, pa, qa); |
---|
487 | c2->s = geom(spherical, eb, pb, qb); |
---|
488 | c2->h = geom(hyperbolic, eb, pb, qb); |
---|
489 | c3->s = geom(spherical, ec, pc, qc); |
---|
490 | c3->h = geom(hyperbolic, ec, pc, qc); |
---|
491 | c4->s = geom(spherical, ed, pd, qd); |
---|
492 | c4->h = geom(hyperbolic, ed, pd, qd); |
---|
493 | } |
---|
494 | } |
---|
495 | } |
---|
496 | } |
---|
497 | } |
---|
498 | } |
---|
499 | |
---|
500 | static void |
---|
501 | p(ModeInfo *mi, circle c) |
---|
502 | { |
---|
503 | apollonianstruct *cp = &apollonians[MI_SCREEN(mi)]; |
---|
504 | char string[10]; |
---|
505 | double g, e; |
---|
506 | int g_width; |
---|
507 | |
---|
508 | #ifdef DEBUG |
---|
509 | (void) printf("c.e=%g c.s=%g c.h=%g c.x=%g c.y=%g\n", |
---|
510 | c.e, c.s, c.h, c.x, c.y); |
---|
511 | #endif |
---|
512 | g = (cp->geometry == spherical) ? c.s : (cp->geometry == hyperbolic) ? |
---|
513 | c.h : c.e; |
---|
514 | if (c.e < 0.0) { |
---|
515 | if (g < 0.0) |
---|
516 | g = -g; |
---|
517 | if (MI_NPIXELS(mi) <= 2) |
---|
518 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), |
---|
519 | MI_WHITE_PIXEL(mi)); |
---|
520 | else |
---|
521 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), |
---|
522 | MI_PIXEL(mi, ((int) ((g + cp->color_offset) * |
---|
523 | g)) % MI_NPIXELS(mi))); |
---|
524 | XDrawArc(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
525 | ((int) (cp->size * (-cp->c1.e) * (c.x - 1.0) / |
---|
526 | (-2.0 * c.e) + cp->size / 2.0 + cp->offset.x)), |
---|
527 | ((int) (cp->size * (-cp->c1.e) * (c.y - 1.0) / |
---|
528 | (-2.0 * c.e) + cp->size / 2.0 + cp->offset.y)), |
---|
529 | (int) (cp->c1.e * cp->size / c.e), |
---|
530 | (int) (cp->c1.e * cp->size / c.e), 0, 23040); |
---|
531 | if (!cp->label) { |
---|
532 | #ifdef DEBUG |
---|
533 | (void) printf("%g\n", -g); |
---|
534 | #endif |
---|
535 | return; |
---|
536 | } |
---|
537 | (void) sprintf(string, "%g", (g == 0.0) ? 0 : -g); |
---|
538 | if (cp->size >= 10 * FONT_WIDTH) { |
---|
539 | /* hard code these to corners */ |
---|
540 | XDrawString(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
541 | ((int) (cp->size * c.x / (2.0 * c.e))) + cp->offset.x, |
---|
542 | ((int) (cp->size * c.y / (2.0 * c.e))) + FONT_HEIGHT, |
---|
543 | string, (g == 0.0) ? 1 : ((g < 10.0) ? 2 : |
---|
544 | ((g < 100.0) ? 3 : 4))); |
---|
545 | } |
---|
546 | if (cp->altgeom && MI_HEIGHT(mi) >= 30 * FONT_WIDTH) { |
---|
547 | XDrawString(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
548 | ((int) (cp->size * c.x / (2.0 * c.e) + cp->offset.x)), |
---|
549 | ((int) (cp->size * c.y / (2.0 * c.e) + MI_HEIGHT(mi) - |
---|
550 | FONT_HEIGHT / 2)), space_string[cp->geometry], |
---|
551 | strlen(space_string[cp->geometry])); |
---|
552 | } |
---|
553 | return; |
---|
554 | } |
---|
555 | if (MI_NPIXELS(mi) <= 2) |
---|
556 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), MI_WHITE_PIXEL(mi)); |
---|
557 | else |
---|
558 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), |
---|
559 | MI_PIXEL(mi, ((int) ((g + cp->color_offset) * g)) % |
---|
560 | MI_NPIXELS(mi))); |
---|
561 | if (c.e == 0.0) { |
---|
562 | if (c.x == 0.0 && c.y != 0.0) { |
---|
563 | XDrawLine(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
564 | 0, (int) ((c.y + 1.0) * cp->size / 2.0 + cp->offset.y), |
---|
565 | MI_WIDTH(mi), |
---|
566 | (int) ((c.y + 1.0) * cp->size / 2.0 + cp->offset.y)); |
---|
567 | } else if (c.y == 0.0 && c.x != 0.0) { |
---|
568 | XDrawLine(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
569 | (int) ((c.x + 1.0) * cp->size / 2.0 + cp->offset.x), 0, |
---|
570 | (int) ((c.x + 1.0) * cp->size / 2.0 + cp->offset.x), |
---|
571 | MI_HEIGHT(mi)); |
---|
572 | } |
---|
573 | return; |
---|
574 | } |
---|
575 | e = (cp->c1.e >= 0.0) ? 1.0 : -cp->c1.e; |
---|
576 | XFillArc(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
577 | ((int) (cp->size * e * (c.x - 1.0) / (2.0 * c.e) + |
---|
578 | cp->size / 2.0 + cp->offset.x)), |
---|
579 | ((int) (cp->size * e * (c.y - 1.0) / (2.0 * c.e) + |
---|
580 | cp->size / 2.0 + cp->offset.y)), |
---|
581 | (int) (e * cp->size / c.e), (int) (e * cp->size / c.e), |
---|
582 | 0, 23040); |
---|
583 | if (!cp->label) { |
---|
584 | #ifdef DEBUG |
---|
585 | (void) printf("%g\n", g); |
---|
586 | #endif |
---|
587 | return; |
---|
588 | } |
---|
589 | if (MI_NPIXELS(mi) <= 2) |
---|
590 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), MI_BLACK_PIXEL(mi)); |
---|
591 | else |
---|
592 | XSetForeground(MI_DISPLAY(mi), MI_GC(mi), |
---|
593 | MI_PIXEL(mi, ((int) ((g + cp->color_offset) * g) + |
---|
594 | MI_NPIXELS(mi) / 2) % MI_NPIXELS(mi))); |
---|
595 | g_width = (g < 10.0) ? 1: ((g < 100.0) ? 2 : 3); |
---|
596 | if (c.e < e * cp->size / (FONT_LENGTH + 5 * g_width) && g < 1000.0) { |
---|
597 | (void) sprintf(string, "%g", g); |
---|
598 | XDrawString(MI_DISPLAY(mi), MI_WINDOW(mi), MI_GC(mi), |
---|
599 | ((int) (cp->size * e * c.x / (2.0 * c.e) + |
---|
600 | cp->size / 2.0 + cp->offset.x)) - |
---|
601 | g_width * FONT_WIDTH / 2, |
---|
602 | ((int) (cp->size * e * c.y / (2.0 * c.e) + |
---|
603 | cp->size / 2.0 + cp->offset.y)) + |
---|
604 | FONT_HEIGHT / 2, |
---|
605 | string, g_width); |
---|
606 | } |
---|
607 | } |
---|
608 | |
---|
609 | #define BIG 7 |
---|
610 | static void |
---|
611 | f(ModeInfo *mi, circle c1, circle c2, circle c3, circle c4) |
---|
612 | { |
---|
613 | apollonianstruct *cp = &apollonians[MI_SCREEN(mi)]; |
---|
614 | int e = (int) ((cp->c1.e >= 0.0) ? 1.0 : -cp->c1.e); |
---|
615 | circle c; |
---|
616 | |
---|
617 | c.e = 2*(c1.e+c2.e+c3.e) - c4.e; |
---|
618 | c.s = 2*(c1.s+c2.s+c3.s) - c4.s; |
---|
619 | c.h = 2*(c1.h+c2.h+c3.h) - c4.h; |
---|
620 | c.x = 2*(c1.x+c2.x+c3.x) - c4.x; |
---|
621 | c.y = 2*(c1.y+c2.y+c3.y) - c4.y; |
---|
622 | if (c.e > cp->size * e || c.x / c.e > BIG || c.y / c.e > BIG || |
---|
623 | c.x / c.e < -BIG || c.y / c.e < -BIG) |
---|
624 | return; |
---|
625 | p(mi, c); |
---|
626 | f(mi, c2, c3, c, c1); |
---|
627 | f(mi, c1, c3, c, c2); |
---|
628 | f(mi, c1, c2, c, c3); |
---|
629 | } |
---|
630 | |
---|
631 | static void |
---|
632 | free_apollonian(Display *display, apollonianstruct *cp) |
---|
633 | { |
---|
634 | if (cp->quad != NULL) { |
---|
635 | (void) free((void *) cp->quad); |
---|
636 | cp->quad = (apollonian_quadruple *) NULL; |
---|
637 | } |
---|
638 | #ifdef DOFONT |
---|
639 | if (cp->gc != None) { |
---|
640 | XFreeGC(display, cp->gc); |
---|
641 | cp->gc = None; |
---|
642 | } |
---|
643 | if (cp->font != None) { |
---|
644 | XFreeFont(display, cp->font); |
---|
645 | cp->font = None; |
---|
646 | } |
---|
647 | #endif |
---|
648 | } |
---|
649 | |
---|
650 | #ifndef DEBUG |
---|
651 | void |
---|
652 | randomize_c(int randomize, circle * c) |
---|
653 | { |
---|
654 | if (randomize / 2) { |
---|
655 | double temp; |
---|
656 | |
---|
657 | temp = c->x; |
---|
658 | c->x = c->y; |
---|
659 | c->y = temp; |
---|
660 | } |
---|
661 | if (randomize % 2) { |
---|
662 | c->x = -c->x; |
---|
663 | c->y = -c->y; |
---|
664 | } |
---|
665 | } |
---|
666 | #endif |
---|
667 | |
---|
668 | void |
---|
669 | init_apollonian(ModeInfo * mi) |
---|
670 | { |
---|
671 | apollonianstruct *cp; |
---|
672 | int i; |
---|
673 | |
---|
674 | if (apollonians == NULL) { |
---|
675 | if ((apollonians = (apollonianstruct *) calloc(MI_NUM_SCREENS(mi), |
---|
676 | sizeof (apollonianstruct))) == NULL) |
---|
677 | return; |
---|
678 | } |
---|
679 | cp = &apollonians[MI_SCREEN(mi)]; |
---|
680 | |
---|
681 | cp->size = MAX(MIN(MI_WIDTH(mi), MI_HEIGHT(mi)) - 1, 1); |
---|
682 | cp->offset.x = (MI_WIDTH(mi) - cp->size) / 2; |
---|
683 | cp->offset.y = (MI_HEIGHT(mi) - cp->size) / 2; |
---|
684 | cp->color_offset = NRAND(MI_NPIXELS(mi)); |
---|
685 | |
---|
686 | #ifdef DOFONT |
---|
687 | if (cp->font == None) { |
---|
688 | if ((cp->font = getFont(MI_DISPLAY(mi))) == None) |
---|
689 | return False; |
---|
690 | } |
---|
691 | #endif |
---|
692 | cp->label = label; |
---|
693 | cp->altgeom = cp->label && altgeom; |
---|
694 | |
---|
695 | if (cp->quad == NULL) { |
---|
696 | cp->count = ABS(MI_COUNT(mi)); |
---|
697 | if ((cp->quad = (apollonian_quadruple *) malloc(cp->count * |
---|
698 | sizeof (apollonian_quadruple))) == NULL) { |
---|
699 | return; |
---|
700 | } |
---|
701 | dquad(cp->count, cp->quad); |
---|
702 | } |
---|
703 | cp->game = NRAND(PREDEF_CIRCLE_GAMES + cp->count); |
---|
704 | cp->geometry = (cp->game && cp->altgeom) ? NRAND(3) : 0; |
---|
705 | |
---|
706 | if (cp->game < PREDEF_CIRCLE_GAMES) { |
---|
707 | cp->c1 = examples[cp->game][0]; |
---|
708 | cp->c2 = examples[cp->game][1]; |
---|
709 | cp->c3 = examples[cp->game][2]; |
---|
710 | cp->c4 = examples[cp->game][3]; |
---|
711 | /* do not label non int */ |
---|
712 | cp->label = cp->label && (cp->c4.e == (int) cp->c4.e); |
---|
713 | } else { /* uses results of dquad, all int */ |
---|
714 | i = cp->game - PREDEF_CIRCLE_GAMES; |
---|
715 | cp->c1.e = cp->quad[i].a; |
---|
716 | cp->c2.e = cp->quad[i].b; |
---|
717 | cp->c3.e = cp->quad[i].c; |
---|
718 | cp->c4.e = cp->quad[i].d; |
---|
719 | if (cp->geometry) |
---|
720 | cquad(&(cp->c1), &(cp->c2), &(cp->c3), &(cp->c4)); |
---|
721 | } |
---|
722 | cp->time = 0; |
---|
723 | MI_CLEARWINDOW(mi); |
---|
724 | if (cp->game != 0) { |
---|
725 | double q123; |
---|
726 | |
---|
727 | if (cp->c1.e == 0.0 || cp->c1.e == -cp->c2.e) |
---|
728 | return; |
---|
729 | cp->c1.x = 0.0; |
---|
730 | cp->c1.y = 0.0; |
---|
731 | cp->c2.x = -(cp->c1.e + cp->c2.e) / cp->c1.e; |
---|
732 | cp->c2.y = 0; |
---|
733 | q123 = sqrt(cp->c1.e * cp->c2.e + cp->c1.e * cp->c3.e + |
---|
734 | cp->c2.e * cp->c3.e); |
---|
735 | #ifdef DEBUG |
---|
736 | (void) printf("q123 = %g, ", q123); |
---|
737 | #endif |
---|
738 | cp->c3.x = (cp->c1.e * cp->c1.e - q123 * q123) / (cp->c1.e * |
---|
739 | (cp->c1.e + cp->c2.e)); |
---|
740 | cp->c3.y = -2.0 * q123 / (cp->c1.e + cp->c2.e); |
---|
741 | q123 = -cp->c1.e - cp->c2.e + q123; |
---|
742 | cp->c4.x = (cp->c1.e * cp->c1.e - q123 * q123) / (cp->c1.e * |
---|
743 | (cp->c1.e + cp->c2.e)); |
---|
744 | cp->c4.y = -2.0 * q123 / (cp->c1.e + cp->c2.e); |
---|
745 | #ifdef DEBUG |
---|
746 | (void) printf("q124 = %g\n", q123); |
---|
747 | (void) printf("%g %g %g %g %g %g %g %g\n", |
---|
748 | cp->c1.x, cp->c1.y, cp->c2.x, cp->c2.y, |
---|
749 | cp->c3.x, cp->c3.y, cp->c4.x, cp->c4.y); |
---|
750 | #endif |
---|
751 | } |
---|
752 | #ifndef DEBUG |
---|
753 | if (LRAND() & 1) { |
---|
754 | cp->c3.y = -cp->c3.y; |
---|
755 | cp->c4.y = -cp->c4.y; |
---|
756 | } |
---|
757 | i = NRAND(4); |
---|
758 | randomize_c(i, &(cp->c1)); |
---|
759 | randomize_c(i, &(cp->c2)); |
---|
760 | randomize_c(i, &(cp->c3)); |
---|
761 | randomize_c(i, &(cp->c4)); |
---|
762 | #endif |
---|
763 | } |
---|
764 | |
---|
765 | void |
---|
766 | draw_apollonian(ModeInfo * mi) |
---|
767 | { |
---|
768 | apollonianstruct *cp; |
---|
769 | |
---|
770 | if (apollonians == NULL) |
---|
771 | return; |
---|
772 | cp = &apollonians[MI_SCREEN(mi)]; |
---|
773 | |
---|
774 | |
---|
775 | MI_IS_DRAWN(mi) = True; |
---|
776 | |
---|
777 | if (cp->time < 5) { |
---|
778 | switch (cp->time) { |
---|
779 | case 0: |
---|
780 | p(mi, cp->c1); |
---|
781 | p(mi, cp->c2); |
---|
782 | p(mi, cp->c3); |
---|
783 | p(mi, cp->c4); |
---|
784 | break; |
---|
785 | case 1: |
---|
786 | f(mi, cp->c1, cp->c2, cp->c3, cp->c4); |
---|
787 | break; |
---|
788 | case 2: |
---|
789 | f(mi, cp->c1, cp->c2, cp->c4, cp->c3); |
---|
790 | break; |
---|
791 | case 3: |
---|
792 | f(mi, cp->c1, cp->c3, cp->c4, cp->c2); |
---|
793 | break; |
---|
794 | case 4: |
---|
795 | f(mi, cp->c2, cp->c3, cp->c4, cp->c1); |
---|
796 | } |
---|
797 | } |
---|
798 | if (++cp->time > MI_CYCLES(mi)) |
---|
799 | { |
---|
800 | #ifdef STANDALONE |
---|
801 | erase_full_window(MI_DISPLAY(mi), MI_WINDOW(mi)); |
---|
802 | #endif /* STANDALONE */ |
---|
803 | init_apollonian(mi); |
---|
804 | } |
---|
805 | } |
---|
806 | |
---|
807 | void |
---|
808 | release_apollonian(ModeInfo * mi) |
---|
809 | { |
---|
810 | if (apollonians != NULL) { |
---|
811 | int screen; |
---|
812 | |
---|
813 | for (screen = 0; screen < MI_NUM_SCREENS(mi); screen++) |
---|
814 | free_apollonian(MI_DISPLAY(mi), &apollonians[screen]); |
---|
815 | (void) free((void *) apollonians); |
---|
816 | apollonians = (apollonianstruct *) NULL; |
---|
817 | } |
---|
818 | } |
---|
819 | |
---|
820 | #endif /* MODE_apollonian */ |
---|