1
2 /* Hello reader!
3
4 * Are you sure you want read this? Its very cryptic and strange code. YOU
5 * HAVE BEEN WARNED! Its purpose is to generate as fast as possible
6 * calculation loops for various formulas/algorithms. It uses lots of
7 * coprocesor magic. It is included from formulas.c
8 */
9
10 #ifndef VARIABLES /*supply defaultd values */
11 #define VARIABLES
12 #endif
13 #ifndef PRETEST
14 #define PRETEST 0
15 #endif
16 #ifndef INIT
17 #define INIT
18 #endif
19 #ifndef POSTCALC
20 #define POSTCALC
21 #endif
22 #ifndef PRESMOOTH
23 #define PRESMOOTH zre = rp + ip
24 #endif
25 #ifndef UFORMULA
26 #define UFORMULA FORMULA
27 #endif
28 #ifndef UEND
29 #define UEND
30 #endif
31 #ifndef SAVE
32 #define SAVE
33 #endif
34 #ifndef SAVEVARIABLES
35 #define SAVEVARIABLES
36 #endif
37 #ifndef RESTORE
38 #define RESTORE
39 #endif
40 #ifndef RANGE
41 #define RANGE 2
42 #endif
43
44 /* Prepare main loop */
45 #ifndef NSFORMULALOOP
46 #define NSFORMULALOOP(iter) \
47 do { /*try first 8 iterations */ \
48 FORMULA; \
49 iter--; \
50 } while (BTEST && iter)
51 #endif
52 #ifndef SFORMULALOOP
53 #define SFORMULALOOP(iter) \
54 do { /*try first 8 iterations */ \
55 SAVEZMAG; \
56 FORMULA; \
57 iter--; \
58 } while (BTEST && iter)
59 #endif
60 #ifndef FORMULALOOP
61 #ifdef SMOOTHMODE
62 #define FORMULALOOP SFORMULALOOP
63 #else
64 #define FORMULALOOP NSFORMULALOOP
65 #endif
66 #endif
67
68 #ifdef SMOOTHMODE
69 #ifdef CUSTOMSAVEZMAG
70 #define SAVEZMAG CUSTOMSAVEZMAG;
71 #else
72 #define SAVEZMAG szmag = rp + ip;
73 #endif
74 #else
75 #define SAVEZMAG
76 #endif
77
78 #ifdef UNCOMPRESS
79 /*uncompressed version of loop */
80 #ifdef SMOOTHMODE
81 static unsigned int SCALC(number_t zre, number_t zim, number_t pre,
82 number_t pim);
83
SCALC(number_t zre,number_t zim,number_t pre,number_t pim)84 static unsigned int SCALC(number_t zre, number_t zim, number_t pre,
85 number_t pim)
86 #else
87 static unsigned int CALC(number_t zre, number_t zim, number_t pre,
88 number_t pim);
89
90 static unsigned CALC(number_t zre, number_t zim, number_t pre, number_t pim)
91 #endif
92 {
93 unsigned int iter = cfractalc.maxiter;
94 number_t szre = 0, szim = 0;
95 #ifdef RPIP
96 number_t rp = 0, ip;
97 #endif
98 #ifdef SMOOTHMODE
99 number_t szmag = 0;
100 #endif
101 SAVEVARIABLES VARIABLES;
102 INIT;
103 if (PRETEST)
104 iter = 0;
105 else {
106 #ifdef RPIP
107 rp = zre * zre;
108 ip = zim * zim;
109 #endif
110 if (iter < 16) {
111
112 /*try first 8 iterations */
113 if (BTEST && iter) {
114 FORMULALOOP(iter);
115 }
116 /*
117 while (BTEST && iter)
118 {
119 SAVEZMAG;
120 FORMULA;
121 iter--;
122 } */
123 } else {
124 iter = 8 + (cfractalc.maxiter & 7);
125
126 /*try first 8 iterations */
127 if (BTEST && iter) {
128 FORMULALOOP(iter);
129 }
130 /*
131 while (BTEST && iter)
132 {
133 SAVEZMAG;
134 FORMULA;
135 iter--;
136 } */
137 if (BTEST) {
138 iter = (cfractalc.maxiter - 8) & (~7);
139 iter >>= 3;
140 /*do next 8 iteration w/o out of bounds checking */
141 do {
142 /*hmm..we are probably in some deep area. */
143 szre = zre; /*save current position */
144 szim = zim;
145 SAVE;
146 UFORMULA;
147 UFORMULA;
148 UFORMULA;
149 UFORMULA;
150 UFORMULA;
151 UFORMULA;
152 UFORMULA;
153 UFORMULA;
154 UEND;
155 iter--;
156 } while (BTEST && iter);
157 if (!(BTEST)) { /*we got out of bounds */
158 iter <<= 3;
159 iter += 8; /*restore saved position */
160 RESTORE;
161 zre = szre;
162 zim = szim;
163 #ifdef RPIP
164 rp = zre * zre;
165 ip = zim * zim;
166 #endif
167 FORMULALOOP(iter);
168 /*
169 do
170 {
171 SAVEZMAG
172 FORMULA;
173 iter--;
174 }
175 while (BTEST && iter); */
176 }
177 } else
178 iter += cfractalc.maxiter - 8 - (cfractalc.maxiter & 7);
179 }
180 }
181 #ifdef SMOOTHMODE
182 if (iter)
183 SMOOTHOUTPUT();
184 POSTCALC;
185 iter = cfractalc.maxiter - iter;
186 INOUTPUT();
187 #else
188 POSTCALC;
189 iter = cfractalc.maxiter - iter;
190 OUTPUT();
191 #endif
192 }
193 #else
194 #ifdef SMOOTHMODE
195 static unsigned int SCALC(number_t zre, number_t zim, number_t pre,
196 number_t pim);
197
SCALC(number_t zre,number_t zim,number_t pre,number_t pim)198 static unsigned int SCALC(number_t zre, number_t zim, number_t pre,
199 number_t pim)
200 #else
201 static unsigned int CALC(number_t zre, number_t zim, number_t pre,
202 number_t pim);
203
204 static unsigned int CALC(number_t zre, number_t zim, number_t pre, number_t pim)
205 #endif
206 {
207 unsigned int iter = cfractalc.maxiter /*& (~(int) 3) */;
208 #ifdef RPIP
209 number_t rp, ip;
210 #endif
211 #ifdef SMOOTHMODE
212 number_t szmag = 0;
213 #endif
214 VARIABLES;
215 INIT;
216 if (PRETEST)
217 iter = 0;
218 else {
219 #ifdef RPIP
220 rp = zre * zre;
221 ip = zim * zim;
222 #endif
223 if (BTEST && iter) {
224 FORMULALOOP(iter);
225 }
226 /*
227 while (BTEST && iter)
228 {
229 SAVEZMAG
230 FORMULA;
231 iter--;
232
233 } */
234 }
235 #ifdef SMOOTHMODE
236 if (iter)
237 SMOOTHOUTPUT();
238 POSTCALC;
239 iter = cfractalc.maxiter - iter;
240 INOUTPUT();
241 #else
242 POSTCALC;
243 iter = cfractalc.maxiter - iter;
244 OUTPUT();
245 #endif
246 }
247 #endif
248
249 /*F. : Periodicity checking routines. (16-02-97)
250 All comments preceded by F. are mine (Fabrice Premel premelfa@etu.utc.fr).
251 Tried to make code as efficient as possible.
252 Next to do is convert lim in a variable that would be updated sometimes
253 I'll try to make here a short explanation on periodicity checking :
254 first, we'll define 2 variables : whentosave and whenincsave, which are,
255 respectively, a measure of when we have to update saved values to be checked,
256 and when to increase interval between 2 updates, as if they're too close,
257 we'll miss large periods. We save Z at the beginning, and then we compare
258 each new iteration with this Z, and if naerly equal, we declare the suite to
259 be periodic. When ( iter mod whentosave ) == 0, we store a new value, and we
260 repeat.
261
262 UNCOMPRESSed form is just an extension, with careful that if we only check
263 whentosave all 8 iterations, number of iterations must be well set at the
264 beginning.This is done by adding a (iter&7) in the while statement preceding
265 then uncompressed calculation. */
266
267 /*F. : This is from then lim factor that depends all periodicity check spped :
268 the bigger it is, the faster we can detect periodicity, but the bigger it is,
269 the more we can introduce errors. I suggest a value of
270 (maxx-minx)/(double)getmaxx for a classic Mandelbrot Set, and maybe a lesser
271 value for an extra power Mandelbrot. But this should be calculated outer
272 from here (ie each frame, for example), to avoid new calculus */
273 #ifdef PERI
274 #define PCHECK \
275 (abs_less_than(r1 - zre, cfractalc.periodicity_limit) && \
276 abs_less_than(s1 - zim, cfractalc.periodicity_limit))
277
278 #ifndef UNCOMPRESS
279
280 #ifdef SMOOTHMODE
281 static unsigned int SPERI(number_t zre, number_t zim, number_t pre,
282 number_t pim);
283
SPERI(number_t zre,number_t zim,number_t pre,number_t pim)284 static unsigned int SPERI(number_t zre, number_t zim, number_t pre,
285 number_t pim)
286 #else
287 static unsigned int PERI(number_t zre, number_t zim, number_t pre,
288 number_t pim);
289
290 static unsigned int PERI(number_t zre, number_t zim, number_t pre, number_t pim)
291 #endif
292 {
293 unsigned int iter = cfractalc.maxiter /*& (~(int) 3) */, iter1 = 8;
294 number_t r1, s1;
295 int whensavenew, whenincsave;
296 #ifdef RPIP
297 number_t rp, ip;
298 #endif
299 #ifdef SMOOTHMODE
300 number_t szmag = 0;
301 #endif
302 VARIABLES;
303 INIT;
304 if (PRETEST)
305 iter = 0;
306 else {
307 #ifdef RPIP
308 rp = zre * zre;
309 ip = zim * zim;
310 #endif
311 if (iter < iter1)
312 iter1 = iter, iter = 8;
313
314 /*H. : do first few iterations w/o checking */
315 if (BTEST && iter1) {
316 FORMULALOOP(iter1);
317 }
318 /*
319 while (BTEST && iter1)
320 {
321 SAVEZMAG;
322 FORMULA;
323 iter1--;
324 } */
325 if (iter1) {
326 if (iter >= 8)
327 iter -= 8 - iter1;
328 goto end;
329 }
330 if (iter <= 8) {
331 iter = iter1;
332 } else {
333 iter -= 8;
334 r1 = zre;
335 s1 = zim;
336 whensavenew = 3; /*You should adapt these values */
337 /*F. : We should always define whensavenew as 2^N-1, so we could use
338 * a AND instead of % */
339
340 whenincsave = 10;
341 /*F. : problem is that after deep zooming, peiodicity is never
342 detected early, cause is is quite slow before going in a periodic
343 loop. So, we should start checking periodicity only after some
344 times */
345 while (BTEST && iter) {
346 SAVEZMAG;
347 FORMULA;
348 if ((iter & whensavenew) == 0) { /*F. : changed % to & */
349 r1 = zre;
350 s1 = zim;
351 whenincsave--;
352 if (!whenincsave) {
353 whensavenew =
354 ((whensavenew + 1) << 1) -
355 1; /*F. : Changed to define a new AND mask */
356 whenincsave = 10;
357 }
358 } else {
359 if (PCHECK) {
360 PERIINOUTPUT();
361 }
362 }
363 iter--;
364 }
365 }
366 }
367 end:
368 #ifdef SMOOTHMODE
369 if (iter)
370 SMOOTHOUTPUT();
371 POSTCALC;
372 iter = cfractalc.maxiter - iter;
373 INOUTPUT();
374 #else
375 POSTCALC;
376 iter = cfractalc.maxiter - iter;
377 OUTPUT();
378 #endif
379 }
380
381 #else
382
383 /*F. : UNCOMPRESSed version. Note that whensavenew+1 should be a multiple of 8,
384 else periodicity won't be able to detect anything. */
385 /*F. : this macros definitions are really strange, but after a while, it's good
386 */
387
388 #ifdef SMOOTHMODE
389 static unsigned int SPERI(number_t zre, number_t zim, number_t pre,
390 number_t pim);
391
SPERI(number_t zre,number_t zim,number_t pre,number_t pim)392 static unsigned int SPERI(number_t zre, number_t zim, number_t pre,
393 number_t pim)
394 #else
395 static unsigned int PERI(number_t zre, number_t zim, number_t pre,
396 number_t pim);
397
398 static unsigned int PERI(number_t zre, number_t zim, number_t pre, number_t pim)
399 #endif
400 {
401 unsigned int iter = cfractalc.maxiter /*& (~(int) 3) */;
402 number_t r1 = zre, s1 = zim;
403 number_t szre = 0,
404 szim =
405 0; /*F. : Didn't declared register, cause they are few used */
406 unsigned int whensavenew, whenincsave;
407 #ifdef RPIP
408 number_t rp = 0, ip;
409 #endif
410 #ifdef SMOOTHMODE
411 number_t szmag = 0;
412 #endif
413 SAVEVARIABLES VARIABLES;
414 INIT;
415 if (PRETEST)
416 iter = 0;
417 else {
418 if (cfractalc.maxiter <= 16) {
419 #ifdef RPIP
420 rp = zre * zre;
421 ip = zim * zim;
422 #endif
423 /*F. : Added iter&7 to be sure we'll be on a 8 multiple */
424 if (BTEST && iter) {
425 FORMULALOOP(iter);
426 }
427 /*
428 while (BTEST && iter)
429 {
430 SAVEZMAG
431 FORMULA;
432 iter--;
433 } */
434 } else {
435 whensavenew = 7; /*You should adapt these values */
436 /*F. : We should always define whensavenew as 2^N-1, so we could use
437 * a AND instead of % */
438
439 whenincsave = 10;
440 #ifdef RPIP
441 rp = zre * zre;
442 ip = zim * zim;
443 #endif
444 /*F. : problem is that after deep zooming, peiodicity is never
445 detected early, cause is is quite slow before going in a periodic
446 loop. So, we should start checking periodicity only after some
447 times */
448 iter = 8 + (cfractalc.maxiter & 7);
449 while (BTEST && iter) { /*F. : Added iter&7 to be sure we'll be on a
450 8 multiple */
451 SAVEZMAG FORMULA;
452 iter--;
453 }
454 if (BTEST) { /*F. : BTEST is calculed two times here, isn't it ? */
455 /*H. : No gcc is clever and adds test to the end :) */
456 iter = (cfractalc.maxiter - 8) & (~7);
457 do {
458 szre = zre, szim = zim;
459 SAVE;
460 SAVEZMAG
461 FORMULA; /*F. : Calculate one time */
462 if (PCHECK)
463 goto periodicity;
464 FORMULA;
465 if (PCHECK)
466 goto periodicity;
467 FORMULA;
468 if (PCHECK)
469 goto periodicity;
470 FORMULA;
471 if (PCHECK)
472 goto periodicity;
473 FORMULA;
474 if (PCHECK)
475 goto periodicity;
476 FORMULA;
477 if (PCHECK)
478 goto periodicity;
479 FORMULA;
480 if (PCHECK)
481 goto periodicity;
482 FORMULA;
483 if (PCHECK)
484 goto periodicity;
485 iter -= 8;
486 /*F. : We only test this now, as it can't be true before */
487 if ((iter & whensavenew) == 0) { /*F. : changed % to & */
488 r1 = zre, s1 = zim; /*F. : Save new values */
489 whenincsave--;
490 if (!whenincsave) {
491 whensavenew =
492 ((whensavenew + 1) << 1) -
493 1; /*F. : Changed to define a new AND mask */
494 whenincsave = 10; /*F. : Start back */
495 }
496 }
497 } while (BTEST && iter);
498 if (!BTEST) { /*we got out of bounds */
499 iter += 8; /*restore saved position */
500 RESTORE;
501 zre = szre;
502 zim = szim;
503 #ifdef RPIP
504 rp = zre * zre;
505 ip = zim * zim;
506 #endif
507 FORMULALOOP(iter);
508 /*
509 do
510 {
511 SAVEZMAG
512 FORMULA;
513 iter--;
514 }
515 while (BTEST && iter); */
516 }
517 } else
518 iter += cfractalc.maxiter - 8 - (cfractalc.maxiter & 7);
519 }
520 }
521 #ifdef SMOOTHMODE
522 if (iter)
523 SMOOTHOUTPUT();
524 POSTCALC;
525 iter = cfractalc.maxiter - iter;
526 INOUTPUT();
527 #else
528 POSTCALC;
529 iter = cfractalc.maxiter - iter;
530 OUTPUT();
531 #endif
532 periodicity:
533 PERIINOUTPUT();
534 }
535
536 /*else uncompress */
537 #endif
538
539 /*endif PERI */
540 #undef PCHECK
541 #endif
542
543 #ifndef SMOOTHMODE
544 #ifdef JULIA
JULIA(struct image * image,number_t pre,number_t pim)545 static void JULIA(struct image *image, number_t pre, number_t pim)
546 {
547 int i, i1, i2, j, x, y;
548 unsigned char iter, itmp2, itmp;
549 number_t rp = 0, ip = 0;
550 number_t zre, zim, im, xdelta, ydelta, range, rangep;
551 number_t xstep, ystep;
552 unsigned char *queue[QMAX];
553 unsigned char **qptr;
554 unsigned char *addr, **addr1 = image->currlines;
555 #ifdef STATISTICS
556 int guessed = 0, unguessed = 0, iters = 0;
557 #endif
558 VARIABLES;
559 range = (number_t)RANGE;
560 rangep = range * range;
561
562 xdelta = image->width / (RMAX - RMIN);
563 ydelta = image->height / (IMAX - IMIN);
564 xstep = (RMAX - RMIN) / image->width;
565 ystep = (IMAX - IMIN) / image->height;
566 init_julia(image, rangep, range, xdelta, ystep);
567 for (i2 = 0; i2 < 2; i2++)
568 for (i1 = 0; i1 < image->height; i1++) {
569 if (i1 % 2)
570 i = image->height / 2 - i1 / 2;
571 else
572 i = image->height / 2 + i1 / 2 + 1;
573 if (i >= image->height)
574 continue;
575 im = IMIN + (i + 0.5) * ystep;
576 for (j = (i + i2) & 1; j < image->width; j += 2) {
577 STAT(total2++);
578 addr = addr1[i] + j;
579 if (*addr != NOT_CALCULATED)
580 continue;
581 x = j;
582 y = i;
583 if (y > 0 && y < image->height - 1 && *(addr + 1) && x > 0 &&
584 x < image->width - 1) {
585 if ((iter = *(addr + 1)) != NOT_CALCULATED &&
586 iter == *(addr - 1) && iter == addr1[y - 1][x] &&
587 iter == addr1[y + 1][x]) {
588 *addr = *(addr + 1);
589 continue;
590 }
591 }
592 zim = im;
593 zre = RMIN + (j + 0.5) * xstep;
594 iter = (unsigned char)0;
595 qptr = queue;
596 ip = (zim * zim);
597 rp = (zre * zre);
598 INIT;
599 while (1) {
600 if (*addr != NOT_CALCULATED
601 #ifdef SAG
602 && (*addr == INPROCESS ||
603 (*addr != (unsigned char)1 &&
604 (itmp2 = *(addr + 1)) != NOT_CALCULATED &&
605 ((itmp2 != (itmp = *(addr - 1)) &&
606 itmp != NOT_CALCULATED) ||
607 (itmp2 != (itmp = *((addr1[y + 1]) + x)) &&
608 itmp != NOT_CALCULATED) ||
609 (itmp2 != (itmp = *((addr1[y - 1]) + x)) &&
610 itmp != NOT_CALCULATED))))
611 #endif
612 ) {
613 if (*addr == INPROCESS || *addr == INSET) {
614 *qptr = addr;
615 qptr++;
616 STAT(guessed++);
617 goto inset;
618 }
619 STAT(guessed++);
620 iter = *addr;
621 goto outset;
622 }
623 #ifdef STATISTICS
624 if (*addr != NOT_CALCULATED)
625 unguessed++;
626 #endif
627 if (*addr != INPROCESS) {
628 *qptr = addr;
629 qptr++;
630 *addr = INPROCESS;
631 if (qptr >= queue + QMAX)
632 goto inset;
633 }
634 STAT(iters++);
635 FORMULA;
636 ip = (zim * zim);
637 rp = (zre * zre);
638 if (greater_than(rp + ip, RANGE) || !(BTEST))
639 goto outset;
640 x = (int)((zre - RMIN) * xdelta);
641 y = (int)((zim - IMIN) * ydelta);
642 addr = addr1[y] + x;
643 if ((itmp = *(addr + 1)) != NOT_CALCULATED &&
644 itmp == *(addr - 1) && itmp == addr1[y - 1][x] &&
645 itmp == addr1[y + 1][x]) {
646 *addr = *(addr + 1);
647 }
648 }
649 inset:
650 while (qptr > queue) {
651 qptr--;
652 **qptr = INSET;
653 }
654 continue;
655 outset:
656 y = image->palette->size;
657 while (qptr > queue) {
658 qptr--;
659 iter++;
660 if ((int)iter >= y)
661 iter = (unsigned char)1;
662 **qptr = iter;
663 }
664 }
665 }
666 #ifdef STATISTICS
667 printf("guessed %i, unguessed %i, iterations %i\n", guessed, unguessed,
668 iters);
669 guessed2 += guessed;
670 unguessed2 += unguessed;
671 iters2 += iters;
672 #endif
673 }
674 #endif
675 #endif
676
677 #undef FORMULALOOP
678 #undef PCHECK
679 #undef SAVEZMAG
680 #ifndef SMOOTHMODE
681 #ifdef SMOOTH
682 #define SMOOTHMODE
683 #include "docalc.h"
684 #endif
685 #endif
686
687 /*cleanup for next formula */
688 #undef NSFORMULALOOP
689 #undef SFORMULALOOP
690 #undef PRESMOOTH
691 #undef SMOOTH
692 #undef SMOOTHMODE
693 #undef RANGE
694 #undef JULIA
695 #undef PERI
696 #undef SPERI
697 #undef INIT
698 #undef VARIABLES
699 #undef PRETEST
700 #undef BTEST
701 #undef FORMULA
702 #undef CALC
703 #undef SCALC
704 #undef RPIP
705 #undef POSTCALC
706 #undef UNCOMPRESS
707 #undef SAVE
708 #undef SAVEVARIABLES
709 #undef RESTORE
710 #undef UFORMULA
711 #undef UEND
712