Skip to content

compute overlap refac #994

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Draft
wants to merge 2 commits into
base: master
Choose a base branch
from
Draft

compute overlap refac #994

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
9de621f
compute overlap
nicolasaunai Apr 10, 2025
8d39efb
fix other overlap
nicolasaunai Apr 10, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
174 changes: 67 additions & 107 deletions pyphare/pyphare/pharesee/geometry.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -149,6 +149,41 @@ def periodicity_shifts(domain_box):
return shifts


def possible_periodic_shifts(box, domain_box):
from pyphare.core.box import shift

boxes = []
dim = domain_box.ndim

if dim == 1:
for ishift in (-1, 0, 1):
ioffset = ishift * domain_box.shape[0]
offset = ioffset

shifted = shift(box, offset)
boxes += [(offset, shifted)]
if dim == 2:
for ishift in (-1, 0, 1):
ioffset = ishift * domain_box.shape[0]
for jshift in (-1, 0, 1):
joffset = jshift * domain_box.shape[1]
offset = [ioffset, joffset]
shifted = shift(box, offset)
boxes += [(offset, shifted)]
if dim == 3:
for ishift in (-1, 0, 1):
ioffset = ishift * domain_box.shape[0]
for jshift in (-1, 0, 1):
joffset = jshift * domain_box.shape[1]
for kshift in (-1, 0, 1):
koffset = kshift * domain_box.shape[2]

offset = [ioffset, joffset, koffset]
shifted = shift(box, offset)
boxes += [(offset, shifted)]
return boxes


def compute_overlaps(patches, domain_box):
"""
returns a list of overlaps for all patch datas in given patches
Expand All @@ -173,119 +208,44 @@ def compute_overlaps(patches, domain_box):
- offset : tuple of two offsets by which the patchData ghost box is shifted
to compute the overlap box.
"""

# sorting the patches per origin allows to
# consider first and last as possible periodic overlaps
dim = domain_box.ndim

if dim == 1:
patches = sorted(patches, key=lambda p: p.origin.all())

zero_offset = [0] * dim if dim > 1 else 0

overlaps = []
zero_offset = [0] * domain_box.ndim if domain_box.ndim > 1 else 0

# first deal with intra domain overlaps
for ip, refPatch in enumerate(patches):
for cmpPatch in patches[ip + 1 :]:
# for two patches, compare patch_datas of the same quantity
for cmpPatch in patches[ip:]:

for ref_pdname, ref_pd in refPatch.patch_datas.items():
cmp_pd = cmpPatch.patch_datas[ref_pdname]

gb1 = ref_pd.ghost_box
gb2 = cmp_pd.ghost_box
overlap = gb1 * gb2

if overlap is not None:
# boxes indexes represent cells
# therefore for fields, we need to
# adjust the box. This essentially
# add 1 to upper in case field is on corners
# because upper corner can only be grabbed
# if box extends to upper+1 cell
# particles don't need that as they are contained
# in cells.
if ref_pd.quantity == "field":
overlap = toFieldBox(overlap, ref_pd)

overlaps.append(
{
"pdatas": (ref_pd, cmp_pd),
"patches": (refPatch, cmpPatch),
"box": overlap,
"offset": (zero_offset, zero_offset),
}
)

def append(ref_pd, cmp_pd, refPatch, cmpPatch, overlap, offset_tuple):
overlaps.append(
{
"pdatas": (ref_pd, cmp_pd),
"patches": (refPatch, cmpPatch),
"box": overlap,
"offset": offset_tuple,
}
)

def borders_per(patch):
return "".join(
[key for key, side in sides.items() if patch.box * side is not None]
)

sides = domain_border_ghost_boxes(domain_box, patches)
shifts = periodicity_shifts(domain_box)

# now dealing with border patches to see their patchdata overlap
if dim == 1 and len(patches) > 0:
patches = [patches[0], patches[-1]]

# filter out patches not near a border
borders_per_patch = {p: borders_per(p) for p in patches}
border_patches = [
p for p, in_sides in borders_per_patch.items() if len(in_sides) > 0
]

for patch_i, ref_patch in enumerate(border_patches):
in_sides = borders_per_patch[ref_patch]
assert in_sides in shifts

for ref_pdname, ref_pd in ref_patch.patch_datas.items():
for shift in shifts[in_sides]:
for cmp_patch in border_patches[
patch_i:
]: # patches can overlap with themselves
for cmp_pdname, cmp_pd in cmp_patch.patch_datas.items():
if cmp_pdname == ref_pdname:
gb1 = ref_pd.ghost_box
gb2 = cmp_pd.ghost_box

offset = np.asarray(shift)
overlap = gb1 * boxm.shift(gb2, -offset)

if overlap is not None:
other_ovrlp = boxm.shift(gb1, offset) * gb2
assert other_ovrlp is not None

if ref_pd.quantity == "field":
overlap = toFieldBox(overlap, ref_pd)
other_ovrlp = toFieldBox(other_ovrlp, ref_pd)

append(
ref_pd,
cmp_pd,
ref_patch,
cmp_patch,
overlap,
(zero_offset, (-offset).tolist()),
)
append(
ref_pd,
cmp_pd,
ref_patch,
cmp_patch,
other_ovrlp,
(offset.tolist(), zero_offset),
)
gb_ref = ref_pd.ghost_box
gb_cmp = cmp_pd.ghost_box

for offset, shifted_cmp in possible_periodic_shifts(gb_cmp, domain_box):
overlap = gb_ref * shifted_cmp
if overlap is not None and not np.all(
overlap.shape == gb_ref.shape
):
if ref_pd.quantity == "field":
overlap = toFieldBox(overlap, ref_pd)

overlaps.append(
{
"pdatas": (ref_pd, cmp_pd),
"patches": (refPatch, cmpPatch),
"box": overlap,
"offset": (zero_offset, offset),
}
)
if offset != zero_offset:
other_overlap = boxm.shift(gb_ref, -offset) * gb_cmp
overlaps.append(
{
"pdatas": (ref_pd, cmp_pd),
"patches": (refPatch, cmpPatch),
"box": other_overlap,
"offset": (-offset, zero_offset),
}
)

return overlaps

Expand Down
20 changes: 11 additions & 9 deletions pyphare/pyphare_tests/test_pharesee/test_geometry.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@ def setup_hierarchy(self, dim, interp_order, nbr_cells, refinement_boxes, **kwar
domain_size=domain_size,
cell_width=domain_size / nbr_cells,
refinement_boxes=refinement_boxes,
**kwargs
**kwargs,
)

# used for tests without ddt hierarchy overrides
Expand Down Expand Up @@ -78,18 +78,20 @@ def test_overlaps(self):
}

overlaps = hierarchy_overlaps(hierarchy)
tested_overlaps = {
ilvl: [
{"box": overlap["box"], "offset": overlap["offset"]}
for overlap in overlaps[ilvl]
]
for ilvl in range(len(hierarchy.patch_levels))
}

for ilvl, lvl in enumerate(hierarchy.patch_levels):
self.assertEqual(len(expected[ilvl]), len(overlaps[ilvl]))

for exp, actual in zip(expected[ilvl], overlaps[ilvl]):
act_box = actual["box"]
act_offset = actual["offset"]
exp_box = exp["box"]
exp_offset = exp["offset"]

self.assertEqual(act_box, exp_box)
self.assertEqual(act_offset, exp_offset)
for actual in tested_overlaps[ilvl]:
if actual not in expected[ilvl]:
raise AssertionError(f"Unexpected overlap: {actual}")

def test_touch_border(self):
hierarchy = self.basic_hierarchy()
Expand Down
30 changes: 18 additions & 12 deletions pyphare/pyphare_tests/test_pharesee/test_geometry_2d.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -114,15 +114,19 @@
)

level_overlaps = hierarchy_overlaps(hierarchy)
tested_overlaps = {
ilvl: [
{"box": overlap["box"], "offset": overlap["offset"]}
for overlap in level_overlaps[ilvl]
]
for ilvl in expected.keys()
}
for ilvl, lvl in enumerate(hierarchy.levels().items()):
if ilvl not in expected:
continue
self.assertEqual(len(expected[ilvl]), len(level_overlaps[ilvl]))
for exp, actual in zip(expected[ilvl], level_overlaps[ilvl]):
self.assertEqual(actual["box"], exp["box"])
self.assertTrue(
(np.asarray(actual["offset"]) == np.asarray(exp["offset"])).all()
)
for actual in tested_overlaps[ilvl]:
self.assertTrue(actual in expected[ilvl])

Check notice

Code scanning / CodeQL

Imprecise assert Note test

assertTrue(a in b) cannot provide an informative message. Using assertIn(a, b) instead will give more informative messages.

if 1 in level_overlaps:
fig = hierarchy.plot_2d_patches(
Expand Down Expand Up @@ -163,14 +167,16 @@
level_overlaps = hierarchy_overlaps(hierarchy)
ilvl = 0
overlap_boxes = []
tested_overlaps = {
0: [
{"box": overlap["box"], "offset": overlap["offset"]}
for overlap in level_overlaps[ilvl]
]
}

self.assertEqual(len(expected), len(level_overlaps[ilvl]))
for exp, actual in zip(expected, level_overlaps[ilvl]):
self.assertEqual(actual["box"], exp["box"])
self.assertTrue(
(np.asarray(actual["offset"]) == np.asarray(exp["offset"])).all()
)
overlap_boxes += [actual["box"]]
self.assertEqual(len(expected), len(tested_overlaps[ilvl]))
for actual in tested_overlaps[ilvl]:
self.assertTrue(actual in expected)

Check notice

Code scanning / CodeQL

Imprecise assert Note test

assertTrue(a in b) cannot provide an informative message. Using assertIn(a, b) instead will give more informative messages.

fig = hierarchy.plot_2d_patches(
ilvl,
Expand Down
Loading