From ac83b40fb9b0bdc880ebb8c840a63520500fefec Mon Sep 17 00:00:00 2001
From: josh <josh@45c1a28c-8058-47b2-ae61-ca45b979098e>
Date: Sun, 7 Dec 2008 21:44:20 +0000
Subject: [PATCH] finished NewtonComputation

git-svn-id: svn://anubis/gvsu@320 45c1a28c-8058-47b2-ae61-ca45b979098e
---
 cs677/final/NewtonComputation.cc | 61 +++++++++++++++++++++++++++++++-
 cs677/final/mpi-fractals.cc      | 31 +++++++++++++---
 2 files changed, 86 insertions(+), 6 deletions(-)

diff --git a/cs677/final/NewtonComputation.cc b/cs677/final/NewtonComputation.cc
index c9bcfc2..dbdc5ad 100644
--- a/cs677/final/NewtonComputation.cc
+++ b/cs677/final/NewtonComputation.cc
@@ -1,7 +1,66 @@
 
+#include <math.h>               /* sqrt() */
+#include <complex.h>            /* complex, I, creal(), cimag() */
 #include "NewtonComputation.h"
 
+#define LEEWAY 0.001
+#define ITERATIONS 32
+#define CLOSE_ENOUGH(z, x, y) \
+    (fabs(creal(z) - (x)) < LEEWAY && fabs(cimag(z) - (y)) < LEEWAY)
+
+/*
+ * This Computation class generates a design in the complex plane
+ *  which shows which complex points go to which of the six roots
+ *  of the equation z^6-1=0 using Newton's method of finding roots
+ */
 unsigned int NewtonComputation::compute(double x, double y)
 {
-    return 0x00FF8800;
+    double complex rootGuess = x + I * y;
+    static double halfRoot3 = sqrt(3) / 2.0;
+
+    for (int iter = 0; iter < ITERATIONS; iter++)
+    {
+        /* percentage of color to illuminate based on current iteration */
+        double pIlum = (double) (ITERATIONS - iter) / (double) ITERATIONS;
+        /*
+         * These if statements check to see if the complex number (the root guess)
+         *  is within LEEWAY distance of a real root. If so, a unique color is returned
+         *  reflecting which root it is close to and how many iterations it took
+         *  for the root guess to get to that root
+         */
+        if (CLOSE_ENOUGH(rootGuess, 1, 0))
+        {
+            return ((unsigned int)(0xFF * pIlum)) << 16;
+        }
+        if (CLOSE_ENOUGH(rootGuess, -1, 0))
+        {
+            return (((unsigned int)(0xFF * pIlum)) << 16) + (((unsigned int)(0xFF * pIlum)) << 8);
+        }
+        if (CLOSE_ENOUGH(rootGuess, 0.5, halfRoot3))
+        {
+            return ((unsigned int)(0xFF * pIlum)) << 8;
+        }
+        if (CLOSE_ENOUGH(rootGuess, -0.5, halfRoot3))
+        {
+            return (unsigned int)(0x88 * pIlum);
+        }
+        if (CLOSE_ENOUGH(rootGuess, 0.5, -halfRoot3))
+        {
+            return (unsigned int)(0xFF * pIlum);
+        }
+        if (CLOSE_ENOUGH(rootGuess, -0.5, -halfRoot3))
+        {
+            return (((unsigned int)(0xFF * pIlum)) << 16) + (unsigned int)(0xFF * pIlum);
+        }
+
+        /* This expression evaluates the next complex number to be tested
+         * for being close to a root. It uses Newton's method for finding
+         * roots of equations according to the following recursive equation:
+         *     x_n+1 = x_n - y_n / dy_n
+         * --> x_n+1 = x_n - x_n^6 / 6x_n^5
+         * More information can be found at:
+         * http://www.willamette.edu/~sekino/fractal/fractal.htm
+         */
+        rootGuess -= (cpow(rootGuess, 6.0) - 1.0) / (6.0 * cpow(rootGuess, 5.0));
+    }
 }
diff --git a/cs677/final/mpi-fractals.cc b/cs677/final/mpi-fractals.cc
index bd5b58b..540ee21 100644
--- a/cs677/final/mpi-fractals.cc
+++ b/cs677/final/mpi-fractals.cc
@@ -11,6 +11,8 @@
 #include "NewtonComputation.h"
 using namespace std;
 
+/* a "task" will be processing TASK_SIZE pixels */
+#define TASK_SIZE 100
 #define PROGNAME "Josh's CS677 Final : MPI Fractal Generator"
 #define getXVirt(x) (((x) - (width >> 1)) * zoom + x_center)
 #define getYVirt(y) ((-((y) - (height >> 1))) * zoom + y_center)
@@ -18,7 +20,7 @@ using namespace std;
 bool createWindow(int width, int height,
                   SDL_Surface ** screen, Uint32 ** pixels);
 void getSizes(int * rank, int * size, int * nprocs);
-void draw(int rank, int size, int nprocs, int width, int height,
+void draw(int rank, int world_size, int nprocs, int width, int height,
     Uint32 * pixels, Computation * computation);
 
 static double x_center = 0.0;
@@ -69,14 +71,19 @@ int main(int argc, char * argv[])
         SDL_Event event;
         bool going = true;
         bool window_success = createWindow(width, height, &screen, &pixels);
+        bool redraw = true;
 
         if (!window_success)
             going = false;
 
         while (going && SDL_WaitEvent(&event) != 0)
         {
-            draw(my_rank, world_size, nprocs, width, height,
-                pixels, computation);
+            if (redraw)
+            {
+                draw(my_rank, world_size, nprocs, width, height,
+                    pixels, computation);
+                redraw = false;
+            }
             SDL_UpdateRect(screen, 0, 0, 0, 0);
             switch (event.type)
             {
@@ -153,10 +160,10 @@ void getSizes(int * rank, int * size, int * nprocs)
     }
 }
 
-void draw(int rank, int size, int nprocs, int width, int height,
+void draw(int rank, int world_size, int nprocs, int width, int height,
     Uint32 * pixels, Computation * computation)
 {
-    if (size == 1)
+    if (world_size == 1)
     {
         for (int y = 0; y < height; y++)
         {
@@ -170,8 +177,22 @@ void draw(int rank, int size, int nprocs, int width, int height,
     }
     else if (rank == 0)
     {
+        int num_pixels = width * height;
+        for (int to_proc = 1; to_proc < world_size; to_proc++)
+        {
+        }
     }
     else
     {
+        for (;;)
+        {
+            /* wait to be told what to do */
+#if 0
+            MPI_Recv();
+            if ()   /* exit if we are done */
+                break;
+            MPI_Send();     /* report results back */
+#endif
+        }
     }
 }