-#include <stdint.h>
-#include <math.h>
#include <stdio.h>
-#include <stdlib.h>
-#include <string.h>
-#include <pthread.h>
-#include <unistd.h>
-#define STB_IMAGE_WRITE_IMPLEMENTATION
-#include "stb_image_write.h"
-
-#define GRID_X 512
-#define GRID_Y GRID_X
-
-#define SHADOW 0
-#if SHADOW
-#define printf(a,...)
-#endif
-
-typedef struct
-{
- double v[3][3];
-
-} Mat3;
-
-typedef struct
-{
- float a;
- float b;
-} FVec2;
-
-typedef struct
-{
- int max_iterations;
- double diff_a;
- double diff_b;
- double feed;
- double kill;
- double delta_t;
-
- char name[32];
-} RD_Opts;
-
-typedef struct
-{
- RD_Opts opts;
- FVec2** source;
- FVec2** dest;
- int start_x;
- int start_y;
- int width;
- int height;
-
- pthread_mutex_t* mutex;
- pthread_mutex_t worker_mutex;
- pthread_cond_t* worker_cond;
- pthread_cond_t* boss_cond;
- pthread_barrier_t* barrier;
-
- int* waiting_workers;
- int worker_count;
- int iterations;
- int max_iterations;
- int worker_id;
-} worker_arg;
-
-FVec2 ** grid;
-FVec2 ** grid_prime;
-FVec2 ** grid_temp;
-int master_iteration_done;
-
-// Connectivity increases as time passes (iterations or time delta)?
-// Seems we like <1.0 time deltas and higher iterations
-// diff values influence tightness, delta between them shouldn't get wider than ~0.5 or narrower than ~0.4
-RD_Opts puffer =
-{
- .max_iterations = 1e4,
- .diff_a = 0.7f,
- .diff_b = 0.3f,
- .feed = 0.055f,
- .kill = 0.062f,
- .delta_t = 1.0f,
-
- .name = "puffer"
-};
-
-RD_Opts worms =
-{
- .max_iterations = 1e4,
- .diff_a = 0.7f,
- .diff_b = 0.3f,
- .feed = 0.082,
- .kill = 0.06f,
- .delta_t = 0.6f,
-
- .name = "worms"
-};
-
-RD_Opts meiosis =
-{
- .max_iterations = 1e4,
- .diff_a = 0.7f,
- .diff_b = 0.3f,
- .feed = 0.0367f,
- .kill = 0.0649f,
- .delta_t = 1.0f,
-
- .name = "meiosis"
-};
-
-const Mat3 laplacian_kernel =
-{
- // oh shit
- // Sum of all members must balance to 0.0 (or outer members have to do 1.0)
- // it looks like any modifications should be on the order of 0.01 or less
- // upper value of 22 <> 23 seems to be a good point
- .v = {
- {0.05f, 0.20f, 0.05f},
- {0.20f, -1.0f, 0.20f},
- {0.05f, 0.20f, 0.05f}
- }
-};
-
-float kernel_sum(Mat3 kernel, Mat3 source)
-{
- float result = 0.0f;
-
- // Convolve kernel at grid[x,y]
- for (int ix = 0; ix < 3; ix++)
- {
- for (int iy = 0; iy < 3; iy++)
- {
- result += (kernel.v[ix][iy] * source.v[ix][iy]);
- }
- }
-
- return result;
-}
-
-float rd_a_prime(FVec2** source_grid, RD_Opts opts, int x, int y, Mat3 kernel, float A, float B)
-{
- float a_prime = 1.0f;
- int x_less = (x-1 < 0) ? x : x-1;
- int y_less = (y-1 < 0) ? y : y-1;
- int x_more = (x+1 == GRID_X) ? x : x+1;
- int y_more = (y+1 == GRID_Y) ? y : y+1;
- // Use species A in the convolution, b_prime will need B
- // For now we won't iterate over edge rows and columns to avoid special case vomit
- Mat3 source =
- {
- .v =
- {
- {grid[x_less][y_less].a, grid[x][y_less].a, grid[x_more][y_less].a},
- {grid[x_less][y+0].a, grid[x][y+0].a, grid[x_more][y+0].a},
- {grid[x_less][y_more].a, grid[x][y_more].a, grid[x_more][y_more].a},
- }
- };
-
- a_prime = A + (opts.diff_a*(kernel_sum(kernel, source)) - (A*B*B) + opts.feed*(1.0f-A)) * opts.delta_t;
- return a_prime;
-}
-
-float rd_b_prime(FVec2** source_grid, RD_Opts opts, int x, int y, Mat3 kernel, float A, float B)
-{
- float b_prime = 1.0f;
- int x_less = (x-1 < 0) ? x : x-1;
- int y_less = (y-1 < 0) ? y : y-1;
- int x_more = (x+1 == GRID_X) ? x : x+1;
- int y_more = (y+1 == GRID_Y) ? y : y+1;
- // Use species A in the convolution, b_prime will need B
- Mat3 source =
- {
- .v =
- {
- {grid[x_less][y_less].b, grid[x][y_less].b, grid[x_more][y_less].b},
- {grid[x_less][y+0].b, grid[x][y+0].b, grid[x_more][y+0].b},
- {grid[x_less][y_more].b, grid[x][y_more].b, grid[x_more][y_more].b},
- }
- };
-
- b_prime = B + (opts.diff_b*(kernel_sum(kernel, source)) + A*(B*B) - (opts.kill + opts.feed)*B) * opts.delta_t;
- return b_prime;
-}
-
-void* iterator(void* _arg)
-{
- worker_arg* warg = (worker_arg*)_arg;
- RD_Opts opts = warg->opts;
- int start_x = warg->start_x;
- int start_y = warg->start_y;
- int w = warg->width;
- int h = warg->height;
-
- for (warg->iterations = 0; warg->iterations < warg->max_iterations; warg->iterations++)
- {
- printf("worker %d: work unit %d/%d\n", warg->worker_id, warg->iterations, warg->max_iterations);
- for (int x = start_x; x < w + start_x && x < GRID_X; x++)
- {
- for (int y = start_y; y < h + start_y && y < GRID_Y; y++)
- {
- FVec2 each = grid[x][y];
- grid_prime[x][y].a = rd_a_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
- grid_prime[x][y].b = rd_b_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
- }
- }
-
- pthread_mutex_lock(warg->mutex);
- if (++*warg->waiting_workers == warg->worker_count)
- {
- grid_temp = grid;
- grid = grid_prime;
- grid_prime = grid_temp;
- *warg->waiting_workers = 0;
- printf("worker-boss %d: completing report\n", warg->worker_id);
- if ((warg->iterations % 100 == 0))
- {
- char buffer[GRID_X][GRID_Y];
- for (int x = 0; x < GRID_X; x++)
- {
- for (int y = 0; y < GRID_Y; y++)
- {
- buffer[x][y] = (uint8_t)(255.0f * grid[x][y].a);
- }
- }
- char name[64] = {0};
- sprintf(name, "img/%s/%d.png", opts.name, warg->iterations);
- stbi_write_png(name, GRID_X, GRID_Y, 1, buffer, GRID_X * sizeof(uint8_t));
- }
- }
- pthread_mutex_unlock(warg->mutex);
-
- pthread_barrier_wait(warg->barrier);
-
-#if 0
- // Last one done should wake up the boss thread
- int last_worker = 0;
- pthread_mutex_lock(warg->mutex);
- if (++*warg->waiting_workers == warg->worker_count)
- {
- last_worker = 1;
- }
- pthread_mutex_unlock(warg->mutex);
-
- if (last_worker == 0)
- {
- pthread_cond_wait(warg->worker_cond, &warg->worker_mutex);
- }
- else
- {
- grid_temp = grid;
- grid = grid_prime;
- grid_prime = grid_temp;
-
- // segfault somewhere in here lmao
- if ((warg->iterations % 100 == 0))
- {
- char buffer[GRID_X][GRID_Y];
- for (int x = 0; x < GRID_X; x++)
- {
- for (int y = 0; y < GRID_Y; y++)
- {
- buffer[x][y] = (uint8_t)(255.0f * grid[x][y].a);
- }
- }
- char name[64] = {0};
- sprintf(name, "img/%s/%d.png", opts.name, warg->iterations);
- stbi_write_png(name, GRID_X, GRID_Y, 1, buffer, GRID_X * sizeof(uint8_t));
- }
- printf("worker-boss %d: completing report\n", warg->worker_id);
-
- *warg->waiting_workers = 0;
- pthread_cond_broadcast(warg->worker_cond);
- }
-#endif
- }
-
- printf("worker %d: exiting\n", warg->worker_id);
- return _arg;
-}
-
-void* iterator_old(void* _arg)
-{
- worker_arg* warg = (worker_arg*)_arg;
- RD_Opts opts = warg->opts;
- int start_x = warg->start_x;
- int start_y = warg->start_y;
- int w = warg->width;
- int h = warg->height;
-
- for (int iterations = 0; iterations < warg->max_iterations; iterations++)
- {
- printf("worker %d: work unit %d/%d\n", warg->worker_id, iterations, warg->max_iterations);
- for (int x = start_x; x < w + start_x; x++)
- {
- for (int y = start_y; y < h + start_y; y++)
- {
- FVec2 each = grid[x][y];
- grid_prime[x][y].a = rd_a_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
- grid_prime[x][y].b = rd_b_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
- }
- }
-
- // Last one done should wake up the boss thread
- pthread_mutex_lock(warg->mutex);
- ++*warg->waiting_workers;
- if (*warg->waiting_workers == warg->worker_count)
- {
- printf("worker %d: waking up boss\n", warg->worker_id);
- pthread_cond_signal(warg->boss_cond);
- }
- pthread_mutex_unlock(warg->mutex);
-
- pthread_cond_wait(warg->worker_cond, &warg->worker_mutex);
-
- }
-
- printf("worker %d: exiting\n", warg->worker_id);
- return _arg;
-}
+#include "structs.h"
int main(int argc, char** argv)
{
- char chars[6] = {' ', ' ', ' ', '+', '#', '@'};
- const int max_chars = sizeof(chars) / sizeof(chars[0]) - 1;
-
- grid = (FVec2**)malloc(GRID_X * sizeof(FVec2*));
- for (int n = 0; n < GRID_X; n++)
- {
- grid[n] = (FVec2*)malloc(GRID_Y * sizeof(FVec2));
- }
-
- grid_prime = (FVec2**)malloc(GRID_X * sizeof(FVec2*));
- for (int n = 0; n < GRID_X; n++)
- {
- grid_prime[n] = (FVec2*)malloc(GRID_Y * sizeof(FVec2));
- }
-
- srand(time(NULL));
-
- RD_Opts opts = puffer;
-
- for (int x = 0; x < GRID_X; x++)
- {
- for (int y = 0; y < GRID_Y; y++)
- {
- grid[x][y].a = 1.0f;
- grid[x][y].b = 0.0f;
- }
- }
- const int max_rand = sqrt(GRID_X);
- for (int n = 0; n < max_rand; n++)
+ RD_Opts meiosis_2 =
{
- grid[rand() % GRID_X][rand() % GRID_Y].b = 1.0f;
- }
+ .max_iterations = 1e4,
+ .diff_a = 0.8f,
+ .diff_b = 0.4f,
+ .feed = 0.0367f,
+ .kill = 0.0649f,
+ .delta_t = 1.0f,
- int waiting_workers = 0;
- worker_arg warg = {
- opts, grid, grid_prime, 0, 0, (GRID_X), (GRID_Y), .worker_count = 8, .waiting_workers = &waiting_workers, .max_iterations = 1e4,
+ .name = "meiosis_2"
};
- pthread_t threads[warg.worker_count];
- pthread_mutex_t mutex;
- pthread_cond_t worker_cond;
- pthread_cond_t boss_cond;
- pthread_barrier_t barrier;
-
- pthread_mutex_init(&mutex, NULL);
- pthread_cond_init(&worker_cond, NULL);
- pthread_cond_init(&boss_cond, NULL);
- pthread_barrier_init(&barrier, NULL, warg.worker_count);
-
- warg.worker_cond = &worker_cond;
- warg.boss_cond = &boss_cond;
- warg.mutex = &mutex;
- warg.barrier = &barrier;
-
- worker_arg wargs[warg.worker_count];
-
- for (int t = 0; t < warg.worker_count; t++)
- {
- wargs[t] = warg;
- wargs[t].worker_id = t;
- wargs[t].width = (GRID_X/warg.worker_count) + ((t == warg.worker_count-1) ? 0 : 4);
- wargs[t].start_x = (wargs[t].width * t);
- pthread_mutex_init(&wargs[t].worker_mutex, NULL);
- printf("worker %d x_span %d, %d\n", t, wargs[t].start_x, wargs[t].width);
- pthread_create(&threads[t], NULL, iterator, &wargs[t]);
- }
-
- int should_quit = 0;
- while(!should_quit)
- {
- pthread_mutex_lock(warg.mutex);
- if (wargs[0].iterations == wargs[0].max_iterations-1)
- should_quit = 1;
- pthread_mutex_unlock(warg.mutex);
- }
-
- printf("boss: exiting loop\n");
- for (int t = 0; t < warg.worker_count; t++)
- {
- pthread_join(threads[t], NULL);
- }
-
- printf("Opts: {\nIterations: %d\nTime delta: %f\nDiffA: %f\nDiffB: %f\nFeed: %f\nKill: %f\n\n}\n", wargs[0].max_iterations, opts.delta_t, opts.diff_a, opts.diff_b, opts.feed, opts.kill);
- // It's end of execution, who cares
- free(grid);
- free(grid_prime);
+ generate_rd(4, puffer);
return 0;
}
+#include <stdio.h>
+#include <stdlib.h>
+#include <stdint.h>
+#include <math.h>
#include "structs.h"
+#define STB_IMAGE_WRITE_IMPLEMENTATION
+#include "stb_image_write.h"
-FVec2 grid[GRID_X][GRID_Y];
-FVec2 grid_prime[GRID_X][GRID_Y];
+#define SHADOW 0
+#if SHADOW
+#define printf(a,...)
+#endif
-float convolve(Mat3 kernel, Mat3 source)
-{
- float result = 0.0f;
-
- // Convolve kernel at grid[x,y]
- for (int ix = 0; ix < 3; ix++)
- {
- for (int iy = 0; iy < 3; iy++)
- {
- result += (kernel.v[ix][iy] * source.v[ix][iy]);
- }
+worker_arg warg;
+int waiting_workers = 0;
+int should_quit = 0;
+
+pthread_t* threads;
+pthread_mutex_t mutex;
+pthread_barrier_t barrier;
+
+FVec2 ** grid;
+FVec2 ** grid_prime;
+FVec2 ** grid_temp;
+
+float kernel_sum(Mat3 kernel, Mat3 source) {
+ float result = 0.0f;
+
+ // Convolve kernel at grid[x,y]
+ for (int ix = 0; ix < 3; ix++) {
+ for (int iy = 0; iy < 3; iy++) {
+ result += (kernel.v[ix][iy] * source.v[ix][iy]);
}
+ }
+
+ return result;
+}
+
+float rd_a_prime(FVec2 **source_grid, RD_Opts opts, int x, int y, Mat3 kernel,
+ float A, float B) {
+ float a_prime = 1.0f;
+ int x_less = (x - 1 < 0) ? x : x - 1;
+ int y_less = (y - 1 < 0) ? y : y - 1;
+ int x_more = (x + 1 == GRID_X) ? x : x + 1;
+ int y_more = (y + 1 == GRID_Y) ? y : y + 1;
+
+ // Use species A in the convolution, b_prime will need B
+ // For now we won't iterate over edge rows and columns to avoid special case
+ // vomit
+ Mat3 source = {
+ .v = {
+ {grid[x_less][y_less].a, grid[x][y_less].a, grid[x_more][y_less].a},
+ {grid[x_less][y + 0].a, grid[x][y + 0].a, grid[x_more][y + 0].a},
+ {grid[x_less][y_more].a, grid[x][y_more].a, grid[x_more][y_more].a},
+ }};
+
+ a_prime = A + (opts.diff_a * (kernel_sum(kernel, source)) - (A * B * B) +
+ opts.feed * (1.0f - A)) *
+ opts.delta_t;
+ return a_prime;
+}
+
+float rd_b_prime(FVec2** source_grid, RD_Opts opts, int x, int y, Mat3 kernel, float A, float B)
+{
+ float b_prime = 1.0f;
+ int x_less = (x-1 < 0) ? x : x-1;
+ int y_less = (y-1 < 0) ? y : y-1;
+ int x_more = (x+1 == GRID_X) ? x : x+1;
+ int y_more = (y+1 == GRID_Y) ? y : y+1;
+
+ // Use species A in the convolution, b_prime will need B
+ Mat3 source =
+ {
+ .v =
+ {
+ {grid[x_less][y_less].b, grid[x][y_less].b, grid[x_more][y_less].b},
+ {grid[x_less][y+0].b, grid[x][y+0].b, grid[x_more][y+0].b},
+ {grid[x_less][y_more].b, grid[x][y_more].b, grid[x_more][y_more].b},
+ }
+ };
+
+ b_prime = B + (opts.diff_b*(kernel_sum(kernel, source)) + A*(B*B) - (opts.kill + opts.feed)*B) * opts.delta_t;
+ return b_prime;
+}
+
+void* iterator(void* _arg)
+{
+ worker_arg* warg = (worker_arg*)_arg;
+ RD_Opts opts = warg->opts;
+ int start_x = warg->start_x;
+ int start_y = warg->start_y;
+ int w = warg->width;
+ int h = warg->height;
+
+ for (warg->iterations = 0; warg->iterations < warg->max_iterations; warg->iterations++)
+ {
+ printf("worker %d: work unit %d/%d\n", warg->worker_id, warg->iterations, warg->max_iterations);
+ for (int x = start_x; x < w + start_x && x < GRID_X; x++)
+ {
+ for (int y = start_y; y < h + start_y && y < GRID_Y; y++)
+ {
+ FVec2 each = grid[x][y];
+ grid_prime[x][y].a = rd_a_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
+ grid_prime[x][y].b = rd_b_prime(grid, opts, x, y, laplacian_kernel, each.a, each.b);
+ }
+ }
- return result;
+ pthread_mutex_lock(&mutex);
+ if (++waiting_workers == warg->worker_count)
+ {
+ grid_temp = grid;
+ grid = grid_prime;
+ grid_prime = grid_temp;
+ waiting_workers = 0;
+ printf("worker-boss %d: completing report\n", warg->worker_id);
+ if ((warg->iterations % 100 == 0))
+ {
+ char buffer[GRID_X][GRID_Y];
+ for (int x = 0; x < GRID_X; x++)
+ {
+ for (int y = 0; y < GRID_Y; y++)
+ {
+ buffer[x][y] = (uint8_t)(255.0f * grid[x][y].a);
+ }
+ }
+ char name[64] = {0};
+ sprintf(name, "img/%s/%d.png", opts.name, warg->iterations);
+ stbi_write_png(name, GRID_X, GRID_Y, 1, buffer, GRID_X * sizeof(uint8_t));
+ }
+ }
+ pthread_mutex_unlock(&mutex);
+ pthread_barrier_wait(&barrier);
+
+ }
+
+ // One last synchronization so boss thread doesn't die early
+ pthread_barrier_wait(&barrier);
+ should_quit = 1;
+
+ printf("worker %d: exiting\n", warg->worker_id);
+ return _arg;
}
-float rd_a_prime(RD_Opts opts, int x, int y, Mat3 kernel, float A, float B)
+int initialize(int worker_count, RD_Opts active_opt)
{
- float a_prime = 1.0f;
- // Use species A in the convolution, b_prime will need B
- // For now we won't iterate over edge rows and columns to avoid special case vomit
- Mat3 source =
- {
- .v =
- {
- {grid[x-1][y-1].a, grid[x][y-1].a, grid[x+1][y-1].a},
- {grid[x-1][y+0].a, grid[x][y+0].a, grid[x+1][y+0].a},
- {grid[x-1][y+1].a, grid[x][y+1].a, grid[x+1][y+1].a},
- }
- };
-
- a_prime = A + (opts.diff_a*(convolve(kernel, source)) - (A*B*B) + opts.feed*(1.0f-A)) * opts.delta_t;
-
- return a_prime;
+ grid = (FVec2**)malloc(GRID_X * sizeof(FVec2*));
+ for (int n = 0; n < GRID_X; n++)
+ {
+ grid[n] = (FVec2*)malloc(GRID_Y * sizeof(FVec2));
+ }
+
+ grid_prime = (FVec2**)malloc(GRID_X * sizeof(FVec2*));
+ for (int n = 0; n < GRID_X; n++)
+ {
+ grid_prime[n] = (FVec2*)malloc(GRID_Y * sizeof(FVec2));
+ }
+
+ srand(time(NULL));
+
+ for (int x = 0; x < GRID_X; x++)
+ {
+ for (int y = 0; y < GRID_Y; y++)
+ {
+ grid[x][y].a = 1.0f;
+ grid[x][y].b = 0.0f;
+ }
+ }
+
+ const int max_rand = sqrt(GRID_X);
+ for (int n = 0; n < max_rand; n++)
+ {
+ grid[rand() % GRID_X][rand() % GRID_Y].b = 1.0f;
+ }
+
+ warg = (worker_arg){
+ active_opt, grid, grid_prime,
+ 0, 0, (GRID_X), (GRID_Y),
+ .worker_count = worker_count,
+ .max_iterations = active_opt.max_iterations
+ };
+
+ threads = (pthread_t*)malloc(sizeof(pthread_t) * worker_count);
+ pthread_mutex_init(&mutex, NULL);
+ pthread_barrier_init(&barrier, NULL, warg.worker_count);
+
+ worker_arg wargs[worker_count];
+
+ for (int t = 0; t < warg.worker_count; t++)
+ {
+ wargs[t] = warg;
+ wargs[t].worker_id = t;
+ wargs[t].width = (GRID_X/worker_count) + ((t == worker_count-1) ? 0 : 4);
+ wargs[t].start_x = (wargs[t].width * t);
+ printf("worker %d x_span %d, %d\n", t, wargs[t].start_x, wargs[t].width);
+ pthread_create(&threads[t], NULL, iterator, &wargs[t]);
+ }
+ return 0;
}
-float rd_b_prime(RD_Opts opts, int x, int y, Mat3 kernel, float A, float B)
+int cleanup()
{
- float b_prime = 1.0f;
- // Use species A in the convolution, b_prime will need B
- Mat3 source =
- {
- .v =
- {
- {grid[x-1][y-1].b, grid[x][y-1].b, grid[x+1][y-1].b},
- {grid[x-1][y+0].b, grid[x][y+0].b, grid[x+1][y+0].b},
- {grid[x-1][y+1].b, grid[x][y+1].b, grid[x+1][y+1].b},
- }
- };
-
- b_prime = B + (opts.diff_b*(convolve(kernel, source)) + A*(B*B) - (opts.kill + opts.feed)*B) * opts.delta_t;
-
- return b_prime;
+ printf("boss: exiting loop\n");
+ for (int t = 0; t < warg.worker_count; t++)
+ {
+ pthread_join(threads[t], NULL);
+ }
+
+ free(grid);
+ free(grid_prime);
+ free(threads);
+
+ return 0;
+}
+
+int generate_rd(int worker_count, RD_Opts active_opt)
+{
+ initialize(worker_count, active_opt);
+
+ // Would be smart to set up a semaphore and move per-workload processing here
+ while(!should_quit) {}
+
+ printf("Opts: {\nIterations: %d\nTime delta: %f\nDiffA: %f\nDiffB: %f\nFeed: %f\nKill: %f\n\n}\n", active_opt.max_iterations, active_opt.delta_t, active_opt.diff_a, active_opt.diff_b, active_opt.feed, active_opt.kill);
+
+ cleanup();
+ return 0;
}
projects / barrow_crawler / commitdiff