star-sim

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

barnes_hut.c (5825B)


      1 #include "barnes_hut.h"
      2 
      3 #include "star.h"
      4 
      5 #include "math.h"
      6 #include "stdlib.h"
      7 
      8 
      9 struct Cell *cell_init(float center_x, float center_y, float distance_x, float distance_y)
     10 {
     11     struct Cell *cell = (struct Cell *)malloc(sizeof(struct Cell));
     12     if (cell)
     13     {
     14         cell->center_x = center_x;
     15         cell->center_y = center_y;
     16         cell->distance_x = distance_x;
     17         cell->distance_y = distance_y;
     18         cell->mass_center_x = center_x;
     19         cell->mass_center_y = center_x;
     20         cell->mass_total = 0;
     21         cell->stars_sum_x = 0;
     22         cell->stars_sum_y = 0;
     23         cell->star = NULL;
     24     }
     25     return cell;
     26 }
     27 
     28 
     29 void cell_free(struct Cell *c)
     30 {
     31     if (c)
     32     {
     33         free(c);
     34     }
     35 }
     36 
     37 
     38 bool cell_contains_star(struct Cell *c, struct Star *s)
     39 {
     40     return (s->x <= c->center_x + c->distance_x
     41         && s->x >= c->center_x - c->distance_x
     42         && s->y <= c->center_y + c->distance_y
     43         && s->y >= c->center_y - c->distance_y);
     44 }
     45 
     46 
     47 bool cell_is_empty(struct Cell *c)
     48 {
     49     return !(c->star);
     50 }
     51 
     52 
     53 struct QuadTreeNode *quad_tree_node_init(float center_x, float center_y, float distance_x, float distance_y)
     54 {
     55     struct QuadTreeNode *node = (struct QuadTreeNode *)malloc(sizeof (struct QuadTreeNode));
     56     if (node)
     57     {
     58         node->cell = cell_init(center_x, center_y, distance_x, distance_y);
     59         node->ne = NULL;
     60         node->nw = NULL;
     61         node->sw = NULL;
     62         node->se = NULL;
     63     }
     64     return node;
     65 }
     66 
     67 
     68 void quad_tree_node_free(struct QuadTreeNode *node)
     69 {
     70     if (node)
     71     {
     72         quad_tree_node_free(node->ne);
     73         quad_tree_node_free(node->nw);
     74         quad_tree_node_free(node->sw);
     75         quad_tree_node_free(node->se);
     76         free(node->cell);
     77         free(node);
     78     }
     79 }
     80 
     81 
     82 bool quad_tree_node_is_leaf(struct QuadTreeNode *node)
     83 {
     84     if (!node)
     85     {
     86         return false;
     87     }
     88     return !(node->ne && node->nw && node->sw && node->se);
     89 }
     90 
     91 
     92 void quad_tree_node_subdivide(struct QuadTreeNode *node)
     93 {
     94     if (node && quad_tree_node_is_leaf(node))
     95     {
     96         float child_distance_x = node->cell->distance_x / 2;
     97         float child_distance_y = node->cell->distance_y / 2;
     98 
     99         node->ne = quad_tree_node_init(node->cell->center_x + child_distance_x, node->cell->center_y + child_distance_y, child_distance_x, child_distance_y);
    100         node->nw = quad_tree_node_init(node->cell->center_x - child_distance_x, node->cell->center_y + child_distance_y, child_distance_x, child_distance_y);
    101         node->sw = quad_tree_node_init(node->cell->center_x - child_distance_x, node->cell->center_y - child_distance_y, child_distance_x, child_distance_y);
    102         node->se = quad_tree_node_init(node->cell->center_x + child_distance_x, node->cell->center_y - child_distance_y, child_distance_x, child_distance_y);
    103     }
    104 }
    105 
    106 
    107 void quad_tree_node_insert_star(struct QuadTreeNode *node, struct Star *star)
    108 {
    109     if (node && cell_contains_star(node->cell, star))
    110     {
    111         if (quad_tree_node_is_leaf(node))
    112         {
    113             if (cell_is_empty(node->cell))
    114             {
    115                 node->cell->star = star;
    116             }
    117             else
    118             {
    119                 quad_tree_node_subdivide(node);
    120 
    121                 quad_tree_node_insert_star(node->ne, node->cell->star);
    122                 quad_tree_node_insert_star(node->nw, node->cell->star);
    123                 quad_tree_node_insert_star(node->sw, node->cell->star);
    124                 quad_tree_node_insert_star(node->se, node->cell->star);
    125 
    126                 quad_tree_node_insert_star(node->ne, star);
    127                 quad_tree_node_insert_star(node->nw, star);
    128                 quad_tree_node_insert_star(node->sw, star);
    129                 quad_tree_node_insert_star(node->se, star);
    130 
    131                 node->cell->star = NULL;
    132             }
    133         }
    134         else
    135         {
    136             node->cell->mass_total += star->mass;
    137             node->cell->stars_sum_x += star->x * star->mass;
    138             node->cell->stars_sum_y += star->y * star->mass;
    139 
    140             float x_avg = node->cell->stars_sum_x / node->cell->mass_total;
    141             float y_avg = node->cell->stars_sum_y / node->cell->mass_total;
    142 
    143             node->cell->mass_center_x = x_avg;
    144             node->cell->mass_center_y = y_avg;
    145 
    146             quad_tree_node_insert_star(node->ne, star);
    147             quad_tree_node_insert_star(node->nw, star);
    148             quad_tree_node_insert_star(node->sw, star);
    149             quad_tree_node_insert_star(node->se, star);
    150         }
    151     }
    152 }
    153 
    154 
    155 struct QuadTree *quad_tree_init(void)
    156 {
    157     struct QuadTree *qt = (struct QuadTree *)malloc(sizeof (struct QuadTree));
    158     if (qt)
    159     {
    160         qt->root = NULL;
    161     }
    162     return qt;
    163 }
    164 
    165 
    166 void quad_tree_free(struct QuadTree *qt)
    167 {
    168     if (qt)
    169     {
    170         quad_tree_node_free(qt->root);
    171         free(qt);
    172     }
    173 }
    174 
    175 
    176 bool node_is_sufficiently_far(struct QuadTreeNode *node, struct Star *star)
    177 {
    178     float s = node->cell->distance_x * 2;
    179     float dx = node->cell->mass_center_x - star->x;
    180     float dy = node->cell->mass_center_y - star->y;
    181     float d = hypotf(dx, dy);
    182 
    183     return s / d < THETA;
    184 }
    185 
    186 
    187 void quad_tree_calc_force_on_star(struct QuadTreeNode *node, struct Star *star)
    188 {
    189     if (!star || !node)
    190     {
    191         return;
    192     }
    193 
    194     if (quad_tree_node_is_leaf(node) && !cell_is_empty(node->cell) && !cell_contains_star(node->cell, star))
    195     {
    196         star_attract(node->cell->star, star);
    197     }
    198     else if (node_is_sufficiently_far(node, star))
    199     {
    200         star_attract_to_mass(star,
    201                 node->cell->mass_total,
    202                 node->cell->mass_center_x,
    203                 node->cell->mass_center_y);
    204     }
    205     else
    206     {
    207         quad_tree_calc_force_on_star(node->ne, star);
    208         quad_tree_calc_force_on_star(node->nw, star);
    209         quad_tree_calc_force_on_star(node->sw, star);
    210         quad_tree_calc_force_on_star(node->se, star);
    211     }
    212 }