star-sim

Barnes-Hut gravity simulation.
git clone git://git.amin.space/star-sim.git
Log | Files | Refs | README | LICENSE

commit 6941f174461494903945cd66df2485c5d914415d
parent 1d3ac75b8034c7ec9568ba65989cb56947ac519d
Author: amin <dev@aminmesbah.com>
Date:   Thu, 20 Apr 2017 20:09:55 +0000

Add Barnes-Hut partition grid overlay.

The Barnes-Hut algorithm reduces the number of comparisons between
particles in an n-body simulation by recursively partitioning the
simulation space into rectangular cells so that distant groups of
particles can be treated as one particle.

I added a basic quadtree to keep track of the cells. The cell borders
are displayed as green lines. MAX_STARS_PER_CELL can be adjusted to
change the threshold at which a cell is subdivided.

FossilOrigin-Name: 7b3641532726c183fdf400c0858776218526384d182fa25f060107d05b4c3dbb
Diffstat:
Abarnes_hut.h | 108+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mmain.cpp | 109++++++++++++++++++++++++++++++++++++++-----------------------------------------
2 files changed, 161 insertions(+), 56 deletions(-)

diff --git a/barnes_hut.h b/barnes_hut.h @@ -0,0 +1,108 @@ +/* Quadtree implementation to keep track of recursive subdivisions of the + * simulation field into rectangular cells. + */ + +#ifndef BARNES_HUT_H +#define BARNES_HUT_H + +// TODO: Why is there a segfault when this is set to 0? +// TODO: Why is there a segfault only very rarely when this is set to 1? +#define MAX_STARS_PER_CELL 1 + +struct Cell +{ + float center_x; + float center_y; + float distance_x; + float distance_y; +}; + +struct QuadTreeNode +{ + Cell *cell; + QuadTreeNode *children[4]; +}; + + +Cell* cell_init(float center_x, float center_y, float distance_x, float distance_y) +{ + Cell *cell = (Cell*)malloc(sizeof(Cell)); + if (cell) + { + cell->center_x = center_x; + cell->center_y = center_y; + cell->distance_x = distance_x; + cell->distance_y = distance_y; + } + return cell; +} + +int cell_count_stars_contained(Cell *c, Star stars[], int num_stars) +{ + int count = 0; + for (int i = 0; i < num_stars; ++i) + { + if (stars[i].x <= c->center_x + c->distance_x + && stars[i].x >= c->center_x - c->distance_x + && stars[i].y <= c->center_y + c->distance_y + && stars[i].y >= c->center_y - c->distance_y) + { + count++; + } + } + return count; +} + +QuadTreeNode* quad_tree_node_init(float center_x, float center_y, float distance_x, float distance_y) +{ + QuadTreeNode *node = (QuadTreeNode*)malloc(sizeof(QuadTreeNode)); + if (node) + { + node->cell = cell_init(center_x, center_y, distance_x, distance_y); + for (int i = 0; i < 4; ++i) + { + node->children[i] = nullptr; + } + } + return node; +} + +void quad_tree_node_free(QuadTreeNode *node) +{ + if (node) + { + for (int i = 0; i < 4; ++i) + { + if (node->children[i]) + { + quad_tree_node_free(node->children[i]); + } + } + free(node->cell); + free(node); + } +} + +void quad_tree_node_subdivide(QuadTreeNode *node, Star stars[], int num_stars) +{ + if (cell_count_stars_contained(node->cell, stars, num_stars) > MAX_STARS_PER_CELL) + { + float child_distance_x = node->cell->distance_x / 2; + float child_distance_y = node->cell->distance_y / 2; + + node->children[0] = quad_tree_node_init(node->cell->center_x + child_distance_x, node->cell->center_y + child_distance_y, child_distance_x, child_distance_y); + node->children[1] = quad_tree_node_init(node->cell->center_x - child_distance_x, node->cell->center_y + child_distance_y, child_distance_x, child_distance_y); + node->children[2] = quad_tree_node_init(node->cell->center_x - child_distance_x, node->cell->center_y - child_distance_y, child_distance_x, child_distance_y); + node->children[3] = quad_tree_node_init(node->cell->center_x + child_distance_x, node->cell->center_y - child_distance_y, child_distance_x, child_distance_y); + + for (int i = 0; i < 4; ++i) + { + if (node->children[i]) + { + quad_tree_node_subdivide(node->children[i], stars, num_stars); + } + } + } +} + +#endif diff --git a/main.cpp b/main.cpp @@ -6,21 +6,21 @@ #include <string.h> #include <time.h> #include "star.h" +#include "barnes_hut.h" #ifndef MAP_ANONYMOUS #define MAP_ANONYMOUS MAP_ANON #endif //#define TRAILS -//#define TOROID #define FULLSCREEN #define TITLE "Stars" #define SCREEN_WIDTH 640 #define SCREEN_HEIGHT 480 #define BYTES_PER_PIXEL 4 -#define NUM_STARS 200 -#define DRAG 0.994 +#define NUM_STARS 100 +#define DRAG 1 #define SECOND 1000.0f #define FPS 60 @@ -178,6 +178,7 @@ void set_pixel(SDLOffscreenBuffer buffer, uint32_t x, uint32_t y, uint32_t color } } + // TODO: change buffer to a pointer void clear_screen(SDLOffscreenBuffer buffer, uint32_t pixel_value) { @@ -187,7 +188,7 @@ void clear_screen(SDLOffscreenBuffer buffer, uint32_t pixel_value) } -void update(Star stars[], int num_stars, SDLWindowDimension* dimension) +void update(Star stars[], int num_stars, QuadTreeNode *cells_root, SDLWindowDimension* dimension) { for (int i = 0; i < num_stars; ++i) { @@ -200,30 +201,50 @@ void update(Star stars[], int num_stars, SDLWindowDimension* dimension) //printf("%d: (%f, %f)\n", i, stars[i].x, stars[i].y); stars[i].speed *= DRAG; -#ifdef TOROID - if (stars[i].x > dimension->width) - { - stars[i].x = 0; - } - else if (stars[i].x < 0) - { - stars[i].x = dimension->width; - } + //printf("%f, %f\n", stars[i].x, stars[i].y); + } - if (stars[i].y > dimension->height) - { - stars[i].y = 0; - } - else if (stars[i].y < 0) + int center_x = global_back_buffer.width / 2; + int center_y = global_back_buffer.height / 2; + int dist_x = global_back_buffer.width / 2; + int dist_y = global_back_buffer.height / 2; + + quad_tree_node_free(cells_root); + cells_root = quad_tree_node_init(center_x, center_y, dist_x, dist_y); + quad_tree_node_subdivide(cells_root, stars, num_stars); +} + +// TODO: pass a pointer to buffer? +void draw_grid(SDLOffscreenBuffer buffer, QuadTreeNode *node, uint32_t color) +{ + if (node) + { + for (int i = 0; i < 4; ++i) { - stars[i].y = dimension->height; + if (node->children[i]) + { + for ( + int x = (node->cell->center_x - node->cell->distance_x); + x <= (node->cell->center_x + node->cell->distance_x); + ++x) + { + set_pixel(buffer, x, node->cell->center_y, color); + } + + for ( + int y = (node->cell->center_y - node->cell->distance_y); + y <= (node->cell->center_y + node->cell->distance_y); + ++y) + { + set_pixel(buffer, node->cell->center_x, y, color); + } + draw_grid(buffer, node->children[i], color); + } } -#endif - //printf("%f, %f\n", stars[i].x, stars[i].y); } } -void render(SDLOffscreenBuffer buffer, float dt, Star stars[], int num_stars) +void render(SDLOffscreenBuffer buffer, float dt, Star stars[], int num_stars, QuadTreeNode *cells_root) { //printf("%f\n", dt); for (int i = 0; i < num_stars; ++i) @@ -234,6 +255,8 @@ void render(SDLOffscreenBuffer buffer, float dt, Star stars[], int num_stars) stars[i].y - (cosf(stars[i].angle) * stars[i].speed) * dt, stars[i].color); } + + draw_grid(buffer, cells_root, COLOR_GREEN); } @@ -280,41 +303,15 @@ int main(int argc, char **argv) stars[i].mass = 5; stars[i].size = star_calc_size(stars[i].mass); stars[i].color = COLOR_WHITE; - - //enum color_t color = (color_t)(rand() % NUM_COLORS); - //switch(color) - //{ - // case YELLOW: - // stars[i].color = COLOR_YELLOW; - // break; - // case ORANGE: - // stars[i].color = COLOR_ORANGE; - // break; - // case RED: - // stars[i].color = COLOR_RED; - // break; - // case MAGENTA: - // stars[i].color = COLOR_MAGENTA; - // break; - // case VIOLET: - // stars[i].color = COLOR_VIOLET; - // break; - // case BLUE: - // stars[i].color = COLOR_BLUE; - // break; - // case CYAN: - // stars[i].color = COLOR_CYAN; - // break; - // case GREEN: - // stars[i].color = COLOR_GREEN; - // break; - // default: - // stars[i].color = COLOR_WHITE; - // break; - //} //printf("%f, %f, %f\n", stars[i].angle, stars[i].speed, stars[i].mass); } + QuadTreeNode *cells_root = quad_tree_node_init( + global_back_buffer.width / 2, + global_back_buffer.height / 2, + global_back_buffer.width / 2, + global_back_buffer.height / 2); + while (running) { uint64_t current_ms = (SDL_GetPerformanceCounter() * SECOND) / SDL_GetPerformanceFrequency(); @@ -336,14 +333,14 @@ int main(int argc, char **argv) while (lag >= MS_PER_UPDATE) { - update(stars, NUM_STARS, &dimension); + update(stars, NUM_STARS, cells_root, &dimension); //printf("\t%" PRIu64 ", %f\n", lag, MS_PER_UPDATE); lag -= MS_PER_UPDATE; } #ifndef TRAILS clear_screen(global_back_buffer, COLOR_BACKGROUND); #endif - render(global_back_buffer, lag/SECOND, stars, NUM_STARS); + render(global_back_buffer, lag/SECOND, stars, NUM_STARS, cells_root); sdl_update_window(window, renderer, global_back_buffer); if (elapsed_ms <= MS_PER_FRAME) {