Commits

Ben Bass  committed 75f2507

add basic optimisations for mset cardioid and period 2 disc

  • Participants
  • Parent commits 7c01543

Comments (0)

Files changed (2)

File mandelmap.cc

-// MandelMap - copyright Ben Bass 2010-2011
+// MandelMap - Ben Bass 2010-2012
 
 #include <cstdio>
 #include <math.h>  // for NAN. Not in the C++ library...
 const FLOAT MAX_Y = 1.5;
 const FLOAT BAILOUT = 4.0;
 
+// mpoint(real, image, ...) -> 1 if mset else 0
 inline long mpoint(FLOAT r,
                   FLOAT i,
                   pinfo* p,
     FLOAT y = p->y;
     int k = MAX_ITER;
 
+#ifdef CARDIOID_OPT
+    // based on 'Optimizations' section in wikipedia
+    // http://en.wikipedia.org/wiki/Mandelbrot_set#Optimizations
+    double roff = (r-0.25);
+    double isq = i*i;
+    double dist = roff*roff + isq;
+    if (dist * (dist + roff) < 0.25*isq)
+    {
+        return 1;
+    }
+
+    roff = r+1.0;
+    if (roff*roff + isq < 1.0/16.0)
+    {
+        return 1;
+    }
+#endif
     // this is the core 'inner loop'.
     // It should be fast, which is why it doesn't
     // have lots of whitespace in it ;)
         x2 = x*x;
         y2 = y*y;
         if ((x2+y2) >= BAILOUT) break;
-        y = 2.0*x*y+r;
-        x = x2-y2+i;
+        y = 2.0*x*y+i;
+        x = x2-y2+r;
     } while (--k);
 
     if (!k)
     }
 }
 
-int INIT_MAX_ITER;
 pinfo* FPTR_START;
 FLOAT STEP_SIZE;
 pthread_mutex_t acc_lock;
         for (FLOAT x=MIN_X; x<MAX_X; x+=STEP_SIZE)
         {
             pinfo* p=local_fptr++;
-            local_inside += mpoint(y, x, p, MAX_ITER);
+            local_inside += mpoint(x, y, p, MAX_ITER);
         }
         // skip ahead the relevant number of lines
         local_fptr += (row_skip * (NUM_THREADS-1));
     void* sourceBuffer = fptr;  // for unmmaping later
 
     printf("mapped file %d at %p %llu bytes\n", fd, fptr, mapsize);
+    printf("Current MAX_ITER: %d\n", fheader->max_iter);
 
     pthread_t t_handle[NUM_THREADS];
     int t_num[NUM_THREADS];
 
     // setup global read-only vars //BADBADBAD (but could be worse)
-    INIT_MAX_ITER = fheader->max_iter;
     FPTR_START = fptr;
     STEP_SIZE = step;
 
 #ifndef MANDEL_MAP_H
 #define MANDEL_MAP_H
 
+#define CARDIOID_OPT
+
 typedef double FLOAT;
 
 const int BINARY_DIGITS = 11;