From 0f66b4a9be474c2c57409cd9f87cbe6fea0f3a85 Mon Sep 17 00:00:00 2001 From: Jackson Urban <121484871+jurb33@users.noreply.github.com> Date: Tue, 25 Aug 2026 22:58:20 -0500 Subject: [PATCH 1/5] add error case in tests for issue #5402 --- yt/data_objects/tests/test_octree.py | 40 +++++++++++++++++++++++++++- 1 file changed, 39 insertions(+), 1 deletion(-) diff --git a/yt/data_objects/tests/test_octree.py b/yt/data_objects/tests/test_octree.py index 0c52ed337e..4e0b6cddee 100644 --- a/yt/data_objects/tests/test_octree.py +++ b/yt/data_objects/tests/test_octree.py @@ -1,8 +1,11 @@ import numpy as np from numpy.testing import assert_almost_equal, assert_equal -from yt.geometry.oct_container import OctreeContainer +from yt.geometry.oct_container import _ORDER_MAX, OctreeContainer +from yt.geometry.particle_oct_container import ParticleOctreeContainer +from yt.geometry.selection_routines import AlwaysSelector from yt.testing import fake_sph_grid_ds +from yt.utilities.lib.geometry_utils import get_morton_indices n_ref = 4 @@ -142,3 +145,38 @@ def test_num_zones_tuple(): assert oct_scalar is not None assert oct_tuple is not None assert oct_nonuniform is not None + + +def test_octcellindex_neighbours_num_zones(): + """ + Regression test for #5402: fill_octcellindex_neighbours hard-coded the + assumption that every oct holds 2x2x2 zones, so both the loop bounds and + the output buffer size were wrong whenever num_zones wasn't (2, 2, 2). + """ + DLE = np.array([0.0, 0.0, 0.0]) + DRE = np.array([8.0, 8.0, 8.0]) + dx = (DRE - DLE) / (2**_ORDER_MAX) + + # One particle per octant of the root oct: refines exactly one level deep. + centers = np.array( + [[x, y, z] for x in (2, 6) for y in (2, 6) for z in (2, 6)], dtype="float64" + ) + morton = get_morton_indices(np.floor((centers - DLE) / dx).astype("uint64")) + morton.sort() + + for nz in ((2, 2, 2), (2, 3, 4)): + octree = ParticleOctreeContainer((1, 1, 1), DLE, DRE, num_zones=nz) + octree.n_ref = 1 + octree.add(morton) + octree.finalize() + + nzones = nz[0] * nz[1] * nz[2] + n_per_oct = (nz[0] + 2) * (nz[1] + 2) * (nz[2] + 2) # +1 ghost zone each side + + selector = AlwaysSelector(None) + num_octs = selector.count_octs(octree, -1) + _, cell_inds = octree.fill_octcellindex_neighbours(selector) + + assert_equal(cell_inds.size, num_octs * n_per_oct) + assert cell_inds.min() >= 0 + assert cell_inds.max() <= nzones From 2f2e4c662801004f57b2f257cf32de89e6a946fd Mon Sep 17 00:00:00 2001 From: Jackson Urban <121484871+jurb33@users.noreply.github.com> Date: Tue, 25 Aug 2026 23:13:40 -0500 Subject: [PATCH 2/5] fix issue #5402 --- yt/data_objects/tests/test_octree.py | 2 ++ yt/geometry/oct_container.pyx | 14 +++++++++++--- yt/geometry/oct_visitors.pyx | 28 ++++++++++++++-------------- 3 files changed, 27 insertions(+), 17 deletions(-) diff --git a/yt/data_objects/tests/test_octree.py b/yt/data_objects/tests/test_octree.py index 4e0b6cddee..a7ae7a706d 100644 --- a/yt/data_objects/tests/test_octree.py +++ b/yt/data_objects/tests/test_octree.py @@ -177,6 +177,8 @@ def test_octcellindex_neighbours_num_zones(): num_octs = selector.count_octs(octree, -1) _, cell_inds = octree.fill_octcellindex_neighbours(selector) + # old code hard-coded 4**3=64 cells/oct; for nz=(2,3,4) it's really + # 4*5*6=120 - this assertion is what catches that assert_equal(cell_inds.size, num_octs * n_per_oct) assert cell_inds.min() >= 0 assert cell_inds.max() <= nzones diff --git a/yt/geometry/oct_container.pyx b/yt/geometry/oct_container.pyx index 981085b15a..410452f662 100644 --- a/yt/geometry/oct_container.pyx +++ b/yt/geometry/oct_container.pyx @@ -818,12 +818,17 @@ cdef class OctreeContainer: num_octs = selector.count_octs(self, domain_id) cdef NeighbourCellIndexVisitor visitor + cdef int n_per_oct cdef np.uint32_t[::1] cell_inds cdef np.int64_t[::1] oct_inds - cell_inds = np.full(num_octs*4**3, self.nz[0] * self.nz[1] * self.nz[2], dtype=np.uint32) - oct_inds = np.full(num_octs*4**3, -1, dtype=np.int64) + # must match the per-oct cell count that NeighbourCellIndexVisitor.visit() writes + n_per_oct = ((self.nz[0] + 2*n_ghost_zones) + * (self.nz[1] + 2*n_ghost_zones) + * (self.nz[2] + 2*n_ghost_zones)) + cell_inds = np.full(num_octs*n_per_oct, self.nz[0] * self.nz[1] * self.nz[2], dtype=np.uint32) + oct_inds = np.full(num_octs*n_per_oct, -1, dtype=np.int64) visitor = NeighbourCellIndexVisitor(self, -1, n_ghost_zones) visitor.cell_inds = cell_inds @@ -929,7 +934,10 @@ cdef class OctreeContainer: cdef int num_octs if num_cells < 0: num_octs = selector.count_octs(self, domain_id) - num_cells = num_octs * 4**3 + # must match the per-oct cell count that NeighbourCellVisitor.visit() writes + num_cells = num_octs * ((self.nz[0] + 2*n_ghost_zones) + * (self.nz[1] + 2*n_ghost_zones) + * (self.nz[2] + 2*n_ghost_zones)) cdef NeighbourCellVisitor visitor cdef np.ndarray[np.uint8_t, ndim=1] levels diff --git a/yt/geometry/oct_visitors.pyx b/yt/geometry/oct_visitors.pyx index e64cb6bff2..87876edf36 100644 --- a/yt/geometry/oct_visitors.pyx +++ b/yt/geometry/oct_visitors.pyx @@ -388,7 +388,7 @@ cdef class BaseNeighbourVisitor(OctVisitor): fcoords[i] += 1 elif fcoords[i] > 1: fcoords[i] -= 1 - local_oct &= (0 <= ishift[i] <= 1) + local_oct &= (0 <= ishift[i] < self.nz[i]) other_oct = not local_oct # Use octree to find neighbour @@ -412,7 +412,7 @@ cdef class BaseNeighbourVisitor(OctVisitor): # Index of neighbouring cell within its oct for i in range(3): - self.neigh_ind[i] = (ishift[i]) % 2 + self.neigh_ind[i] = (ishift[i]) % self.nz[i] self.other_oct = other_oct if other_oct: @@ -460,15 +460,14 @@ cdef class NeighbourCellIndexVisitor(BaseNeighbourVisitor): self.last = o.domain_ind - cdef int i0, i1 + cdef int i0 i0 = -self.n_ghost_zones - i1 = 2 + self.n_ghost_zones # Loop over cells in and directly around oct - for i in range(i0, i1): + for i in range(i0, self.nz[0] + self.n_ghost_zones): ishift[0] = i - for j in range(i0, i1): + for j in range(i0, self.nz[1] + self.n_ghost_zones): ishift[1] = j - for k in range(i0, i1): + for k in range(i0, self.nz[2] + self.n_ghost_zones): ishift[2] = k self.set_neighbour_info(o, ishift) @@ -480,7 +479,8 @@ cdef class NeighbourCellIndexVisitor(BaseNeighbourVisitor): neigh_cell_ind = self.neighbour_rind() else: neigh_domain_ind = -1 - neigh_cell_ind = 8 + # sentinel: one past the last valid cell index, i.e. "no such cell" + neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] self.cell_inds[self.index] = neigh_cell_ind self.domain_inds[self.index] = neigh_domain_ind @@ -506,15 +506,14 @@ cdef class NeighbourCellVisitor(BaseNeighbourVisitor): self.last = o.domain_ind - cdef int i0, i1 + cdef int i0 i0 = -self.n_ghost_zones - i1 = 2 + self.n_ghost_zones # Loop over cells in and directly around oct - for i in range(i0, i1): + for i in range(i0, self.nz[0] + self.n_ghost_zones): ishift[0] = i - for j in range(i0, i1): + for j in range(i0, self.nz[1] + self.n_ghost_zones): ishift[1] = j - for k in range(i0, i1): + for k in range(i0, self.nz[2] + self.n_ghost_zones): ishift[2] = k self.set_neighbour_info(o, ishift) @@ -532,7 +531,8 @@ cdef class NeighbourCellVisitor(BaseNeighbourVisitor): neigh_level = 255 neigh_domain = -1 neigh_file_ind = -1 - neigh_cell_ind = 8 + # sentinel: one past the last valid cell index, i.e. "no such cell" + neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] self.levels[self.index] = neigh_level self.file_inds[self.index] = neigh_file_ind From 7e320fd07a762112a28ea742bcd8148ed2599bbf Mon Sep 17 00:00:00 2001 From: Jackson Urban <121484871+jurb33@users.noreply.github.com> Date: Wed, 26 Aug 2026 21:26:30 -0500 Subject: [PATCH 3/5] add test commenting and comment function references for tracing --- yt/data_objects/tests/test_octree.py | 4 ++-- yt/geometry/oct_container.pyx | 4 ++-- yt/geometry/oct_visitors.pyx | 4 ++-- 3 files changed, 6 insertions(+), 6 deletions(-) diff --git a/yt/data_objects/tests/test_octree.py b/yt/data_objects/tests/test_octree.py index a7ae7a706d..de2b9a09ec 100644 --- a/yt/data_objects/tests/test_octree.py +++ b/yt/data_objects/tests/test_octree.py @@ -177,8 +177,8 @@ def test_octcellindex_neighbours_num_zones(): num_octs = selector.count_octs(octree, -1) _, cell_inds = octree.fill_octcellindex_neighbours(selector) - # old code hard-coded 4**3=64 cells/oct; for nz=(2,3,4) it's really - # 4*5*6=120 - this assertion is what catches that + #old oct_visitors/containers hardcoded 4**3=64 cells/oct; for nz=(2,3,4) it's really + # 4*5*6=120. this assertion catches that assert_equal(cell_inds.size, num_octs * n_per_oct) assert cell_inds.min() >= 0 assert cell_inds.max() <= nzones diff --git a/yt/geometry/oct_container.pyx b/yt/geometry/oct_container.pyx index 410452f662..3b4ef3e774 100644 --- a/yt/geometry/oct_container.pyx +++ b/yt/geometry/oct_container.pyx @@ -823,7 +823,7 @@ cdef class OctreeContainer: cdef np.uint32_t[::1] cell_inds cdef np.int64_t[::1] oct_inds - # must match the per-oct cell count that NeighbourCellIndexVisitor.visit() writes + #match per-oct cell count in NeighbourCellIndexVisitor.visit() n_per_oct = ((self.nz[0] + 2*n_ghost_zones) * (self.nz[1] + 2*n_ghost_zones) * (self.nz[2] + 2*n_ghost_zones)) @@ -934,7 +934,7 @@ cdef class OctreeContainer: cdef int num_octs if num_cells < 0: num_octs = selector.count_octs(self, domain_id) - # must match the per-oct cell count that NeighbourCellVisitor.visit() writes + #match per-oct cell count in NeighbourCellVisitor.visit() num_cells = num_octs * ((self.nz[0] + 2*n_ghost_zones) * (self.nz[1] + 2*n_ghost_zones) * (self.nz[2] + 2*n_ghost_zones)) diff --git a/yt/geometry/oct_visitors.pyx b/yt/geometry/oct_visitors.pyx index 87876edf36..e1a7822e0c 100644 --- a/yt/geometry/oct_visitors.pyx +++ b/yt/geometry/oct_visitors.pyx @@ -479,7 +479,7 @@ cdef class NeighbourCellIndexVisitor(BaseNeighbourVisitor): neigh_cell_ind = self.neighbour_rind() else: neigh_domain_ind = -1 - # sentinel: one past the last valid cell index, i.e. "no such cell" + #one past the last valid cell index, i.e. "no such cell" neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] self.cell_inds[self.index] = neigh_cell_ind @@ -531,7 +531,7 @@ cdef class NeighbourCellVisitor(BaseNeighbourVisitor): neigh_level = 255 neigh_domain = -1 neigh_file_ind = -1 - # sentinel: one past the last valid cell index, i.e. "no such cell" + #one past the last valid cell index, i.e. "no such cell" neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] self.levels[self.index] = neigh_level From a7b5bbe7fc69567c8c3d37670f8924ac05f18a26 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 2 Sep 2026 08:19:00 +0000 Subject: [PATCH 4/5] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- yt/data_objects/tests/test_octree.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/yt/data_objects/tests/test_octree.py b/yt/data_objects/tests/test_octree.py index de2b9a09ec..ae113562c9 100644 --- a/yt/data_objects/tests/test_octree.py +++ b/yt/data_objects/tests/test_octree.py @@ -177,7 +177,7 @@ def test_octcellindex_neighbours_num_zones(): num_octs = selector.count_octs(octree, -1) _, cell_inds = octree.fill_octcellindex_neighbours(selector) - #old oct_visitors/containers hardcoded 4**3=64 cells/oct; for nz=(2,3,4) it's really + # old oct_visitors/containers hardcoded 4**3=64 cells/oct; for nz=(2,3,4) it's really # 4*5*6=120. this assertion catches that assert_equal(cell_inds.size, num_octs * n_per_oct) assert cell_inds.min() >= 0 From 6fcff742ab312ca96e4fa5192706a0ee3df3dd92 Mon Sep 17 00:00:00 2001 From: Corentin Cadiou Date: Wed, 2 Sep 2026 11:20:47 +0200 Subject: [PATCH 5/5] [STY] Cosmetic changes --- yt/geometry/oct_container.pyx | 20 ++++++++++++-------- yt/geometry/oct_visitors.pyx | 8 ++++---- 2 files changed, 16 insertions(+), 12 deletions(-) diff --git a/yt/geometry/oct_container.pyx b/yt/geometry/oct_container.pyx index 3b4ef3e774..74063ab4f5 100644 --- a/yt/geometry/oct_container.pyx +++ b/yt/geometry/oct_container.pyx @@ -823,10 +823,12 @@ cdef class OctreeContainer: cdef np.uint32_t[::1] cell_inds cdef np.int64_t[::1] oct_inds - #match per-oct cell count in NeighbourCellIndexVisitor.visit() - n_per_oct = ((self.nz[0] + 2*n_ghost_zones) - * (self.nz[1] + 2*n_ghost_zones) - * (self.nz[2] + 2*n_ghost_zones)) + # Match per-oct cell count in NeighbourCellIndexVisitor.visit() + n_per_oct = ( + (self.nz[0] + 2 * n_ghost_zones) + * (self.nz[1] + 2 * n_ghost_zones) + * (self.nz[2] + 2 * n_ghost_zones) + ) cell_inds = np.full(num_octs*n_per_oct, self.nz[0] * self.nz[1] * self.nz[2], dtype=np.uint32) oct_inds = np.full(num_octs*n_per_oct, -1, dtype=np.int64) @@ -934,10 +936,12 @@ cdef class OctreeContainer: cdef int num_octs if num_cells < 0: num_octs = selector.count_octs(self, domain_id) - #match per-oct cell count in NeighbourCellVisitor.visit() - num_cells = num_octs * ((self.nz[0] + 2*n_ghost_zones) - * (self.nz[1] + 2*n_ghost_zones) - * (self.nz[2] + 2*n_ghost_zones)) + # Match per-oct cell count in NeighbourCellVisitor.visit() + num_cells = num_octs * ( + (self.nz[0] + 2 * n_ghost_zones) + * (self.nz[1] + 2 * n_ghost_zones) + * (self.nz[2] + 2 * n_ghost_zones) + ) cdef NeighbourCellVisitor visitor cdef np.ndarray[np.uint8_t, ndim=1] levels diff --git a/yt/geometry/oct_visitors.pyx b/yt/geometry/oct_visitors.pyx index e1a7822e0c..91b0f9e300 100644 --- a/yt/geometry/oct_visitors.pyx +++ b/yt/geometry/oct_visitors.pyx @@ -461,6 +461,7 @@ cdef class NeighbourCellIndexVisitor(BaseNeighbourVisitor): self.last = o.domain_ind cdef int i0 + cdef int out_of_bound_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] i0 = -self.n_ghost_zones # Loop over cells in and directly around oct for i in range(i0, self.nz[0] + self.n_ghost_zones): @@ -479,8 +480,7 @@ cdef class NeighbourCellIndexVisitor(BaseNeighbourVisitor): neigh_cell_ind = self.neighbour_rind() else: neigh_domain_ind = -1 - #one past the last valid cell index, i.e. "no such cell" - neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] + neigh_cell_ind = out_of_bound_cell_ind self.cell_inds[self.index] = neigh_cell_ind self.domain_inds[self.index] = neigh_domain_ind @@ -507,6 +507,7 @@ cdef class NeighbourCellVisitor(BaseNeighbourVisitor): self.last = o.domain_ind cdef int i0 + cdef int out_of_bound_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] i0 = -self.n_ghost_zones # Loop over cells in and directly around oct for i in range(i0, self.nz[0] + self.n_ghost_zones): @@ -531,8 +532,7 @@ cdef class NeighbourCellVisitor(BaseNeighbourVisitor): neigh_level = 255 neigh_domain = -1 neigh_file_ind = -1 - #one past the last valid cell index, i.e. "no such cell" - neigh_cell_ind = self.nz[0] * self.nz[1] * self.nz[2] + neigh_cell_ind = out_of_bound_cell_ind self.levels[self.index] = neigh_level self.file_inds[self.index] = neigh_file_ind