1
\$\begingroup\$

This is supposed to subdivide a n-dimensional aabb (axis aligned bounding box) into 2ⁿ equal pieces as if sliced by n planes that are perpendicular to only one axis and parallel to the rest of the axes at a time.

The code:

CGL_bool CGL_aabb_subdivide_ndim(CGL_int n, CGL_float* aabb_min, CGL_float* aabb_max, CGL_float* aabbs_min_out, CGL_float* aabbs_max_out)
{
    static CGL_float aabb_mid_pos[1024], aabb_axis_deltas[1024]; // I don't think we'll ever need more than 1024 dimensions
    CGL_float nth_dim_val_holder = 0.0f, *pm_a[2];
    CGL_int result_count = 1 << n, nth_dim_val_holder_index = 0;
    for (CGL_int i = 0; i < n; i++)
    {
        aabb_mid_pos[i] = (aabb_min[i] + aabb_max[i]) * 0.5f;
        aabb_axis_deltas[i] = aabb_mid_pos[i] - aabb_min[i];
    }

    for (CGL_int j = 0; j < n; j++)
    {
        pm_a[0] = &aabb_min[j]; pm_a[1] = &aabb_mid_pos[j];
        nth_dim_val_holder = *pm_a[nth_dim_val_holder_index = 0];
        for (CGL_int i = 0; i < result_count; i++)
        {
            aabbs_min_out[i * n + (n - j - 1)] = nth_dim_val_holder;
            if ((i - 1) % (1 << j) == 0) nth_dim_val_holder = *pm_a[nth_dim_val_holder_index = (nth_dim_val_holder_index + 1) % 2];
        }
    }

    for (CGL_int i = 0; i < result_count; i++)
    {
        for (CGL_int j = 0; j < n; j++)
        {
            aabbs_max_out[i * n + j] = aabbs_min_out[i * n + j] + aabb_axis_deltas[j];
        }
    }

    return CGL_TRUE;
}

Here CGL_float, CGL_int, etc are just standard C types renamed with typedefs.

aabb_min is the min point of the n-dim box (it's an array of n elements) Similar for aabb_max. aabb_min_out and aabb_max_out are pre-allocated float arrays of size (2ⁿ * n) to store the result aabb min and max values for all the resulting aabb.

I tested this using:

CGL_float test_for_n_dim(CGL_int n)
{
    CGL_float* aabb_minn = (CGL_float*)malloc(sizeof(CGL_float) * n); if (!aabb_minn) return;
    CGL_float* aabb_maxn = (CGL_float*)malloc(sizeof(CGL_float) * n); if (!aabb_maxn) return;
    for (CGL_int i = 0; i < n; i++) { aabb_minn[i] = 0; aabb_maxn[i] = 1; }
    CGL_float* aabb_out_minn = (CGL_float*)malloc(sizeof(CGL_float) * n * (1 << n)); if (!aabb_out_minn) return;
    CGL_float* aabb_out_maxn = (CGL_float*)malloc(sizeof(CGL_float) * n * (1 << n)); if (!aabb_out_maxn) return;
    CGL_float time_before = CGL_utils_get_time();
    CGL_aabb_subdivide_ndim(n, aabb_minn, aabb_maxn, aabb_out_minn, aabb_out_maxn);
    CGL_float time_after = CGL_utils_get_time();
    if (n <= 3)
    {
        CGL_info("For %dD", n);
        for (CGL_int i = 0; i < (1 << n); i++)
        {
            CGL_printf_green("{(");
            for (CGL_int j = 0; j < n; j++)
            {
                CGL_printf_green("%f", aabb_out_minn[i * n + j]);
                if (j != n - 1) CGL_printf_green(", ");
            }
            CGL_printf_green(") -> (");
            for (CGL_int j = 0; j < n; j++)
            {
                CGL_printf_green("%f", aabb_out_maxn[i * n + j]);
                if (j != n - 1) CGL_printf_green(", ");
            }
            CGL_printf_green(")}\n");
        }
    }
    free(aabb_minn); free(aabb_maxn);free(aabb_out_minn); free(aabb_out_maxn);
    return time_after - time_before;
}

Results:

[INTERNAL] [Tue May  9 22:04:28 202] : Started Logger Session
[INFO] [Tue May  9 22:04:28 202] : Testing NDIM AABB Subdivision
[INFO] [Tue May  9 22:04:28 202] : For 1D
{(0.000000) -> (0.500000)}
{(0.500000) -> (1.000000)}
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000002 (1D)
[INFO] [Tue May  9 22:04:28 202] : For 2D
{(0.000000, 0.000000) -> (0.500000, 0.500000)}
{(0.000000, 0.500000) -> (0.500000, 1.000000)}
{(0.500000, 0.000000) -> (1.000000, 0.500000)}
{(0.500000, 0.500000) -> (1.000000, 1.000000)}
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000000 (2D)
[INFO] [Tue May  9 22:04:28 202] : For 3D
{(0.000000, 0.000000, 0.000000) -> (0.500000, 0.500000, 0.500000)}
{(0.000000, 0.000000, 0.500000) -> (0.500000, 0.500000, 1.000000)}
{(0.500000, 0.500000, 0.000000) -> (1.000000, 1.000000, 0.500000)}
{(0.500000, 0.500000, 0.500000) -> (1.000000, 1.000000, 1.000000)}
{(0.500000, 0.000000, 0.000000) -> (1.000000, 0.500000, 0.500000)}
{(0.500000, 0.000000, 0.500000) -> (1.000000, 0.500000, 1.000000)}
{(0.000000, 0.500000, 0.000000) -> (0.500000, 1.000000, 0.500000)}
{(0.000000, 0.500000, 0.500000) -> (0.500000, 1.000000, 1.000000)}
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000000 (3D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000001 (4D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000001 (5D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000001 (6D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000002 (7D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000005 (8D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000010 (9D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000030 (10D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000091 (11D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000207 (12D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.000394 (13D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.001144 (14D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.001890 (15D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.004921 (16D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.009553 (17D)
[INFO] [Tue May  9 22:04:28 202] : Time taken: 0.092833 (18D)
[INFO] [Tue May  9 22:04:29 202] : Time taken: 0.201626 (19D)
[INFO] [Tue May  9 22:04:29 202] : Time taken: 0.657279 (20D)
[INFO] [Tue May  9 22:04:30 202] : Time taken: 1.104113 (21D)
[INFO] [Tue May  9 22:04:32 202] : Time taken: 1.859703 (22D)
[INFO] [Tue May  9 22:04:36 202] : Time taken: 3.357090 (23D)
[INFO] [Tue May  9 22:04:43 202] : Time taken: 7.311569 (24D)
[INFO] [Tue May  9 22:05:09 202] : Time taken: 25.260529 (25D)
\$\endgroup\$
0

2 Answers 2

2
\$\begingroup\$

Please don't hide important code off the right-hand side of the screen like this:

    CGL_float* aabb_minn = (CGL_float*)malloc(sizeof(CGL_float) * n); if (!aabb_minn) return;
    CGL_float* aabb_maxn = (CGL_float*)malloc(sizeof(CGL_float) * n); if (!aabb_maxn) return;

One line per statement, please! And make the return value consistent with the function's declaration - getting away with this is a clear sign that you didn't enable enough warning options when you compiled.

There's no need to cast the void* which malloc() returns. Also, it's generally clearer to use the assigned-to pointer in the sizeof expression rather than its type (this matters more when the declaration is far away):

    CGL_float* aabb_minn = malloc((sizeof *aabb_min) * n);

Did you see the bug in these lines? When the allocation for aabb_maxn fails, we leak the memory pointed to by aabb_minn. That's why we don't hide the clean-up code off the right-hand margin!

We really need something more like

    CGL_float *const aabb_minn = malloc((sizeof *aabb_minn) * n);
    CGL_float *const aabb_maxn = malloc((sizeof *aabb_maxn) * n);
    if (!aabb_minn || !aabb_maxn) {
        free(aabb_minn);
        free(aabb_maxn);
        return 0.0/0.0;  /* NaN */
    }
\$\endgroup\$
6
  • \$\begingroup\$ The cast is so that this can be used from c++ too. Also i was wondering if there are any possible optimization in th algorithm itself \$\endgroup\$ Commented May 10, 2023 at 8:04
  • \$\begingroup\$ Also i didn't bother checking if malloc failed properly as I was just testing the algorithm with some thow away code and that's not anything important \$\endgroup\$ Commented May 10, 2023 at 8:14
  • \$\begingroup\$ In C++, prefer not to use malloc() - use a smart pointer instead (so you never need to explicitly free()). It's usually a mistake to share source code between languages like that (only do so for declarations that are accessed from both languages). \$\endgroup\$ Commented May 10, 2023 at 8:21
  • \$\begingroup\$ This code is originally for c. But i want it to be able to compile with c++ too so. I guess a extra typecast doesn't hurt much. \$\endgroup\$ Commented May 10, 2023 at 8:31
  • \$\begingroup\$ Still a bad idea. If you need it for programs in both languages, write it in C, and use extern "C" to link it into your C++ program. IMO, that's better than writing code that's simultaneously poor C and poor C++. \$\endgroup\$ Commented May 10, 2023 at 8:54
2
\$\begingroup\$

In The code, there is no indication what CGL_aabb_subdivide_ndim() is supposed to do beyond function and parameter names.
This makes it impossible to tell happenstance from requirement.
(Example: Is the order of result coordinates of any significance?)

// I don't think we'll ever need more than 1024 dimensions
So #define MAX_DIMENSIONS 1024, code a check and document behaviour in case n < 0, n = 0, n = 1, MAX_DIMENSIONS < n - as well as CGL_aabb_subdivide_ndim() being non-reentrant.
(1024 dimensions would mean about \$2^{1037}\$ bytes output.)

You assign 1 << n (not 1LL << n) to a variable of circumscribed type CGL_int.
That will not work from n = #bits in int - 1 up.

In the "main loop" assigning values to aabbs_min_out elements, why take and dereference addresses of values instead of just manipulating the values?

    for (int dim = 0 ; dim < n ; dim++) {
        CGL_float 
            lower[] = { aabb_min[dim], aabb_min[dim] + aabb_max[dim]) * 0.5f },
            upper[] = { lower[1], aabb_max[dim] };
        for (int box = 0 ; box < result_count ; box++) {
            aabbs_min_out[box * n + (n - dim - 1)] = lower[(box >> dim) & 1];
            aabbs_max_out[box * n + (n - dim - 1)] = upper[(box >> dim) & 1];
        }
    }

Walking aabbs_min_out in steps of n isn't cache friendly; populating aabbs_max_out in this same loop doubles the pressure, worse if the memory system can't handle two address streams concurrently.

Your last loop handles coordinates for a single box in order of dimensions, while your main loop assigned mins in reverse order.
I can't imagine that to work - test with different sizes and offsets for each dimension.

Maybe better

    static const int MAX_DIMENSIONS = sizeof(size_t) * NBBY - sizeof(CGL_float);
    if (n <= 0 || MAX_DIMENSIONS < n)
        handle()
    CGL_float aabb_mid_pos[MAX_DIMENSIONS];
    for (int dim = 0 ; dim < n ; dim++)
            aabb_mid_pos[dim] = (aabb_min[dim] + aabb_max[dim]) * 0.5f;
    CGL_float 
        *lower[] = { aabb_min, aabb_mid_pos },
        *upper[] = { aabb_mid_pos, aabb_max };
    for (size_t box = 0, coordinate_index = 0 ; box < result_count ; box++)
        for (int dim = 0, bit = box ; dim < n 
             ; dim++, coordinate_index++, bit >>= 1) {
            aabbs_min_out[coordinate_index] = lower[bit & 1][dim];
            aabbs_max_out[coordinate_index] = upper[bit & 1][dim];
        }
    }
\$\endgroup\$
1
  • \$\begingroup\$ sizeof(size_t) * NBBY - 1 - logb(sizeof(CGL_float))? \$\endgroup\$
    – greybeard
    Commented May 11, 2023 at 9:54

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.