updated documentation/removed TODOs

5 files changed, 122 insertions, 57 deletions

diff --git a/makefile b/makefile
index 103b693..6aa3232 100644
--- a/makefile
+++ b/makefile
@@ -67,7 +67,7 @@ presentation/presentation.html: presentation/presentation.md
 
 analysis: analysis/analysis.html
 
-analysis/analysis.html: analysis/analysis.Rmd # TODO: add compile command
+analysis/analysis.html: analysis/analysis.Rmd analysis/data.csv
 	Rscript -e "setwd('analysis'); rmarkdown::render('analysis.Rmd')"
 
 clean:
diff --git a/src/cuda-fractals.cu b/src/cuda-fractals.cu
index ef94a3f..52e7428 100644
--- a/src/cuda-fractals.cu
+++ b/src/cuda-fractals.cu
@@ -1,10 +1,10 @@
+/*
+ * Compute fractals using CUDA
+ */
 #include <cuda_runtime.h>
 #include <thrust/complex.h>
 #include "fractals.h"
-
-//These won't work/ need to be modifiied for the cuda version
-//#include "grids.h"
-//#include "fractals.h"
+#include "grids.h"
 
 /*
  * Macro for checking CUDA errors
@@ -22,8 +22,13 @@
 #define CBASE double
 #endif
 
+#ifndef BLOCK_SIZE_X
 #define BLOCK_SIZE_X 16
+#endif
+
+#ifndef BLOCK_SIZE_Y
 #define BLOCK_SIZE_Y 16
+#endif
 
 // this ain't pretty, but this logic should never change for kernels
 #define SET_ROW_COL \
@@ -32,6 +37,9 @@
     if(row >= rows || col >= cols) return
 
 
+/*
+ * Device helper function to cconvert a grid point into its corresponding complex number
+ */
 __device__
 thrust::complex<CBASE> grid_to_complex(const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t row, const size_t col, const size_t rows, const size_t cols){
     const CBASE x_min = lower_left.real();
@@ -48,6 +56,9 @@ thrust::complex<CBASE> grid_to_complex(const thrust::complex<CBASE> lower_left,
     return z;
 }
 
+/*
+ * Device function to compute if a point is or is not in the mandelbrot set
+ */
 __device__
 byte mandelbrot(const thrust::complex<CBASE> z0, const byte max_iterations){
     thrust::complex<CBASE> z = z0;
@@ -59,6 +70,9 @@ byte mandelbrot(const thrust::complex<CBASE> z0, const byte max_iterations){
     return iteration;
 }
 
+/*
+ * Kernel to run mandelbrot on an entire grid
+ */
 __global__
 void mandelbrot_kernel(byte* grid_data, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -68,6 +82,9 @@ void mandelbrot_kernel(byte* grid_data, const byte max_iterations, const thrust:
     grid_data[row*cols + col] = mandelbrot(z, max_iterations);
 }
 
+/*
+ * Device function to compute if a point is or is not in the tricorn set
+ */
 __device__ 
 byte tricorn(const thrust::complex<CBASE> z0, const byte max_iterations){
     thrust::complex<CBASE> z = z0;
@@ -79,6 +96,9 @@ byte tricorn(const thrust::complex<CBASE> z0, const byte max_iterations){
     return iteration;
 }
 
+/*
+ * Kernel to run tricorn on an entire grid
+ */
 __global__
 void tricorn_kernel(byte* grid_data, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -88,6 +108,9 @@ void tricorn_kernel(byte* grid_data, const byte max_iterations, const thrust::co
     grid_data[row*cols + col] = tricorn(z, max_iterations);
 }
 
+/*
+ * Device function to compute if a point is or is not in the burning ship set
+ */
 __device__
 byte burning_ship(const thrust::complex<CBASE> z0, const byte max_iterations){
     thrust::complex<CBASE> z = z0;
@@ -101,6 +124,9 @@ byte burning_ship(const thrust::complex<CBASE> z0, const byte max_iterations){
     return iteration;
 }
 
+/*
+ * Kernel to run burning_ship on an entire grid
+ */
 __global__
 void burning_ship_kernel(byte* grid_data, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -110,6 +136,9 @@ void burning_ship_kernel(byte* grid_data, const byte max_iterations, const thrus
     grid_data[row*cols + col] = burning_ship(z, max_iterations);
 }
 
+/*
+ * Device function to compute if a point is or is not in the multibrot set
+ */
 __device__
 byte multibrot(const thrust::complex<CBASE> z0, const byte max_iterations, const double d){
     thrust::complex<CBASE> z = z0;
@@ -121,6 +150,9 @@ byte multibrot(const thrust::complex<CBASE> z0, const byte max_iterations, const
     return iteration;
 }
 
+/*
+ * Kernel to run multibrot on an entire grid
+ */
 __global__
 void multibrot_kernel(byte* grid_data, const double degree, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -130,6 +162,9 @@ void multibrot_kernel(byte* grid_data, const double degree, const byte max_itera
     grid_data[row*cols + col] = multibrot(z, max_iterations, degree);
 }
 
+/*
+ * Device function to compute if a point is or is not in the multicorn set
+ */
 __device__
 byte multicorn(const thrust::complex<CBASE> z0, const byte max_iterations, const double d){
     thrust::complex<CBASE> z = z0;
@@ -141,6 +176,9 @@ byte multicorn(const thrust::complex<CBASE> z0, const byte max_iterations, const
     return iteration;
 }
 
+/*
+ * Kernel to run multicorn on an entire grid
+ */
 __global__
 void multicorn_kernel(byte* grid_data, const double degree, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -150,6 +188,9 @@ void multicorn_kernel(byte* grid_data, const double degree, const byte max_itera
     grid_data[row*cols + col] = multicorn(z, max_iterations, degree);
 }
 
+/*
+ * Device function to compute if a point is or is not in the julia set
+ */
 __device__
 byte julia(const thrust::complex<CBASE> z0, const byte max_iterations, const thrust::complex<CBASE> c, const double R){
     thrust::complex<CBASE> z = z0;
@@ -161,6 +202,9 @@ byte julia(const thrust::complex<CBASE> z0, const byte max_iterations, const thr
     return iteration;
 }
 
+/*
+ * Kernel to run julia on an entire grid
+ */
 __global__
 void julia_kernel(byte* grid_data, const thrust::complex<CBASE> constant, const double radius, const byte max_iterations, const thrust::complex<CBASE> lower_left, const thrust::complex<CBASE> upper_right, const size_t rows, const size_t cols){
     SET_ROW_COL;
@@ -170,6 +214,7 @@ void julia_kernel(byte* grid_data, const thrust::complex<CBASE> constant, const
     grid_data[row*cols + col] = julia(z, max_iterations, constant, radius);
 }
 
+// prevent c++ name mangling so the grid functions can be linked against in standard c code
 extern "C" {
 void mandelbrot_grid(grid_t* grid, const grid_gen_params* params){
     const size_t size = grid->size;
@@ -181,9 +226,7 @@ void mandelbrot_grid(grid_t* grid, const grid_gen_params* params){
 
     byte* d_grid_data;
     CHECK(cudaMalloc(&d_grid_data, size*sizeof(byte)));
-    //TODO: find good sizes
     dim3 block_size(BLOCK_SIZE_X, BLOCK_SIZE_Y);
-    //dim3 grid_size(0,0);
     dim3 grid_size((cols + block_size.x - 1) / block_size.x, (rows + block_size.y - 1) / block_size.y);
     mandelbrot_kernel<<<grid_size, block_size>>>(d_grid_data, max_iterations, lower_left, upper_right, rows, cols);
     CHECK(cudaDeviceSynchronize());
diff --git a/src/fractals.c b/src/fractals.c
index 898f30d..1fcce3d 100644
--- a/src/fractals.c
+++ b/src/fractals.c
@@ -16,10 +16,16 @@
 #define NUM_RUNS 5
 #endif
 
+/*
+ * Prints out usage information for the program
+ */
 void print_usage(FILE* file, const char* program_name){
     fprintf(file, "Usage: %s [-v] [-i iterations] [-x x_res] [-y y_res] [-z magnification] [-d degree] [-c constant] [-r radius] [-l lower_left] [-u upper_right] [-o output_grid] -f fractal\n", program_name);
 }
 
+/*
+ * Print out the full help message for the program
+ */
 void print_help(){
     printf("Options:\n"
             "  -i, --iterations <value>        the number of iterations (default: 25, max 255)\n"
@@ -41,6 +47,9 @@ void print_help(){
             "\nExits with a status code of 1 if the program encounters an error, exits with 2 if an argument is incorrect\n");
 }
 
+/*
+ * Prints out additional information about the program, such as the precision it was compiled with
+ */
 void print_info(const char* program_name){
 #ifdef EXTENDED_PRECISION
     printf("Compiled with long double float precision\n");
@@ -50,6 +59,9 @@ void print_info(const char* program_name){
 #endif
 }
 
+/*
+ * Runs a fractal generator NUM_RUNS times and returns the average of those runs
+ */
 double time_fractal(fractal_generator generator, grid_t* grid, grid_gen_params* params){
     struct timespec start, end;
     clock_gettime(CLOCK_MONOTONIC, &start);
@@ -67,6 +79,47 @@ static inline void parse_complex(const char* string, complex_t* z){
     };
 }
 
+/*
+ * Parses a fractal generator form a string, exitting if it is supplied a parameter which it doesn't support
+ */
+fractal_generator parse_fractal_generator(const char* argument, const bool param_is_degree, const bool param_is_cr){
+    if(strncmp(argument, "mandelbrot", strlen("mandelbrot")) == 0) {
+        return mandelbrot_grid;
+    }
+    else if(strncmp(argument, "tricorn", strlen("tricorn")) == 0) {
+        return tricorn_grid;
+    }
+    else if(strncmp(argument, "multibrot", strlen("multibrot")) == 0) {
+        if(param_is_cr){
+            fprintf(stderr, "multibrot requires a degree, not constant and radius, exitting\n");
+            exit(EXIT_BAD_ARGUMENT);
+        }
+        return multibrot_grid;
+    }
+    else if(strncmp(argument, "multicorn", strlen("multicorn")) == 0) {
+        if(param_is_cr){
+            fprintf(stderr, "multicorn requires a degree, not constant and radius, exitting\n");
+            exit(EXIT_BAD_ARGUMENT);
+        }
+        return multicorn_grid;
+    }
+    else if(strncmp(argument, "burning_ship", strlen("burning_ship")) == 0) {
+        return burning_ship_grid;
+    }
+    else if(strncmp(argument, "julia", strlen("julia")) == 0) {
+        if(param_is_degree){
+            fprintf(stderr, "julia requires a constant and a radius, not a degree, exitting\n");
+            exit(EXIT_BAD_ARGUMENT);
+        }
+        return julia_grid;
+    }
+    else {
+        fprintf(stderr, "Invalid fractal type: %s, see --help for a list of supported fractals\n", argument);
+        exit(EXIT_BAD_ARGUMENT);
+    }
+}
+
+
 int main(const int argc, char *argv[]) {
     struct winsize w;
     ioctl(STDOUT_FILENO, TIOCGWINSZ, &w);
@@ -82,6 +135,7 @@ int main(const int argc, char *argv[]) {
     bool performance = false;
 
     //degree is mutually exclusive with constant and radius
+    //this could be simplified if grid_gen_params was a struct instead of a union but ¯\_(ツ)_/¯
     bool param_is_degree = false;
     bool param_is_cr = false;
     CBASE degree = 2;
@@ -177,46 +231,8 @@ int main(const int argc, char *argv[]) {
                 param_is_cr = true;
                 break;
             case 'f':
-                if(strncmp(optarg, "mandelbrot", strlen("mandelbrot")) == 0) {
-                    fractal_name = "mandelbrot";
-                    generator = mandelbrot_grid;
-                }
-                else if(strncmp(optarg, "tricorn", strlen("tricorn")) == 0) {
-                    fractal_name = "tricorn";
-                    generator = tricorn_grid;
-                }
-                else if(strncmp(optarg, "multibrot", strlen("multibrot")) == 0) {
-                    if(param_is_cr){
-                        fprintf(stderr, "multibrot requires a degree, not constant and radius, exitting\n");
-                        exit(EXIT_BAD_ARGUMENT);
-                    }
-                    fractal_name = "multibrot";
-                    generator = multibrot_grid;
-                }
-                else if(strncmp(optarg, "multicorn", strlen("multicorn")) == 0) {
-                    if(param_is_cr){
-                        fprintf(stderr, "multicorn requires a degree, not constant and radius, exitting\n");
-                        exit(EXIT_BAD_ARGUMENT);
-                    }
-                    fractal_name = "multicorn";
-                    generator = multicorn_grid;
-                }
-                else if(strncmp(optarg, "burning_ship", strlen("burning_ship")) == 0) {
-                    fractal_name = "burning ship";
-                    generator = burning_ship_grid;
-                }
-                else if(strncmp(optarg, "julia", strlen("julia")) == 0) {
-                    if(param_is_degree){
-                        fprintf(stderr, "julia requires a constant and a radius, not a degree, exitting\n");
-                        exit(EXIT_BAD_ARGUMENT);
-                    }
-                    fractal_name = "julia";
-                    generator = julia_grid;
-                }
-                else {
-                    fprintf(stderr, "Invalid fractal type: %s, see --help for a list of supported fractals\n", optarg);
-                    exit(EXIT_BAD_ARGUMENT);
-                }
+                fractal_name = optarg;
+                generator = parse_fractal_generator(optarg, param_is_degree, param_is_cr);
                 break;
             case 'z':
                 if(sscanf(optarg, CFORMAT, &magnification) != 1){
diff --git a/src/grids.c b/src/grids.c
index e65ba79..279fdb5 100644
--- a/src/grids.c
+++ b/src/grids.c
@@ -1,3 +1,6 @@
+/*
+ * Functions for the creation, manipulation, and evaluation of grids.
+ */
 #include <complex.h>
 #include <stdio.h>
 #include <stdlib.h>
@@ -92,6 +95,8 @@ bool grid_equal(const grid_t* grid1_p, const grid_t* grid2_p){
 
 /*
  * Checks if two grids have a given maximum difference
+ * This function has no few if any uses, experimentall the grids are either exact or very different,
+ * indicating an error in their generation
  */
 bool grid_allclose(const grid_t* restrict grid1, const grid_t* restrict grid2, const byte max_error){
     if(grid1->x != grid2->x || grid1->y != grid2->y ||
@@ -101,8 +106,6 @@ bool grid_allclose(const grid_t* restrict grid1, const grid_t* restrict grid2, c
     }
     const size_t size = grid1->size;
     for(size_t i = 0; i < size; i++){
-        //FIXME: figure out how to handle difference between two unsigned values
-        //       possibly cast them to sized longs?
         if(abs(grid1->data[i] - grid2->data[i]) >= max_error) return false;
     }
     return true;
@@ -172,12 +175,11 @@ void zoom_grid(grid_t* restrict grid, const CBASE magnification){
  *
  * Returns 0 on success
  *
- * FIXME: add dynamic bounding box resolution by storing sizeof points in the file
  * The .grid format is a binary file format
- * The first 3 bytes of the file are a magic number defined in grids.h
+ * The first 3 bytes of the file are a magic number defined in grids.h (a goose !)
  * The next 16 bytes are the grid dimensions (x then y)
- * The next 8 bytes is the max_iterations
- * The next 8 bytes are the size of a grid point in bytes
+ * The next 1 byte is the max_iterations
+ * The next 8 bytes are the precision for the complex numbers
  * The next 2*precision bytes are the lower left and upper right corners
  * The rest of the file is the data for the grid, which should be exactly x*y*8 bytes
  */
@@ -239,14 +241,13 @@ void print_grid(FILE* file, const grid_t* grid){
     const byte iterations = grid->max_iterations;
     const byte* data = grid->data;
 
-    //TODO: set values in output buffer rather than multiple printf calls
-    //      the buffer needs to be larger to hold newlines
     char* output_buffer = malloc(size + grid->y-1);
     if(!output_buffer){
         fprintf(stderr, "Failed to allocate output buffer for %zu points\n", size);
         return;
     }
 
+    // sets an output buffer to minimize number of syscalls
     setvbuf(file, output_buffer, _IOFBF, size + grid->y - 1);
 
     const char point_types[] = { ' ', '.', '*', '%', '#'};
@@ -323,8 +324,10 @@ grid_t* read_grid(FILE* restrict file){
     if(fread(&lower_left, sizeof(complex_t), 1, file) != 1){ longjmp(file_read_error, 1); }
     if(fread(&upper_right, sizeof(complex_t), 1, file) != 1){ longjmp(file_read_error, 1); }
 
-    //TODO: look into mmaping the file to data, offseting by the bounding and resolution information
-    //      this would likely require an alloc_grid function, similar to jeff's implementation in hw03
+    // NOTE: for very large sizes it might make sense to memmap the file
+    //       this would work on the CPU but break the cuda implmentation
+    //       a potential fix would to just make multiple grids, aand stich them in
+    //       after image rendering
     grid_t* grid = create_grid(x, y, max_iterations, lower_left, upper_right);
     if(!grid){
         return NULL;
diff --git a/src/grids.h b/src/grids.h
index 6bf96f4..129cfd5 100644
--- a/src/grids.h
+++ b/src/grids.h
@@ -20,6 +20,9 @@ typedef struct {
     CBASE im;
 } complex_t;
 
+// I've gone back and forth on the size for a grid point, eventually I settled on a byte
+// ideally it should just be its own type which can be configured similarly to precision.h
+// but it does not seem worth
 typedef unsigned char byte;
 
 typedef struct {