Add a bounding volume hierarchy

aabb.c

@@ -0,0 +1,64 @@
+#include "aabb.h"
+
+#include <math.h>
+
+bool aabb_hit(const AABB *aabb, Ray r, double t_min, double t_max) {
+  /* X */
+  {
+    double inv_d = 1.0 / r.direction.x;
+    double t0 = (aabb->min.x - r.origin.x) * inv_d;
+    double t1 = (aabb->max.x - r.origin.x) * inv_d;
+    if (inv_d < 0.0) {
+      double tmp = t0;
+      t0 = t1;
+      t1 = tmp;
+    }
+    t_min = t0 > t_min ? t0 : t_min;
+    t_max = t1 < t_max ? t1 : t_max;
+    if (t_max <= t_min)
+      return false;
+  }
+  /* Y */
+  {
+    double inv_d = 1.0 / r.direction.y;
+    double t0 = (aabb->min.y - r.origin.y) * inv_d;
+    double t1 = (aabb->max.y - r.origin.y) * inv_d;
+    if (inv_d < 0.0) {
+      double tmp = t0;
+      t0 = t1;
+      t1 = tmp;
+    }
+    t_min = t0 > t_min ? t0 : t_min;
+    t_max = t1 < t_max ? t1 : t_max;
+    if (t_max <= t_min)
+      return false;
+  }
+  /* Z */
+  {
+    double inv_d = 1.0 / r.direction.z;
+    double t0 = (aabb->min.z - r.origin.z) * inv_d;
+    double t1 = (aabb->max.z - r.origin.z) * inv_d;
+    if (inv_d < 0.0) {
+      double tmp = t0;
+      t0 = t1;
+      t1 = tmp;
+    }
+    t_min = t0 > t_min ? t0 : t_min;
+    t_max = t1 < t_max ? t1 : t_max;
+    if (t_max <= t_min)
+      return false;
+  }
+
+  return true;
+}
+
+AABB aabb_surrounding_box(const AABB *first, const AABB *second) {
+  return (AABB){
+      .min = {fmin(first->min.x, second->min.x),
+              fmin(first->min.y, second->min.y),
+              fmin(first->min.z, second->min.z)},
+      .max = {fmax(first->max.x, second->max.x),
+              fmax(first->max.y, second->max.y),
+              fmax(first->max.z, second->max.z)},
+  };
+}

aabb.h

@@ -0,0 +1,17 @@
+#ifndef INCLUDED_AABB_H
+#define INCLUDED_AABB_H
+
+#include "point3.h"
+#include "ray.h"
+
+#include <stdbool.h>
+
+typedef struct AABB {
+  Point3 min;
+  Point3 max;
+} AABB;
+
+bool aabb_hit(const AABB *aabb, Ray r, double t_min, double t_max);
+AABB aabb_surrounding_box(const AABB *first, const AABB *second);
+
+#endif /* INCLUDED_AABB_H */

hittable.c

@@ -1,6 +1,8 @@
 #include "hittable.h"
+#include "aabb.h"
 #include "arena.h"
 #include "point3.h"
+#include "utils.h"
 #include "vec3.h"
 
 #include <assert.h>
@@ -25,10 +27,37 @@ bool hittable_hit(const Hittable *hittable, Ray r, double t_min, double t_max,
   case HITTABLE_MOVING_SPHERE:
     return moving_sphere_hit((const MovingSphere *)hittable, r, t_min, t_max,
                              record);
+  case HITTABLE_BVH_NODE:
+    return bvh_node_hit((const BVHNode *)hittable, r, t_min, t_max, record);
   }
   return false;
 }
 
+bool hittable_bounding_box(const Hittable *hittable, double time_start,
+                           double time_end, AABB *bounding_box) {
+  switch (hittable->type) {
+  case HITTABLE_LIST:
+    return hittable_list_bounding_box((const HittableList *)hittable,
+                                      time_start, time_end, bounding_box);
+  case HITTABLE_SPHERE:
+    return sphere_bounding_box((const Sphere *)hittable, time_start, time_end,
+                               bounding_box);
+  case HITTABLE_MOVING_SPHERE:
+    return moving_sphere_bounding_box((const MovingSphere *)hittable,
+                                      time_start, time_end, bounding_box);
+  case HITTABLE_BVH_NODE:
+    return bvh_node_bounding_box((const BVHNode *)hittable, time_start,
+                                 time_end, bounding_box);
+  }
+  return false;
+}
+
+HittableList *hittable_list_create(Arena *arena) {
+  HittableList *list = arena_alloc(arena, sizeof(HittableList));
+  list->type = HITTABLE_LIST;
+  return list;
+}
+
 static void hittable_list_grow(HittableList *list, size_t n, Arena *arena) {
   if (list->objects) {
     const Hittable **new_objects =
@@ -65,6 +94,26 @@ bool hittable_list_hit(const HittableList *list, Ray r, double t_min,
   return hit_anything;
 }
 
+bool hittable_list_bounding_box(const HittableList *list, double time_start,
+                                double time_end, AABB *bounding_box) {
+  if (list->size == 0)
+    return false;
+
+  AABB temp_box;
+  bool first_box = true;
+
+  for (size_t i = 0; i < list->size; ++i) {
+    if (!hittable_bounding_box(list->objects[i], time_start, time_end,
+                               &temp_box))
+      return false;
+    *bounding_box =
+        first_box ? temp_box : aabb_surrounding_box(bounding_box, &temp_box);
+    first_box = false;
+  }
+
+  return true;
+}
+
 Sphere *sphere_create(Point3 center, double radius, const Material *material,
                       Arena *arena) {
   Sphere *sphere = arena_alloc(arena, sizeof(Sphere));
@@ -102,6 +151,18 @@ bool sphere_hit(const Sphere *sphere, Ray r, double t_min, double t_max,
   return true;
 }
 
+bool sphere_bounding_box(const Sphere *sphere, double time_start,
+                         double time_end, AABB *bounding_box) {
+  (void)time_start, (void)time_end;
+  *bounding_box = (AABB){
+      .min = point3_add(sphere->center, (Vec3){-sphere->radius, -sphere->radius,
+                                               -sphere->radius}),
+      .max = point3_add(sphere->center,
+                        (Vec3){sphere->radius, sphere->radius, sphere->radius}),
+  };
+  return true;
+}
+
 MovingSphere *moving_sphere_create(Point3 center_start, Point3 center_end,
                                    double start, double end, double radius,
                                    const Material *material, Arena *arena) {
@@ -150,3 +211,102 @@ bool moving_sphere_hit(const MovingSphere *sphere, Ray r, double t_min,
   record->material = sphere->material;
   return true;
 }
+
+bool moving_sphere_bounding_box(const MovingSphere *sphere, double time_start,
+                                double time_end, AABB *bounding_box) {
+  AABB box_start = {
+      .min =
+          point3_add(moving_sphere_center(sphere, time_start),
+                     (Vec3){-sphere->radius, -sphere->radius, -sphere->radius}),
+      .max = point3_add(moving_sphere_center(sphere, time_start),
+                        (Vec3){sphere->radius, sphere->radius, sphere->radius}),
+  };
+  AABB box_end = {
+      .min =
+          point3_add(moving_sphere_center(sphere, time_end),
+                     (Vec3){-sphere->radius, -sphere->radius, -sphere->radius}),
+      .max = point3_add(moving_sphere_center(sphere, time_end),
+                        (Vec3){sphere->radius, sphere->radius, sphere->radius}),
+  };
+  *bounding_box = aabb_surrounding_box(&box_start, &box_end);
+  return true;
+}
+
+typedef int BoxCompareFunc(const void *lhs, const void *rhs);
+
+#define BOX_COMPARATOR(axis)                                                   \
+  static int box_##axis##_compare(const void *lhs, const void *rhs) {          \
+    AABB lhs_box, rhs_box;                                                     \
+    if (!hittable_bounding_box(*(const Hittable **)lhs, 0, 0, &lhs_box) ||     \
+        !hittable_bounding_box(*(const Hittable **)rhs, 0, 0, &rhs_box)) {     \
+      fprintf(stderr, "No bounding-box in BVH node");                          \
+      exit(1);                                                                 \
+    }                                                                          \
+                                                                               \
+    return lhs_box.min.axis < rhs_box.min.axis;                                \
+  }
+BOX_COMPARATOR(x)
+BOX_COMPARATOR(y)
+BOX_COMPARATOR(z)
+#undef BOX_COMPARATOR
+
+BVHNode *bvh_node_create(const Hittable **objects, size_t start, size_t end,
+                         double time_start, double time_end, Arena *arena) {
+  BVHNode *node = arena_alloc(arena, sizeof(BVHNode));
+  node->type = HITTABLE_BVH_NODE;
+
+  int axis = random_int_in_range(0, 2);
+  BoxCompareFunc *comparator = (axis == 0)   ? box_x_compare
+                               : (axis == 1) ? box_y_compare
+                                             : box_z_compare;
+
+  size_t object_span = end - start;
+  if (object_span == 1) {
+    node->left = node->right = objects[start];
+  } else if (object_span == 2) {
+    if (comparator(&objects[start], &objects[start + 1])) {
+      node->left = objects[start];
+      node->right = objects[start + 1];
+    } else {
+      node->left = objects[start + 1];
+      node->right = objects[start];
+    }
+  } else {
+    qsort(objects + start, object_span, sizeof(const Hittable *), comparator);
+    size_t mid = start + object_span / 2;
+    node->left = (const Hittable *)bvh_node_create(objects, start, mid,
+                                                   time_start, time_end, arena);
+    node->right = (const Hittable *)bvh_node_create(
+        objects, mid, end, time_start, time_end, arena);
+  }
+
+  AABB left_box, right_box;
+
+  if (!hittable_bounding_box(node->left, time_start, time_end, &left_box) ||
+      !hittable_bounding_box(node->right, time_start, time_end, &right_box)) {
+    fprintf(stderr, "No bounding-box in BVH node");
+    exit(1);
+  }
+
+  node->box = aabb_surrounding_box(&left_box, &right_box);
+  return node;
+}
+
+bool bvh_node_hit(const BVHNode *node, Ray r, double t_min, double t_max,
+                  HitRecord *record) {
+  if (!aabb_hit(&node->box, r, t_min, t_max))
+    return false;
+
+  bool hit_left = hittable_hit(node->left, r, t_min, t_max, record);
+  bool hit_right =
+      hittable_hit(node->right, r, t_min, hit_left ? record->t : t_max, record);
+
+  return hit_left || hit_right;
+}
+
+bool bvh_node_bounding_box(const BVHNode *node, double time_start,
+                           double time_end, AABB *bounding_box) {
+  (void)time_start, (void)time_end;
+  *bounding_box = node->box;
+  return true;
+}

hittable.h

@@ -1,6 +1,7 @@
 #ifndef INCLUDED_HITTABLE_H
 #define INCLUDED_HITTABLE_H
 
+#include "aabb.h"
 #include "arena.h"
 #include "material.h"
 #include "point3.h"
@@ -24,6 +25,7 @@ typedef enum HittableType {
   HITTABLE_LIST,
   HITTABLE_SPHERE,
   HITTABLE_MOVING_SPHERE,
+  HITTABLE_BVH_NODE,
 } HittableType;
 
 typedef struct Hittable {
@@ -32,6 +34,8 @@ typedef struct Hittable {
 
 bool hittable_hit(const Hittable *hittable, Ray r, double t_min, double t_max,
                   HitRecord *record);
+bool hittable_bounding_box(const Hittable *hittable, double time_start,
+                           double time_end, AABB *bounding_box);
 
 typedef struct HittableList {
   HittableType type;
@@ -40,10 +44,13 @@ typedef struct HittableList {
   size_t capacity;
 } HittableList;
 
+HittableList *hittable_list_create(Arena *arena);
 void hittable_list_add(HittableList *list, const Hittable *hittable,
                        Arena *arena);
 bool hittable_list_hit(const HittableList *list, Ray r, double t_min,
                        double t_max, HitRecord *record);
+bool hittable_list_bounding_box(const HittableList *list, double time_start,
+                                double time_end, AABB *bounding_box);
 
 typedef struct Sphere {
   HittableType type;
@@ -56,6 +63,8 @@ Sphere *sphere_create(Point3 center, double radius, const Material *material,
                       Arena *arena);
 bool sphere_hit(const Sphere *sphere, Ray r, double t_min, double t_max,
                 HitRecord *record);
+bool sphere_bounding_box(const Sphere *sphere, double time_start,
+                         double time_end, AABB *bounding_box);
 
 typedef struct MovingSphere {
   HittableType type;
@@ -71,5 +80,21 @@ MovingSphere *moving_sphere_create(Point3 center_start, Point3 center_end,
 Point3 moving_sphere_center(const MovingSphere *sphere, double t);
 bool moving_sphere_hit(const MovingSphere *sphere, Ray r, double t_min,
                        double t_max, HitRecord *record);
+bool moving_sphere_bounding_box(const MovingSphere *sphere, double time_start,
+                                double time_end, AABB *bounding_box);
+
+typedef struct BVHNode {
+  HittableType type;
+  const Hittable *left;
+  const Hittable *right;
+  AABB box;
+} BVHNode;
+
+BVHNode *bvh_node_create(const Hittable **objects, size_t start, size_t end,
+                         double time_start, double time_end, Arena *arena);
+bool bvh_node_hit(const BVHNode *node, Ray r, double t_min, double t_max,
+                  HitRecord *record);
+bool bvh_node_bounding_box(const BVHNode *node, double time_start,
+                           double time_end, AABB *bounding_box);
 
 #endif /* INCLUDED_HITTABLE_H */

main.c

@@ -35,14 +35,14 @@ static Color ray_color(Ray r, const Hittable *world, int depth) {
 }
 
 static const Hittable *generate_random_scene(Arena *arena) {
-  static HittableList world = {.type = HITTABLE_LIST};
+  HittableList *world = hittable_list_create(arena);
 
   const Lambertian *ground_material =
       lambertian_create((Color){0.5, 0.5, 0.5}, arena);
   const Sphere *ground_sphere =
       sphere_create((Point3){0.0, -1000.0, 0.0}, 1000.0,
                     (const Material *)ground_material, arena);
-  hittable_list_add(&world, (const Hittable *)ground_sphere, arena);
+  hittable_list_add(world, (const Hittable *)ground_sphere, arena);
 
   const Dielectric *glass = dielectric_create(1.5, arena);
 
@@ -61,18 +61,18 @@ static const Hittable *generate_random_scene(Arena *arena) {
           const MovingSphere *sphere =
               moving_sphere_create(center, center2, 0.0, 1.0, 0.2,
                                    (const Material *)material, arena);
-          hittable_list_add(&world, (const Hittable *)sphere, arena);
+          hittable_list_add(world, (const Hittable *)sphere, arena);
         } else if (choose_material < 0.95) {
           const Metal *material =
               metal_create(color_random_in_range(0.5, 1),
                            random_double_in_range(0.5, 1.0), arena);
           const Sphere *sphere =
               sphere_create(center, 0.2, (const Material *)material, arena);
-          hittable_list_add(&world, (const Hittable *)sphere, arena);
+          hittable_list_add(world, (const Hittable *)sphere, arena);
         } else {
           const Sphere *sphere =
               sphere_create(center, 0.2, (const Material *)glass, arena);
-          hittable_list_add(&world, (const Hittable *)sphere, arena);
+          hittable_list_add(world, (const Hittable *)sphere, arena);
         }
       }
     }
@@ -84,18 +84,22 @@ static const Hittable *generate_random_scene(Arena *arena) {
 
   const Sphere *sphere1 = sphere_create((Point3){0.0, 1.0, 0.0}, 1.0,
                                         (const Material *)glass, arena);
-  hittable_list_add(&world, (const Hittable *)sphere1, arena);
+  hittable_list_add(world, (const Hittable *)sphere1, arena);
   const Sphere *sphere2 = sphere_create((Point3){-4.0, 1.0, 0.0}, 1.0,
                                         (const Material *)lambertian, arena);
-  hittable_list_add(&world, (const Hittable *)sphere2, arena);
+  hittable_list_add(world, (const Hittable *)sphere2, arena);
   const Sphere *sphere3 = sphere_create((Point3){4.0, 1.0, 0.0}, 1.0,
                                         (const Material *)metal, arena);
-  hittable_list_add(&world, (const Hittable *)sphere3, arena);
+  hittable_list_add(world, (const Hittable *)sphere3, arena);
 
-  return (const Hittable *)&world;
+  BVHNode *bvh_root =
+      bvh_node_create(world->objects, 0, world->size, 0.0, 1.0, arena);
+  return (const Hittable *)bvh_root;
 }
 
 int main(void) {
+  srand(42);
+
   /* Memory management */
 
   const unsigned int buffer_size = 1 * 1024 * 1024;
@@ -152,6 +156,7 @@ int main(void) {
   return 0;
 }
 
+#include "aabb.c"
 #include "arena.c"
 #include "camera.c"
 #include "color.c"

utils.c

@@ -11,6 +11,10 @@ double random_double_in_range(double min, double max) {
   return min + (max - min) * random_double();
 }
 
+int random_int_in_range(int min, int max) {
+  return (int)(random_double_in_range(min, max + 1));
+}
+
 double clamp(double x, double min, double max) {
   if (x < min)
     return min;

utils.h

@@ -5,6 +5,7 @@ double degrees_to_radians(double degrees);
 
 double random_double(void);
 double random_double_in_range(double min, double max);
+int random_int_in_range(int min, int max);
 
 double clamp(double x, double min, double max);