From cfbb058289835140d9ddb3c3ae91bee8fdf67037 Mon Sep 17 00:00:00 2001 From: Randy McShandy Date: Thu, 7 Dec 2023 19:48:36 -0600 Subject: [PATCH 1/1] Hey, it kinda works! Single-threaded --- main.c | 191 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 191 insertions(+) create mode 100755 main.c diff --git a/main.c b/main.c new file mode 100755 index 0000000..58c65fb --- /dev/null +++ b/main.c @@ -0,0 +1,191 @@ +#include +#include +#include +#define GRID_X 256 +#define GRID_Y 256 + +typedef struct +{ + float v[3][3]; + +} Mat3; + +typedef struct +{ + float a; + float b; +} FVec2; + +typedef struct +{ + float diff_a; + float diff_b; + float feed; + float kill; + float delta_t; +} RD_Opts; + +FVec2 grid[GRID_X][GRID_Y]; +FVec2 grid_prime[GRID_X][GRID_Y]; + +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 + .v = { + {0.05f, 0.2f, 0.05f}, + {0.20, -1.0f, 0.20}, + {0.05f, 0.2f, 0.05f} + } +}; + +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]); + } + } + + return result; +} + +float rd_a_prime(RD_Opts opts, int x, int y, Mat3 kernel, float A, float B) +{ + 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; + if(isnan(a_prime)) + { + int bp = 0; + printf("nan!\n"); + } + return a_prime; +} + +float rd_b_prime(RD_Opts opts, int x, int y, Mat3 kernel, float A, float B) +{ + 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; + + if(isnan(b_prime)) + { + int bp = 0; + printf("nan!\n"); + } + return b_prime; +} + +int main(int argc, char** argv) +{ + char chars[6] = {' ', ' ', ' ', '+', '#', '@'}; + const int max_chars = sizeof(chars) / sizeof(chars[0]) - 1; + // Connectivity increases as time passes (iterations or time delta)? + // Seems we like <1.0 time deltas and higher iterations + + RD_Opts puffer_opts = + { + .diff_a = 1.0f, + .diff_b = 0.5f, + .feed = 0.055f, + .kill = 0.062f, + .delta_t = 0.6f, + }; + + RD_Opts coral_opts = + { + .diff_a = 0.7f, + .diff_b = 0.3f, + .feed = 0.0545f, + .kill = 0.062f, + .delta_t = 1.0f, + }; + + RD_Opts opts = coral_opts; + + + 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; + } + } + + grid[GRID_X/2][GRID_Y/2].b = 1.0f; + grid[GRID_X/4][GRID_Y/2].b = 1.0f; + grid[(GRID_X/6) * 5][4 * (GRID_Y/9)].b = 1.0f; + grid[GRID_X/4 * 3][GRID_Y/3].b = 1.0f; + grid[GRID_X-1][GRID_Y/5].b = 1.0f; + + const int max_iterations = 1 * 1e4; + for (int iterations = 0; iterations < max_iterations; iterations++) + { + for (int x = 1; x < GRID_X-1; x++) + { + for (int y = 1; y < GRID_Y-1; y++) + { + FVec2 each = grid[x][y]; + grid_prime[x][y].a = rd_a_prime(opts, x, y, laplacian_kernel, each.a, each.b); + grid_prime[x][y].b = rd_b_prime(opts, x, y, laplacian_kernel, each.a, each.b); + } + } + if (iterations % (int)1e3 == 0) + { + for (int x = 1; x < GRID_X-1; x++) + { + for (int y = 1; y < GRID_Y-1; y++) + { + int idx = floor(grid[x][y].a * (max_chars)); + printf("%c", chars[idx]); + } + printf("\n"); + } + } + memcpy(grid, grid_prime, sizeof(FVec2) * GRID_X * GRID_Y); + } + + printf("Opts: {\nIterations: %d\nTime delta: %f\nDiffA: %f\nDiffB: %f\nFeed: %f\nKill: %f\n\n}\n", max_iterations, opts.delta_t, opts.diff_a, opts.diff_b, opts.feed, opts.kill); + + for (int x = 1; x < GRID_X-1; x++) + { + for (int y = 1; y < GRID_Y-1; y++) + { + int idx = floor(grid[x][y].a * (max_chars)); + printf("%c", chars[idx]); + } + printf("\n"); + } + + return 0; +} + -- 2.49.0