Compare commits

...
Author SHA1 Message Date
ClaudeandClaude Opus 5.5 5f69ebaa90 run_refined_dfg: pick the coarse table builder by the imported source's API
jet_scenario=none host tests run against the current checkout (no
blocked_base_tables), while historical runs import the pinned source (no
blocked_coarse_tables). Select whichever single name the source defines;
fail if neither or both. No src alias; pins and manifests unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 05:11:34 +00:00
ClaudeandClaude Opus 5.5 9d0588dcad Merge codex/src-coherence-callers (7b652b2) into claude/src-coherence (ae53fad)
Conflicts (6 files): tests took the caller-migration side, then the
358bed1 renames (_level_boundaries -> _device_boundaries); run_refined_dfg.py
kept ae53fad (same jet_interval level= repair plus blocked_coarse_tables on the
current-checkout path). Follow-up renames in 5 auto-merged tests;
blocked_fine_tables imported from refinement_tables. Frozen drivers untouched
(run_jet_circle_cv.py, check_projected_pairs.py remap_child, historical DFG path).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 04:45:48 +00:00
ClaudeandClaude Opus 5.5 ae53fade13 Contract: name the in-EDGE half-cell projection; repair DFG runner current-checkout Jet path
- No NODE<->EDGE input conversion still leaves the root-edge -> fine-local
  half-cell projection (prepare_fine_boundary, fine_hits_to_edge_frame); say so.
- run_refined_dfg.py (jet_scenario != none, current checkout only): pass
  level=1 to FixedPatchStepper.set_jet_amplitude on reversal; use
  blocked_coarse_tables there. jet_scenario=none keeps the pinned source's
  blocked_base_tables; manifests and results unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 02:45:30 +00:00
ClaudeandClaude Opus 5.5 358bed1962 Naming/contract delta: _root_cfg, _device_boundaries, blocked_coarse_tables; NODE vs EDGE frame
- Coarse/fine name the two members of any pair (Hierarchy docstring); only the
  misleading blocked_base_tables -> blocked_coarse_tables. prepare_fine_boundary
  and blocked_fine_tables keep their pair-relative names.
- Depth-1 static circles use the integer-NODE frame (unchanged behavior); the
  root EDGE frame applies to refined levels only. Contract text corrected.
- _kernel_configs documents why it stays distinct from Hierarchy.configs.
- Drop the stale re-export comment; two tests import table builders from
  refinement_tables. Mechanical caller renames for this delta only.
No numerical code or config derivation changed.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 02:34:07 +00:00
23 changed files with 140 additions and 102 deletions

No files matched your search

+19 -4
View File
@@ -8,8 +8,16 @@ pending a GPU run).
- Level `l = 0 .. depth-1` (0 = root). Pair `hierarchy.pairs[l-1]` builds level `l`. - Level `l = 0 .. depth-1` (0 = root). Pair `hierarchy.pairs[l-1]` builds level `l`.
- Array shapes are `(ny, nx)`; rectangles, origins and centers are `(x, y)`. - Array shapes are `(ny, nx)`; rectangles, origins and centers are `(x, y)`.
- FixedPatchStepper bodies/Jets use the root EDGE frame (root cell `i` = `[i, i+1)`); - Coarse/fine name the two members (parent/child) of ANY pair, not absolute L0/L1:
uniform `Simulation` objects use the integer-NODE frame. Nothing converts between them. `fine_shape`, `fine_origin`, `prepare_fine_boundary`, `blocked_fine_tables`,
`blocked_coarse_tables`.
- Refined levels (l >= 1) of FixedPatchStepper place bodies/Jets in the root EDGE frame
(root cell `i` = `[i, i+1)`). Uniform `Simulation` objects AND the static circles of a
single-level (depth 1) FixedPatchStepper use the integer-NODE frame. Nothing converts
between them. This concerns input placement only: within the EDGE frame, refined
levels still project root-edge coordinates to fine-local cell indices with the
half-cell shift (`prepare_fine_boundary`: `(v - o)/dx - .5`; hits back via `fine_hits_to_edge_frame`
`origin + dx*(hit + .5)`), using `edge_frame`/`fine_origin` below.
- `Hierarchy.edge_frame(l)` -> `((ox, oy), h)`: cell `(i, j)` of level `l`'s allocated - `Hierarchy.edge_frame(l)` -> `((ox, oy), h)`: cell `(i, j)` of level `l`'s allocated
array (ghost band and holes included) covers `[o + h*i, o + h*(i+1))`. A scalar array (ghost band and holes included) covers `[o + h*i, o + h*(i+1))`. A scalar
FixedPatch2D window carries the half-coarse-cell ghost shift `bounds - HALO/2` in FixedPatch2D window carries the half-coarse-cell ghost shift `bounds - HALO/2` in
@@ -26,8 +34,11 @@ pending a GPU run).
| `base_layout`, `fine_layout`, `child_layout` (stepper properties) | `layouts[0]`, `layouts[1]`, `layouts[2]` (`FixedPatchSet2D.fine_layout(block)` is unchanged) | | `base_layout`, `fine_layout`, `child_layout` (stepper properties) | `layouts[0]`, `layouts[1]`, `layouts[2]` (`FixedPatchSet2D.fine_layout(block)` is unchanged) |
| `fine_boundary`, `child_boundary` | `boundaries[1]`, `boundaries[2]` (host geometry; `None` = no bodies/Jets) | | `fine_boundary`, `child_boundary` | `boundaries[1]`, `boundaries[2]` (host geometry; `None` = no bodies/Jets) |
| `_jets` | `jets[1]` | | `_jets` | `jets[1]` |
| `_boundaries[l]` (dict from L1, level keys) | `_level_boundaries[l]` (list from L0, `None` at L0 unless single-level static circles) | | `_boundaries[l]` (dict from L1, level keys) | `_device_boundaries[l]` (list from L0, `None` at L0 unless single-level static circles) |
| `_fine_tables`, `_actuated`, `_jet_cache` | `_level_boundaries[1].tables` / `.actuated` / `.jet_cache` | | `_level_boundaries` (private, b4eacf7) | `_device_boundaries` |
| `_fine_tables`, `_actuated`, `_jet_cache` | `_device_boundaries[1].tables` / `.actuated` / `.jet_cache` |
| `_remap_cfg` (private) | `_root_cfg` |
| `refinement_tables.blocked_base_tables` (also via `refinement_stepper`) | `refinement_tables.blocked_coarse_tables` (import from `refinement_tables`; the stepper no longer re-exports) |
| `_fine_leaf`, `_child_leaf` | `leaf[1]`, `leaf[2]` (copies) | | `_fine_leaf`, `_child_leaf` | `leaf[1]`, `leaf[2]` (copies) |
| `_middle_kind` | `hierarchy.kinds[1]` | | `_middle_kind` | `hierarchy.kinds[1]` |
| `_boundary_hosts` (private) | `_boundaries` (private list behind `boundaries`) | | `_boundary_hosts` (private) | `_boundaries` (private list behind `boundaries`) |
@@ -56,6 +67,10 @@ the replaced level relative to its parent).
## Not done / open ## Not done / open
- Two per-level config derivations stay distinct on purpose: `Hierarchy.configs`
(`_middle_config`, wall fields reset on every refined level, geometry validation) vs.
`FixedPatchStepper._kernel_configs` (reset only on level 1 of a slot_reservoir root;
generates each level's config.h). Folding them needs a config.h/PTX diff gate.
- No single-`Simulation` facade over FixedPatchStepper: the two own different state models - No single-`Simulation` facade over FixedPatchStepper: the two own different state models
(uniform node frame, body sync, forcing, checkpoint vs. immutable blocked hierarchy, edge (uniform node frame, body sync, forcing, checkpoint vs. immutable blocked hierarchy, edge
frame, remap-as-new-instance). A facade would only forward calls; no concrete need found. frame, remap-as-new-instance). A facade would only forward calls; no concrete need found.
+9 -2
View File
@@ -26,8 +26,15 @@ class Hierarchy:
Index conventions (the single place they are defined): Index conventions (the single place they are defined):
* level l = 0 (root) .. depth-1; pair index p = l-1 builds level l from level p; * level l = 0 (root) .. depth-1; pair index p = l-1 builds level l from level p;
* array shapes are (ny, nx); rectangles, origins and centers are (x, y); * array shapes are (ny, nx); rectangles, origins and centers are (x, y);
* the root EDGE frame puts root cell i at [i, i+1), centers at i+.5 (bodies and * a pair's two members are its COARSE (parent) and FINE (child) level, for every
Jets of FixedPatchStepper); uniform Simulation objects use the integer-NODE frame; pair, not absolute L0/L1 (fine_shape, fine_origin, blocked_coarse_tables, ...);
* refined levels (l >= 1) of FixedPatchStepper place bodies and Jets in the root
EDGE frame: root cell i covers [i, i+1), center i+.5. Uniform Simulation objects
and the static circles of a single-level (depth 1) FixedPatchStepper use the
integer-NODE frame (CircleGeometry.build_flag_mask at integer (x, y)). No
NODE<->EDGE input conversion is applied anywhere; inside the EDGE frame, refined
boundaries still project root-edge coordinates to fine-local cell indices with
the half-cell shift ((v - o)/h - .5, see prepare_fine_boundary);
* cell (i, j) of level l's allocated array (ghost band and holes included) covers * cell (i, j) of level l's allocated array (ghost band and holes included) covers
[o + h*i, o + h*(i+1)) with (o, h) = edge_frame(l); level 1 is built from the root [o + h*i, o + h*(i+1)) with (o, h) = edge_frame(l); level 1 is built from the root
config directly, deeper levels from their parent frame (see edge_frame). config directly, deeper levels from their parent frame (see edge_frame).
+27 -22
View File
@@ -22,7 +22,7 @@ from .channel import channel_masks, validate_channel_config, validate_channel_tr
from ..body.jets import jet_link_velocity from ..body.jets import jet_link_velocity
from .hierarchy import build_hierarchy from .hierarchy import build_hierarchy
from .refinement import ACTIVE, HALO, FixedPatch2D, FixedPatchSet2D, prepare_fine_boundary from .refinement import ACTIVE, HALO, FixedPatch2D, FixedPatchSet2D, prepare_fine_boundary
from .refinement_tables import (blocked_base_tables, blocked_fine_tables, # noqa: F401 (re-export) from .refinement_tables import (blocked_coarse_tables, blocked_fine_tables,
blocked_kind, projected_explosion_tables, _prepare_blocked_boundary, blocked_kind, projected_explosion_tables, _prepare_blocked_boundary,
_BlockedBoundary, _BlockedBoundary,
BLOCKED_BOUNDARY_DEVICE_COLUMNS) BLOCKED_BOUNDARY_DEVICE_COLUMNS)
@@ -66,8 +66,12 @@ class _Level:
def _kernel_configs(cfg, pairs, wall_aware): def _kernel_configs(cfg, pairs, wall_aware):
"""Per-level kernel LBMConfig (window shape, viscosity doubling per level). """Per-level kernel LBMConfig (window shape, viscosity doubling per level).
Geometry configs live in Hierarchy.configs; level 1 of a slot_reservoir root Deliberately distinct from Hierarchy.configs (refinement._middle_config), which
replaces the root wall fields by the internal bounce-back/velocity ghost setup. resets the wall fields on EVERY refined level for geometry validation. Here only
level 1 of a slot_reservoir root resets them (internal bounce-back/velocity ghost
setup); other levels inherit their parent's fields. These configs generate each
level's compiled config.h, so folding the two derivations is a numerical-parity
change (config.h/PTX diff gate), not a style one.
""" """
configs = [replace(cfg)] configs = [replace(cfg)]
for index, pair in enumerate(pairs): for index, pair in enumerate(pairs):
@@ -119,7 +123,8 @@ class FixedPatchStepper:
supported). supported).
``refinements=[]`` is one blocked level: periodic SRT root without LES, ``refinements=[]`` is one blocked level: periodic SRT root without LES,
Jets or exchange observation; its ``level_objects[0]`` are stationary Jets or exchange observation; its ``level_objects[0]`` are stationary
passive circles (static BFL). passive circles (static BFL) in the integer-NODE frame, like Simulation
objects (see Hierarchy for frames).
``coarse_boundary`` is ``periodic`` (default) or ``slot_reservoir``: ``coarse_boundary`` is ``periodic`` (default) or ``slot_reservoir``:
a bottom on-site slot and top fixed-density reservoir from ``cfg``. a bottom on-site slot and top fixed-density reservoir from ``cfg``.
@@ -131,7 +136,7 @@ class FixedPatchStepper:
Periodic mode does not apply the wall/slot boundary fields in ``cfg``. Periodic mode does not apply the wall/slot boundary fields in ``cfg``.
Outer walls remain coarse; fine halo donors must be ordinary fluid. Outer walls remain coarse; fine halo donors must be ordinary fluid.
Bodies are fixed circles in root (L0) physical coordinates; velocities are Bodies of refined levels are fixed circles in root (L0) physical coordinates; velocities are
in root lattice units (dx=dt=1). Independent local copies are projected at in root lattice units (dx=dt=1). Independent local copies are projected at
fine cell centers; their complete stencil must be active and disjoint from fine cell centers; their complete stencil must be active and disjoint from
the actual two-hop transfer footprint. These use the edge/cell-center frame the actual two-hop transfer footprint. These use the edge/cell-center frame
@@ -258,7 +263,7 @@ class FixedPatchStepper:
raise ValueError(f"level {level} blocks exceed launch limits") raise ValueError(f"level {level} blocks exceed launch limits")
self._layouts.append(layout) self._layouts.append(layout)
if self._channel: if self._channel:
base = blocked_base_tables(pairs[0], self._layouts[0], np.argsort(pairs[0].coalesce_coarse)) base = blocked_coarse_tables(pairs[0], self._layouts[0], np.argsort(pairs[0].coalesce_coarse))
validate_channel_transfers(pairs[0], self._layouts[0], base[1]) validate_channel_transfers(pairs[0], self._layouts[0], base[1])
self._level_cfgs = _kernel_configs(cfg, pairs, self._wall_aware) self._level_cfgs = _kernel_configs(cfg, pairs, self._wall_aware)
# Host boundary geometry per level in the root edge frame (Hierarchy.edge_frame). # Host boundary geometry per level in the root edge frame (Hierarchy.edge_frame).
@@ -294,7 +299,7 @@ class FixedPatchStepper:
if depth == 1 and level_objects[0]: if depth == 1 and level_objects[0]:
tables = static_circle_tables(self._layouts[0], level_objects[0]) tables = static_circle_tables(self._layouts[0], level_objects[0])
root = _BlockedBoundary(None, tables, (), False) root = _BlockedBoundary(None, tables, (), False)
self._level_boundaries = [root] + [ self._device_boundaries = [root] + [
_prepare_blocked_boundary(self._boundaries[level], self._layouts[level], _prepare_blocked_boundary(self._boundaries[level], self._layouts[level],
self._level_cfgs[level], level_jets[level], self._level_cfgs[level], level_jets[level],
_shifted_blocked=self._ddf_shifting) _shifted_blocked=self._ddf_shifting)
@@ -310,7 +315,7 @@ class FixedPatchStepper:
if self._channel: if self._channel:
self._leaf[0] &= channel_masks((cfg.ny, cfg.nx))[0] self._leaf[0] &= channel_masks((cfg.ny, cfg.nx))[0]
# Private constructor inputs for remap_finest. # Private constructor inputs for remap_finest.
self._remap_cfg = replace(cfg) self._root_cfg = replace(cfg)
self._level_blocks = level_blocks self._level_blocks = level_blocks
self.ctx = context self.ctx = context
self._lock = RLock() self._lock = RLock()
@@ -366,7 +371,7 @@ class FixedPatchStepper:
for index, layout in enumerate(self._layouts): for index, layout in enumerate(self._layouts):
total += 2 * 4 * int(np.prod(layout.blocked_shape)) total += 2 * 4 * int(np.prod(layout.blocked_shape))
total += layout.halo_dst.nbytes + layout.halo_src.nbytes total += layout.halo_dst.nbytes + layout.halo_src.nbytes
b = self._level_boundaries[index] b = self._device_boundaries[index]
if b is not None: if b is not None:
total += sum(np.asarray(v).nbytes for k, v in b.tables.items() total += sum(np.asarray(v).nbytes for k, v in b.tables.items()
if k == "flags" or k in BLOCKED_BOUNDARY_DEVICE_COLUMNS) if k == "flags" or k in BLOCKED_BOUNDARY_DEVICE_COLUMNS)
@@ -389,19 +394,19 @@ class FixedPatchStepper:
grid = ((config.nx + self._block[0] - 1) // self._block[0], config.ny, layout.nblocks) grid = ((config.nx + self._block[0] - 1) // self._block[0], config.ny, layout.nblocks)
level = _Level(shape, module, self._allocate(size), self._allocate(size), level = _Level(shape, module, self._allocate(size), self._allocate(size),
module.get_function(name), grid, layout) module.get_function(name), grid, layout)
level.boundary = self._level_boundaries[index] level.boundary = self._device_boundaries[index]
level.halo_copy = module.get_function("BlockHaloCopy") level.halo_copy = module.get_function("BlockHaloCopy")
if len(layout.halo_dst): if len(layout.halo_dst):
level.halo_dst = self._upload(layout.halo_dst) level.halo_dst = self._upload(layout.halo_dst)
level.halo_src = self._upload(layout.halo_src) level.halo_src = self._upload(layout.halo_src)
if index == 0 and self._channel: if index == 0 and self._channel:
cfg = self._remap_cfg cfg = self._root_cfg
ptr, size = module.get_global("d_block_channel") ptr, size = module.get_global("d_block_channel")
if size != struct.calcsize("=IIf"): if size != struct.calcsize("=IIf"):
raise RuntimeError("Block channel host/device layout mismatch") raise RuntimeError("Block channel host/device layout mismatch")
cuda.memcpy_htod(ptr, struct.pack("=IIf", cfg.nx, cfg.ny, cfg.velocity)) cuda.memcpy_htod(ptr, struct.pack("=IIf", cfg.nx, cfg.ny, cfg.velocity))
if index == 0 and self._wall_aware: if index == 0 and self._wall_aware:
cfg = self._remap_cfg cfg = self._root_cfg
ptr, size = module.get_global("d_block_wall") ptr, size = module.get_global("d_block_wall")
if size != struct.calcsize("=IIIf"): if size != struct.calcsize("=IIIf"):
raise RuntimeError("Block wall host/device layout mismatch") raise RuntimeError("Block wall host/device layout mismatch")
@@ -415,11 +420,11 @@ class FixedPatchStepper:
order = np.argsort(pair.coalesce_coarse) order = np.argsort(pair.coalesce_coarse)
ghosts, children, child_kind = blocked_fine_tables(pair, child.layout, order) ghosts, children, child_kind = blocked_fine_tables(pair, child.layout, order)
if index == 0: if index == 0:
donors, targets, parent_kind = blocked_base_tables(pair, parent.layout, order, donors, targets, parent_kind = blocked_coarse_tables(pair, parent.layout, order,
self._wall_aware) self._wall_aware)
else: else:
# Nested pairs keep inherited INNER/OUTER and sparse holes of the parent. # Nested pairs keep inherited INNER/OUTER and sparse holes of the parent.
donors, targets, parent_kind = blocked_base_tables( donors, targets, parent_kind = blocked_coarse_tables(
pair, parent.layout, order, composed_kind=self.hierarchy.kinds[index]) pair, parent.layout, order, composed_kind=self.hierarchy.kinds[index])
result = _Pair(child.module.get_function("BlockExplode"), result = _Pair(child.module.get_function("BlockExplode"),
parent.module.get_function("BlockCoalesce"), parent.module.get_function("BlockCoalesce"),
@@ -552,12 +557,12 @@ class FixedPatchStepper:
if (self._wall_aware or self._channel or self.hierarchy.depth != 3 if (self._wall_aware or self._channel or self.hierarchy.depth != 3
or not isinstance(pairs[1], FixedPatchSet2D) or not isinstance(pairs[1], FixedPatchSet2D)
or self._boundaries[1] is not None or self._boundaries[1] is not None
or any(b is not None for b in self._level_boundaries[:2])): or any(b is not None for b in self._device_boundaries[:2])):
raise ValueError("remap_finest requires a periodic three-level hierarchy with a " raise ValueError("remap_finest requires a periodic three-level hierarchy with a "
"level-2 set and no level-0/1 bodies") "level-2 set and no level-0/1 bodies")
cfg = self._remap_cfg cfg = self._root_cfg
objects = tuple(objects) objects = tuple(objects)
actual = self._level_boundaries[2] actual = self._device_boundaries[2]
jets = () if actual is None else actual.jets jets = () if actual is None else actual.jets
old = RemapTopology.from_config(cfg, pairs[0], pairs[1].bounds, old = RemapTopology.from_config(cfg, pairs[0], pairs[1].bounds,
child_objects=objects, child_jets=jets, child_objects=objects, child_jets=jets,
@@ -628,7 +633,7 @@ class FixedPatchStepper:
if (isinstance(level, bool) or not isinstance(level, Integral) if (isinstance(level, bool) or not isinstance(level, Integral)
or not 1 <= level < self.hierarchy.depth): or not 1 <= level < self.hierarchy.depth):
raise ValueError("Jet level must name an existing refined level (1 .. depth-1)") raise ValueError("Jet level must name an existing refined level (1 .. depth-1)")
boundary = self._level_boundaries[level] boundary = self._device_boundaries[level]
jets = list(self._level_jets[level]) jets = list(self._level_jets[level])
jets[index] = jets[index].with_amplitude(amplitude) # validates first jets[index] = jets[index].with_amplitude(amplitude) # validates first
velocity = None velocity = None
@@ -735,7 +740,7 @@ class FixedPatchStepper:
self._exchange_valid = False self._exchange_valid = False
for _ in range(n): for _ in range(n):
if self._observe_raw_exchange: # stream-ordered: holds the last interval only if self._observe_raw_exchange: # stream-ordered: holds the last interval only
for b in self._level_boundaries: for b in self._device_boundaries:
if b is not None: if b is not None:
cuda.memset_d8_async(b.interval_integral, 0, cuda.memset_d8_async(b.interval_integral, 0,
len(b.tables["cell"]) * 32, self._stream) len(b.tables["cell"]) * 32, self._stream)
@@ -758,7 +763,7 @@ class FixedPatchStepper:
raise RuntimeError("window exchange observation is disabled") raise RuntimeError("window exchange observation is disabled")
self._window_valid = False self._window_valid = False
try: try:
for b in self._level_boundaries: for b in self._device_boundaries:
if b is not None and b.window_integral is not None: if b is not None and b.window_integral is not None:
cuda.memset_d8(b.window_integral, 0, len(b.tables["cell"]) * 32) cuda.memset_d8(b.window_integral, 0, len(b.tables["cell"]) * 32)
except BaseException: except BaseException:
@@ -771,7 +776,7 @@ class FixedPatchStepper:
if (isinstance(level, bool) or not isinstance(level, Integral) if (isinstance(level, bool) or not isinstance(level, Integral)
or not 0 <= level < len(self._levels)): or not 0 <= level < len(self._levels)):
raise ValueError("level must name an existing level") raise ValueError("level must name an existing level")
b = self._level_boundaries[level] b = self._device_boundaries[level]
if (b is None or isinstance(body_id, bool) or not isinstance(body_id, Integral) if (b is None or isinstance(body_id, bool) or not isinstance(body_id, Integral)
or not 0 <= body_id < len(b.host.objects)): or not 0 <= body_id < len(b.host.objects)):
raise ValueError("body_id must name a body on this level") raise ValueError("body_id must name a body on this level")
@@ -868,7 +873,7 @@ class FixedPatchStepper:
self._stream.synchronize() self._stream.synchronize()
except Exception as exc: except Exception as exc:
error = exc error = exc
for b in getattr(self, "_level_boundaries", ()): for b in getattr(self, "_device_boundaries", ()):
if b is not None: if b is not None:
b.flags = b.device = b.kernel = None b.flags = b.device = b.kernel = None
b.last_raw = b.interval_integral = b.accumulate = None b.last_raw = b.interval_integral = b.accumulate = None
+1 -1
View File
@@ -36,7 +36,7 @@ def blocked_kind(layout, dense_kind):
kind[owner[owned]] = np.asarray(dense_kind).ravel()[owned] kind[owner[owned]] = np.asarray(dense_kind).ravel()[owned]
return kind return kind
def blocked_base_tables(patch, layout, order, wall_aware=False, *, composed_kind=None): def blocked_coarse_tables(patch, layout, order, wall_aware=False, *, composed_kind=None):
"""Parent donor cells and receiver-local arrival slots, in transfer row order. """Parent donor cells and receiver-local arrival slots, in transfer row order.
Stage arrivals only after child reads and parent halo gather, immediately Stage arrivals only after child reads and parent halo gather, immediately
+2 -2
View File
@@ -7,7 +7,7 @@ from CelerisLab.config import LBMConfig
from CelerisLab.lbm.blocks import BlockLayout2D from CelerisLab.lbm.blocks import BlockLayout2D
from CelerisLab.lbm.channel import channel_masks, validate_channel_config, validate_channel_transfers from CelerisLab.lbm.channel import channel_masks, validate_channel_config, validate_channel_transfers
from CelerisLab.lbm.refinement import FixedPatch2D from CelerisLab.lbm.refinement import FixedPatch2D
from CelerisLab.lbm.refinement_tables import blocked_base_tables from CelerisLab.lbm.refinement_tables import blocked_coarse_tables
from .test_numerical_operators import compile_run from .test_numerical_operators import compile_run
@@ -31,7 +31,7 @@ def transfer(bounds, block=(4, 6)):
cfg = LBMConfig(nx=24, ny=20, viscosity=.1) cfg = LBMConfig(nx=24, ny=20, viscosity=.1)
pair = FixedPatch2D.from_config(cfg, bounds) pair = FixedPatch2D.from_config(cfg, bounds)
layout = BlockLayout2D.build(20, 24, *block) layout = BlockLayout2D.build(20, 24, *block)
targets = blocked_base_tables(pair, layout, np.argsort(pair.coalesce_coarse))[1] targets = blocked_coarse_tables(pair, layout, np.argsort(pair.coalesce_coarse))[1]
return pair, layout, targets return pair, layout, targets
+5 -5
View File
@@ -102,7 +102,7 @@ def test_simulation_before_initialize_commits_host_only(monkeypatch):
def stepper_state(run): def stepper_state(run):
boundary = run._level_boundaries[2] boundary = run._device_boundaries[2]
return boundary.tables["surface_velocity"].copy(), boundary.jets, run.jets[2] return boundary.tables["surface_velocity"].copy(), boundary.jets, run.jets[2]
@@ -116,8 +116,8 @@ def test_stepper_invalid_input_writes_nothing_and_stays_ready():
with pytest.raises(error): with pytest.raises(error):
run.set_jet_amplitude(*args, **kw) run.set_jet_amplitude(*args, **kw)
module.cuda.memcpy_htod.assert_not_called() module.cuda.memcpy_htod.assert_not_called()
np.testing.assert_array_equal(run._level_boundaries[2].tables["surface_velocity"], velocity) np.testing.assert_array_equal(run._device_boundaries[2].tables["surface_velocity"], velocity)
assert (run._level_boundaries[2].jets, run.jets[2]) == (jets, level_jets) assert (run._device_boundaries[2].jets, run.jets[2]) == (jets, level_jets)
assert run._status == "ready" assert run._status == "ready"
run.set_jet_amplitude(0, .017, level=2) run.set_jet_amplitude(0, .017, level=2)
assert run.jets[2][0].amplitude == .017 assert run.jets[2][0].amplitude == .017
@@ -136,8 +136,8 @@ def test_stepper_device_failure_keeps_host_and_blocks_reuse(where):
run._stream.synchronize.side_effect = RuntimeError("upload failed") run._stream.synchronize.side_effect = RuntimeError("upload failed")
with pytest.raises(RuntimeError, match="upload failed"): with pytest.raises(RuntimeError, match="upload failed"):
run.set_jet_amplitude(0, .017, level=2) run.set_jet_amplitude(0, .017, level=2)
np.testing.assert_array_equal(run._level_boundaries[2].tables["surface_velocity"], velocity) np.testing.assert_array_equal(run._device_boundaries[2].tables["surface_velocity"], velocity)
assert (run._level_boundaries[2].jets, run.jets[2]) == (jets, level_jets) assert (run._device_boundaries[2].jets, run.jets[2]) == (jets, level_jets)
assert run._status == "failed" assert run._status == "failed"
module.cuda.memcpy_htod.side_effect = None module.cuda.memcpy_htod.side_effect = None
run._stream.synchronize.side_effect = None run._stream.synchronize.side_effect = None
+4 -4
View File
@@ -96,11 +96,11 @@ class MultiplePinnedRemap(unittest.TestCase):
before = run.download() before = run.download()
clone = run.remap_finest((PINS+MOVED+ADDED)[::-1], objects=bodies) clone = run.remap_finest((PINS+MOVED+ADDED)[::-1], objects=bodies)
try: try:
self.assertEqual(tuple(j.amplitude for j in clone._level_boundaries[2].jets), (.017, -.011)) self.assertEqual(tuple(j.amplitude for j in clone._device_boundaries[2].jets), (.017, -.011))
np.testing.assert_array_equal(clone._level_boundaries[2].tables['surface_velocity'], np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
run._level_boundaries[2].tables['surface_velocity']) run._device_boundaries[2].tables['surface_velocity'])
clone.set_jet_amplitude(1, .005, level=2) clone.set_jet_amplitude(1, .005, level=2)
self.assertEqual(run._level_boundaries[2].jets[1].amplitude, -.011) self.assertEqual(run._device_boundaries[2].jets[1].amplitude, -.011)
clone.step() clone.step()
self.assertEqual((run.step_count, clone.step_count), (7, 8)) self.assertEqual((run.step_count, clone.step_count), (7, 8))
for got, want in zip(run.download(), before): for got, want in zip(run.download(), before):
+5 -5
View File
@@ -62,9 +62,9 @@ class PinnedBodyRemap(unittest.TestCase):
run.set_jet_amplitude(0,.017,level=2) run.set_jet_amplitude(0,.017,level=2)
before=run.download() before=run.download()
clone=run.remap_finest(B,objects=objects) clone=run.remap_finest(B,objects=objects)
self.assertEqual(clone._level_boundaries[2].jets[0].amplitude,.017) self.assertEqual(clone._device_boundaries[2].jets[0].amplitude,.017)
np.testing.assert_array_equal(clone._level_boundaries[2].tables['surface_velocity'], np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
run._level_boundaries[2].tables['surface_velocity']) run._device_boundaries[2].tables['surface_velocity'])
self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2])) self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2]))
clone.step() clone.step()
self.assertEqual((run.step_count,clone.step_count),(7,8)) self.assertEqual((run.step_count,clone.step_count),(7,8))
@@ -74,8 +74,8 @@ class PinnedBodyRemap(unittest.TestCase):
def test_static_body_without_actuation(self): def test_static_body_without_actuation(self):
run,module,allocations,storage,objects=ready(jets=()) run,module,allocations,storage,objects=ready(jets=())
clone=run.remap_finest(A,objects=objects) clone=run.remap_finest(A,objects=objects)
self.assertFalse(clone._level_boundaries[2].actuated) self.assertFalse(clone._device_boundaries[2].actuated)
self.assertNotIn('surface_velocity',clone._level_boundaries[2].tables) self.assertNotIn('surface_velocity',clone._device_boundaries[2].tables)
for a,b in zip(run.download(),clone.download()): np.testing.assert_array_equal(a,b) for a,b in zip(run.download(),clone.download()): np.testing.assert_array_equal(a,b)
clone.close();run.close() clone.close();run.close()
+2 -2
View File
@@ -85,7 +85,7 @@ class RawExchange(unittest.TestCase):
def test_four_level_schedule_publication_and_fixed_storage(self): def test_four_level_schedule_publication_and_fixed_storage(self):
run, host, calls, arrays = backend() run, host, calls, arrays = backend()
b = run._level_boundaries[3] b = run._device_boundaries[3]
rows = len(b.tables['cell']) rows = len(b.tables['cell'])
self.assertEqual(b.last_raw.size, rows * 16) self.assertEqual(b.last_raw.size, rows * 16)
self.assertEqual(b.interval_integral.size, rows * 32) self.assertEqual(b.interval_integral.size, rows * 32)
@@ -128,7 +128,7 @@ class RawExchange(unittest.TestCase):
def test_disabled_has_nulls_and_no_observation_allocation_or_launch(self): def test_disabled_has_nulls_and_no_observation_allocation_or_launch(self):
run, host, calls, _ = backend(False) run, host, calls, _ = backend(False)
b = run._level_boundaries[3] b = run._device_boundaries[3]
self.assertIs(b.last_raw, b.interval_integral) self.assertIs(b.last_raw, b.interval_integral)
self.assertIs(b.interval_integral, None) self.assertIs(b.interval_integral, None)
run.step() run.step()
@@ -141,7 +141,7 @@ class Gate(unittest.TestCase):
def test_sparse_parent_tables_preserve_inherited_kinds(self): def test_sparse_parent_tables_preserve_inherited_kinds(self):
child,kind=r.nest_patch_set(CFG,self.parent,CHILDREN) child,kind=r.nest_patch_set(CFG,self.parent,CHILDREN)
layout=self.parent.fine_layout((4,4));order=np.argsort(child.coalesce_coarse) layout=self.parent.fine_layout((4,4));order=np.argsort(child.coalesce_coarse)
donors,targets,bkind=host.blocked_base_tables(child,layout,order,composed_kind=kind) donors,targets,bkind=host.blocked_coarse_tables(child,layout,order,composed_kind=kind)
self.assertTrue((donors < layout.hole).all()) self.assertTrue((donors < layout.hole).all())
self.assertTrue((targets < np.prod(layout.blocked_shape)).all()) self.assertTrue((targets < np.prod(layout.blocked_shape)).all())
np.testing.assert_array_equal(bkind[layout.owner[layout.owned]],kind[layout.owned]) np.testing.assert_array_equal(bkind[layout.owner[layout.owned]],kind[layout.owned])
@@ -158,10 +158,10 @@ class Gate(unittest.TestCase):
hole=int(np.flatnonzero(~layout.owned.ravel())[0]) hole=int(np.flatnonzero(~layout.owned.ravel())[0])
bad=child.ghost_coarse.copy();bad[0]=hole bad=child.ghost_coarse.copy();bad[0]=hole
with self.assertRaisesRegex(ValueError,'owned uncovered'): with self.assertRaisesRegex(ValueError,'owned uncovered'):
host.blocked_base_tables(replace(child,ghost_coarse=bad),layout,order,composed_kind=kind) host.blocked_coarse_tables(replace(child,ghost_coarse=bad),layout,order,composed_kind=kind)
bad=child.coalesce_coarse.copy();bad[0]=hole bad=child.coalesce_coarse.copy();bad[0]=hole
with self.assertRaisesRegex(ValueError,'owned uncovered'): with self.assertRaisesRegex(ValueError,'owned uncovered'):
host.blocked_base_tables(replace(child,coalesce_coarse=bad),layout,order,composed_kind=kind) host.blocked_coarse_tables(replace(child,coalesce_coarse=bad),layout,order,composed_kind=kind)
def bodies_and_jets(speed=.002): def bodies_and_jets(speed=.002):
+4 -4
View File
@@ -67,7 +67,7 @@ class ShiftedJets(unittest.TestCase):
def test_switch_upload_only_then_raw_remap_keeps_live_velocity(self): def test_switch_upload_only_then_raw_remap_keeps_live_velocity(self):
run,module,allocations,storage,objects=ready();clone=None run,module,allocations,storage,objects=ready();clone=None
try: try:
boundary=run._level_boundaries[2] boundary=run._device_boundaries[2]
before=run._download_raw() before=run._download_raw()
original=boundary.tables['surface_velocity'].copy() original=boundary.tables['surface_velocity'].copy()
module.cuda.memcpy_htod.reset_mock() module.cuda.memcpy_htod.reset_mock()
@@ -79,7 +79,7 @@ class ShiftedJets(unittest.TestCase):
for a,b in zip(before.values,run._download_raw().values):self.assertEqual(a.tobytes(),b.tobytes()) for a,b in zip(before.values,run._download_raw().values):self.assertEqual(a.tobytes(),b.tobytes())
clone=run.remap_finest(B,objects=objects) clone=run.remap_finest(B,objects=objects)
self.assertEqual(clone.jets[2][0].amplitude,-jet().amplitude) self.assertEqual(clone.jets[2][0].amplitude,-jet().amplitude)
np.testing.assert_array_equal(clone._level_boundaries[2].tables['surface_velocity'], np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
boundary.tables['surface_velocity']) boundary.tables['surface_velocity'])
self.assertEqual(clone.step_count,7) self.assertEqual(clone.step_count,7)
saved=run._download_raw() saved=run._download_raw()
@@ -92,12 +92,12 @@ class ShiftedJets(unittest.TestCase):
def test_invalid_amplitude_does_not_upload_or_commit(self): def test_invalid_amplitude_does_not_upload_or_commit(self):
run,module,_,_,_=ready() run,module,_,_,_=ready()
try: try:
before=run._level_boundaries[2].tables['surface_velocity'].copy() before=run._device_boundaries[2].tables['surface_velocity'].copy()
module.cuda.memcpy_htod.reset_mock() module.cuda.memcpy_htod.reset_mock()
for amplitude in (float('nan'),float('inf'),1.): for amplitude in (float('nan'),float('inf'),1.):
with self.assertRaises(ValueError):run.set_jet_amplitude(0,amplitude,level=2) with self.assertRaises(ValueError):run.set_jet_amplitude(0,amplitude,level=2)
module.cuda.memcpy_htod.assert_not_called() module.cuda.memcpy_htod.assert_not_called()
np.testing.assert_array_equal(before,run._level_boundaries[2].tables['surface_velocity']) np.testing.assert_array_equal(before,run._device_boundaries[2].tables['surface_velocity'])
finally:run.close() finally:run.close()
+1 -1
View File
@@ -73,7 +73,7 @@ class SingleLevel(unittest.TestCase):
def test_circle_solids_are_not_state(self): def test_circle_solids_are_not_state(self):
run, _, _ = backend((circle(),)) run, _, _ = backend((circle(),))
solid = run._level_boundaries[0].tables["solid"] solid = run._device_boundaries[0].tables["solid"]
self.assertTrue(solid.any()) self.assertTrue(solid.any())
np.testing.assert_array_equal(run._state[0], ~solid) np.testing.assert_array_equal(run._state[0], ~solid)
+1 -1
View File
@@ -37,7 +37,7 @@ class UnionBodySupport(unittest.TestCase):
run.set_jet_amplitude(0,.017,level=2) run.set_jet_amplitude(0,.017,level=2)
clone=run.remap_finest(((16,16,24,40),(24,16,32,40),A[1]),objects=objects) clone=run.remap_finest(((16,16,24,40),(24,16,32,40),A[1]),objects=objects)
try: try:
self.assertEqual(clone._level_boundaries[2].jets[0].amplitude,.017) self.assertEqual(clone._device_boundaries[2].jets[0].amplitude,.017)
self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2])) self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2]))
for a,b in zip(run.download(),clone.download()):np.testing.assert_array_equal(a,b) for a,b in zip(run.download(),clone.download()):np.testing.assert_array_equal(a,b)
finally:clone.close() finally:clone.close()
@@ -240,7 +240,7 @@ class SingleLevelGPU(unittest.TestCase):
with CudaContext(0) as context: with CudaContext(0) as context:
for split in SPLITS: for split in SPLITS:
with single(context=context, block=split, objects=[circle()]) as backend: with single(context=context, block=split, objects=[circle()]) as backend:
tables = backend._level_boundaries[0].tables tables = backend._device_boundaries[0].tables
for key, values, dtype in (("q", [l.q for l in links], np.float32), for key, values, dtype in (("q", [l.q for l in links], np.float32),
("dir", [l.dir for l in links], np.uint8), ("dir", [l.dir for l in links], np.uint8),
("fallback", [l.fallback for l in links], np.uint8)): ("fallback", [l.fallback for l in links], np.uint8)):
@@ -59,7 +59,7 @@ class ExchangeWindowCPU(unittest.TestCase):
stack.enter_context(patch.object(cuda, "memcpy_dtoh", stack.enter_context(patch.object(cuda, "memcpy_dtoh",
side_effect=lambda out, src: out.fill(0))) side_effect=lambda out, src: out.fill(0)))
backend = build(m) backend = build(m)
b = backend._level_boundaries[1] b = backend._device_boundaries[1]
self.assertIsNotNone(b.window_integral) self.assertIsNotNone(b.window_integral)
with self.assertRaises(RuntimeError): with self.assertRaises(RuntimeError):
backend.start_exchange_window() # not initialized backend.start_exchange_window() # not initialized
@@ -96,7 +96,7 @@ class ExchangeWindowCPU(unittest.TestCase):
with mocked_backend() as m, ExitStack() as stack: with mocked_backend() as m, ExitStack() as stack:
stack.enter_context(patch.object(cuda, "memset_d8_async")) stack.enter_context(patch.object(cuda, "memset_d8_async"))
backend = build(m, window=False) backend = build(m, window=False)
self.assertIsNone(backend._level_boundaries[1].window_integral) self.assertIsNone(backend._device_boundaries[1].window_integral)
backend.initialize(*[np.full(l.shape, 1 / 9, np.float32) for l in backend._levels]) backend.initialize(*[np.full(l.shape, 1 / 9, np.float32) for l in backend._levels])
backend.step(1) backend.step(1)
self.assertEqual(len(accumulations(m)), 2) self.assertEqual(len(accumulations(m)), 2)
@@ -120,7 +120,7 @@ class ExchangeWindowGPU(unittest.TestCase):
windowed.step(2) windowed.step(2)
plain.step(2) plain.step(2)
windowed.start_exchange_window() windowed.start_exchange_window()
b = windowed._level_boundaries[1] b = windowed._device_boundaries[1]
count = len(b.tables["cell"]) count = len(b.tables["cell"])
factors = np.array([.25, .25, .125, .25]) # h*dt, h*h*dt, h*h for level-1 h=.5 factors = np.array([.25, .25, .125, .25]) # h*dt, h*h*dt, h*h for level-1 h=.5
scale = np.zeros(4) scale = np.zeros(4)
@@ -212,7 +212,7 @@ class ExchangeWindowShiftedMRTGPU(unittest.TestCase):
run.step(2) run.step(2)
if run is not plain: if run is not plain:
run.start_exchange_window() run.start_exchange_window()
b = windowed._level_boundaries[1] b = windowed._device_boundaries[1]
count = len(b.tables["cell"]) count = len(b.tables["cell"])
record = dict(scale=np.zeros(4), substeps=0) record = dict(scale=np.zeros(4), substeps=0)
audit_substeps(windowed, b, count, record) audit_substeps(windowed, b, count, record)
@@ -228,7 +228,7 @@ class ExchangeWindowShiftedMRTGPU(unittest.TestCase):
same_bits(self, a, e) same_bits(self, a, e)
rows = [np.empty((count, 4), np.float64) for _ in range(2)] rows = [np.empty((count, 4), np.float64) for _ in range(2)]
cuda.memcpy_dtoh(rows[0], b.window_integral) cuda.memcpy_dtoh(rows[0], b.window_integral)
cuda.memcpy_dtoh(rows[1], split._level_boundaries[1].window_integral) cuda.memcpy_dtoh(rows[1], split._device_boundaries[1].window_integral)
np.testing.assert_array_equal(rows[0].view(np.uint64), rows[1].view(np.uint64)) np.testing.assert_array_equal(rows[0].view(np.uint64), rows[1].view(np.uint64))
out = windowed.read_exchange_window(0, level=1) out = windowed.read_exchange_window(0, level=1)
@@ -57,7 +57,7 @@ def make(context):
def topology(run, rectangles, objects): def topology(run, rectangles, objects):
return RemapTopology.from_config( return RemapTopology.from_config(
CFG, run.hierarchy.pairs[0], rectangles, child_objects=objects, CFG, run.hierarchy.pairs[0], rectangles, child_objects=objects,
child_jets=run._level_boundaries[2].jets, child_jets=run._device_boundaries[2].jets,
) )
@@ -74,7 +74,7 @@ class MultiplePinnedRemapGPU(unittest.TestCase):
def change_amplitudes(self, run, amplitudes): def change_amplitudes(self, run, amplitudes):
# Each jet selects a different body. Updating one must preserve the # Each jet selects a different body. Updating one must preserve the
# other body's current per-link velocities, including its driven links. # other body's current per-link velocities, including its driven links.
boundary = run._level_boundaries[2] boundary = run._device_boundaries[2]
body_ids = run.boundaries[2].curved_host[3] body_ids = run.boundaries[2].curved_host[3]
for index, amplitude in enumerate(amplitudes): for index, amplitude in enumerate(amplitudes):
before = boundary.tables['surface_velocity'].copy() before = boundary.tables['surface_velocity'].copy()
@@ -91,7 +91,7 @@ class MultiplePinnedRemapGPU(unittest.TestCase):
self.assertTrue(same_physical_boundary(new.boundaries[2], mesh.boundary)) self.assertTrue(same_physical_boundary(new.boundaries[2], mesh.boundary))
# Compare physical cut-link rows, never layout-dependent blocked addresses. # Compare physical cut-link rows, never layout-dependent blocked addresses.
for run in (old, new): for run in (old, new):
boundary = run._level_boundaries[2] boundary = run._device_boundaries[2]
self.assertEqual(tuple(jet.amplitude for jet in boundary.jets), amplitudes) self.assertEqual(tuple(jet.amplitude for jet in boundary.jets), amplitudes)
exact(self, boundary.tables['surface_velocity'], mesh.surface_velocity, exact(self, boundary.tables['surface_velocity'], mesh.surface_velocity,
'latest physical per-link velocities') 'latest physical per-link velocities')
@@ -39,9 +39,9 @@ class PinnedBodyRemapGPU(unittest.TestCase):
exact(self,getattr(after.boundaries[2],name),getattr(before.boundaries[2],name),name) exact(self,getattr(after.boundaries[2],name),getattr(before.boundaries[2],name),name)
for a,b in zip(after.boundaries[2].curved_host,before.boundaries[2].curved_host): for a,b in zip(after.boundaries[2].curved_host,before.boundaries[2].curved_host):
exact(self,a,b,'physical cut links') exact(self,a,b,'physical cut links')
exact(self,after._level_boundaries[2].tables['surface_velocity'], exact(self,after._device_boundaries[2].tables['surface_velocity'],
before._level_boundaries[2].tables['surface_velocity'],'current per-link velocities') before._device_boundaries[2].tables['surface_velocity'],'current per-link velocities')
self.assertEqual(after._level_boundaries[2].jets[0].amplitude,.017) self.assertEqual(after._device_boundaries[2].jets[0].amplitude,.017)
def test_noop_next_step_bitwise(self): def test_noop_next_step_bitwise(self):
with CudaContext(0) as ctx, ExitStack() as owned: with CudaContext(0) as ctx, ExitStack() as owned:
@@ -63,7 +63,7 @@ class PinnedBodyRemapGPU(unittest.TestCase):
current.set_jet_amplitude(0,.017,level=2) current.set_jet_amplitude(0,.017,level=2)
for rects in (B,A): for rects in (B,A):
previous=current;before=previous.download();count=previous.step_count previous=current;before=previous.download();count=previous.step_count
jets=previous._level_boundaries[2].jets jets=previous._device_boundaries[2].jets
old=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],previous.hierarchy.pairs[1].bounds, old=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],previous.hierarchy.pairs[1].bounds,
child_objects=objects,child_jets=jets) child_objects=objects,child_jets=jets)
new=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],rects, new=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],rects,
@@ -105,7 +105,7 @@ class RawExchangeGPU(unittest.TestCase):
arrays = [state(l.shape) for l in measured._levels] arrays = [state(l.shape) for l in measured._levels]
measured.initialize(*arrays) measured.initialize(*arrays)
plain.initialize(*arrays) plain.initialize(*arrays)
b = measured._level_boundaries[1] b = measured._device_boundaries[1]
count = len(b.tables['cell']) count = len(b.tables['cell'])
kernel = measured._levels[1].kernel kernel = measured._levels[1].kernel
samples = [] samples = []
@@ -158,7 +158,7 @@ class RawExchangeGPU(unittest.TestCase):
FixedPatchStepper(cfg, context=ctx, **settings) as plain: FixedPatchStepper(cfg, context=ctx, **settings) as plain:
arrays = [state(l.shape, ux=0.) for l in run._levels] arrays = [state(l.shape, ux=0.) for l in run._levels]
run.initialize(*arrays); plain.initialize(*arrays) run.initialize(*arrays); plain.initialize(*arrays)
b = run._level_boundaries[3] b = run._device_boundaries[3]
count = len(b.tables['cell']) count = len(b.tables['cell'])
kernel = run._levels[3].kernel kernel = run._levels[3].kernel
samples, independent_samples = [], [] samples, independent_samples = [], []
@@ -213,7 +213,7 @@ class RawExchangeGPU(unittest.TestCase):
FixedPatchStepper(CFG, context=ctx, **fine((16, 20)), observe_raw_exchange=True) as amr: FixedPatchStepper(CFG, context=ctx, **fine((16, 20)), observe_raw_exchange=True) as amr:
amr.initialize(*(state(l.shape) for l in amr._levels)) amr.initialize(*(state(l.shape) for l in amr._levels))
uniform.initialize(state((9, 48, 64))) uniform.initialize(state((9, 48, 64)))
b = amr._level_boundaries[1] b = amr._device_boundaries[1]
count = len(b.tables['cell']) count = len(b.tables['cell'])
links = obj.geometry.build_cut_links(64, 48) links = obj.geometry.build_cut_links(64, 48)
np.testing.assert_array_equal([v.dir for v in links], b.tables['dir']) np.testing.assert_array_equal([v.dir for v in links], b.tables['dir'])
@@ -221,7 +221,7 @@ class RawExchangeGPU(unittest.TestCase):
# Same fine-grid force arms after translation to the same center: the # Same fine-grid force arms after translation to the same center: the
# single level's static BFL gets the AMR Jet column and force arms; the raw # single level's static BFL gets the AMR Jet column and force arms; the raw
# output pointer is substituted at launch (no stepper-side accumulation). # output pointer is substituted at launch (no stepper-side accumulation).
root = uniform._level_boundaries[0] root = uniform._device_boundaries[0]
velocity = uniform._upload(b.tables['surface_velocity']) velocity = uniform._upload(b.tables['surface_velocity'])
root.device.update(surface_velocity=velocity, rx=uniform._upload(b.tables['rx']), root.device.update(surface_velocity=velocity, rx=uniform._upload(b.tables['rx']),
ry=uniform._upload(b.tables['ry'])) ry=uniform._upload(b.tables['ry']))
@@ -24,7 +24,8 @@ from CelerisLab.lbm import refinement_stepper
from CelerisLab.lbm.blocks import BlockLayout2D from CelerisLab.lbm.blocks import BlockLayout2D
from CelerisLab.lbm.descriptors import D2Q9_EX, D2Q9_EY from CelerisLab.lbm.descriptors import D2Q9_EX, D2Q9_EY
from CelerisLab.lbm.refinement import ACTIVE, OUTER_GHOST, FixedPatch2D from CelerisLab.lbm.refinement import ACTIVE, OUTER_GHOST, FixedPatch2D
from CelerisLab.lbm.refinement_stepper import FixedPatchStepper, blocked_base_tables from CelerisLab.lbm.refinement_stepper import FixedPatchStepper
from CelerisLab.lbm.refinement_tables import blocked_coarse_tables
from test_refinement_blocked_fine import NX, NY, BOUNDS, config, state, circle, actuation from test_refinement_blocked_fine import NX, NY, BOUNDS, config, state, circle, actuation
BASE_BLOCKS = ((24, 32), (12, 8), (3, 4), (1, 1)) BASE_BLOCKS = ((24, 32), (12, 8), (3, 4), (1, 1))
@@ -108,7 +109,7 @@ class BlockedBaseCPU(unittest.TestCase):
leaves[[0, -1]] = False leaves[[0, -1]] = False
for by, bx in BASE_BLOCKS: for by, bx in BASE_BLOCKS:
layout = BlockLayout2D.build(NY, NX, by, bx, periodic=(True, mode == "periodic")) layout = BlockLayout2D.build(NY, NX, by, bx, periodic=(True, mode == "periodic"))
donors, targets, kind = blocked_base_tables(p, layout, order, mode != "periodic") donors, targets, kind = blocked_coarse_tables(p, layout, order, mode != "periodic")
ly, lx = layout.local_shape ly, lx = layout.local_shape
L = ly * lx L = ly * lx
expected, source_cells = [], [] expected, source_cells = [], []
@@ -138,7 +139,7 @@ class BlockedBaseCPU(unittest.TestCase):
p = FixedPatch2D.from_config(config("periodic"), BOUNDS) p = FixedPatch2D.from_config(config("periodic"), BOUNDS)
order = np.argsort(p.coalesce_coarse) order = np.argsort(p.coalesce_coarse)
layout = BlockLayout2D.build(NY, NX, 3, 4) layout = BlockLayout2D.build(NY, NX, 3, 4)
donors, targets, _ = blocked_base_tables(p, layout, order) donors, targets, _ = blocked_coarse_tables(p, layout, order)
dense = np.arange(9 * NY * NX, dtype=np.float32).reshape(9, NY, NX) dense = np.arange(9 * NY * NX, dtype=np.float32).reshape(9, NY, NX)
initial = layout.to_blocked(dense) initial = layout.to_blocked(dense)
arrivals = np.arange(len(targets), dtype=np.float32) + 10000 arrivals = np.arange(len(targets), dtype=np.float32) + 10000
@@ -178,7 +179,7 @@ class BlockedBaseCPU(unittest.TestCase):
# With one parent block, all eight grazing sources are owner cores. # With one parent block, all eight grazing sources are owner cores.
# Each is read for four fine ghost children: late explosion corrupts32. # Each is read for four fine ghost children: late explosion corrupts32.
single = BlockLayout2D.build(NY, NX, NY, NX) single = BlockLayout2D.build(NY, NX, NY, NX)
donors, targets, _ = blocked_base_tables(p, single, order) donors, targets, _ = blocked_coarse_tables(p, single, order)
initial = single.to_blocked(dense) initial = single.to_blocked(dense)
staged = initial.copy() staged = initial.copy()
staged.ravel()[targets.astype(np.intp)] = arrivals staged.ravel()[targets.astype(np.intp)] = arrivals
@@ -21,7 +21,8 @@ from CelerisLab.lbm.descriptors import D2Q9_EX, D2Q9_EY
from CelerisLab.lbm.descriptors import SOLID from CelerisLab.lbm.descriptors import SOLID
from CelerisLab.lbm.refinement import (ACTIVE, INNER_GHOST, OUTER_GHOST, FixedPatch2D, from CelerisLab.lbm.refinement import (ACTIVE, INNER_GHOST, OUTER_GHOST, FixedPatch2D,
prepare_fine_boundary) prepare_fine_boundary)
from CelerisLab.lbm.refinement_stepper import FixedPatchStepper, blocked_fine_tables from CelerisLab.lbm.refinement_stepper import FixedPatchStepper
from CelerisLab.lbm.refinement_tables import blocked_fine_tables
NX, NY, BOUNDS, STEPS = 32, 24, (8, 6, 24, 18), 8 NX, NY, BOUNDS, STEPS = 32, 24, (8, 6, 24, 18), 8
FINE_BLOCKS = ((32, 40), (16, 20), (16, 10)) # fine window 32x40: 1x1, 2x2, 4x2 FINE_BLOCKS = ((32, 40), (16, 20), (16, 10)) # fine window 32x40: 1x1, 2x2, 4x2
@@ -186,8 +187,8 @@ class BlockedFineCPU(unittest.TestCase):
backend, calls = self.mocked_backend([circle(spec)]) backend, calls = self.mocked_backend([circle(spec)])
c, f = backend._levels c, f = backend._levels
coarse, cur, spare = c.current, f.current, f.spare coarse, cur, spare = c.current, f.current, f.spare
flags, kind = backend._level_boundaries[1].flags, backend._levels[1].kind flags, kind = backend._device_boundaries[1].flags, backend._levels[1].kind
velocity = backend._level_boundaries[1].device.get("surface_velocity") velocity = backend._device_boundaries[1].device.get("surface_velocity")
self.assertEqual(velocity is None, spec is None) self.assertEqual(velocity is None, spec is None)
backend.step(1) backend.step(1)
backend.close() backend.close()
@@ -212,8 +213,8 @@ class BlockedFineCPU(unittest.TestCase):
def test_actuation_without_velocity_buffer_fails_before_bfl(self): def test_actuation_without_velocity_buffer_fails_before_bfl(self):
for speed in (0.0, .002): # any spec, even zero amplitude for speed in (0.0, .002): # any spec, even zero amplitude
backend, calls = self.mocked_backend([circle(actuation(speed))]) backend, calls = self.mocked_backend([circle(actuation(speed))])
self.assertIsNotNone(backend._level_boundaries[1].tables["surface_velocity"]) self.assertIsNotNone(backend._device_boundaries[1].tables["surface_velocity"])
del backend._level_boundaries[1].device["surface_velocity"] del backend._device_boundaries[1].device["surface_velocity"]
with self.assertRaisesRegex(RuntimeError, "Surface velocity device buffer is missing"): with self.assertRaisesRegex(RuntimeError, "Surface velocity device buffer is missing"):
backend.step(1) backend.step(1)
names = [name for name, _ in calls] names = [name for name, _ in calls]
@@ -62,8 +62,8 @@ class JetStepperCPU(unittest.TestCase):
C = (16.0, 12.0) C = (16.0, 12.0)
whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.) whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.)
backend, _ = mocked([circle()], [whole], FINE_BLOCKS[1]) backend, _ = mocked([circle()], [whole], FINE_BLOCKS[1])
_, jet_id = backend._level_boundaries[1].jet_cache _, jet_id = backend._device_boundaries[1].jet_cache
self.assertGreater(len(np.unique(backend._level_boundaries[1].tables["block"][jet_id == 0])), 1) self.assertGreater(len(np.unique(backend._device_boundaries[1].tables["block"][jet_id == 0])), 1)
backend.close() backend.close()
def test_fine_ownership_uses_remapped_physical_owner(self): def test_fine_ownership_uses_remapped_physical_owner(self):
@@ -76,10 +76,10 @@ class JetStepperCPU(unittest.TestCase):
source = SimObject(99, geometry=CircleGeometry(*C, 2.15), center=C, radius=2.15) source = SimObject(99, geometry=CircleGeometry(*C, 2.15), center=C, radius=2.15)
whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.) whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.)
backend, _ = mocked([source], [whole], FINE_BLOCKS[1]) backend, _ = mocked([source], [whole], FINE_BLOCKS[1])
_, jet_id = backend._level_boundaries[1].jet_cache _, jet_id = backend._device_boundaries[1].jet_cache
owners = backend.boundaries[1].curved_host[3][jet_id == 0] owners = backend.boundaries[1].curved_host[3][jet_id == 0]
self.assertEqual(np.unique(owners).tolist(), [0]) self.assertEqual(np.unique(owners).tolist(), [0])
self.assertGreater(len(np.unique(backend._level_boundaries[1].tables["block"][jet_id == 0])), 1) self.assertGreater(len(np.unique(backend._device_boundaries[1].tables["block"][jet_id == 0])), 1)
self.assertEqual(source.obj_id, 99) self.assertEqual(source.obj_id, 99)
backend.close() backend.close()
cfg = config("periodic") cfg = config("periodic")
@@ -97,8 +97,8 @@ class JetStepperCPU(unittest.TestCase):
def test_blocked_jet_velocity_matches_legacy_rows(self): def test_blocked_jet_velocity_matches_legacy_rows(self):
legacy, _ = mocked([circle(actuation(SPEED))], (), FINE_BLOCKS[1]) legacy, _ = mocked([circle(actuation(SPEED))], (), FINE_BLOCKS[1])
driven, _ = mocked([circle()], [jet()], FINE_BLOCKS[1]) driven, _ = mocked([circle()], [jet()], FINE_BLOCKS[1])
a = legacy._level_boundaries[1].tables["surface_velocity"] a = legacy._device_boundaries[1].tables["surface_velocity"]
b = driven._level_boundaries[1].tables["surface_velocity"] b = driven._device_boundaries[1].tables["surface_velocity"]
diff = np.abs(a - b).max() diff = np.abs(a - b).max()
print(f"fine jet vs legacy velocity rows: max |dU| = {diff:.3e} (exact={diff == 0})") print(f"fine jet vs legacy velocity rows: max |dU| = {diff:.3e} (exact={diff == 0})")
self.assertLessEqual(diff, 4 * np.finfo(np.float32).eps * SPEED) self.assertLessEqual(diff, 4 * np.finfo(np.float32).eps * SPEED)
@@ -106,8 +106,8 @@ class JetStepperCPU(unittest.TestCase):
def test_blocked_amplitude_update_rewrites_only_velocity_column(self): def test_blocked_amplitude_update_rewrites_only_velocity_column(self):
backend, htod = mocked([circle()], [jet()], FINE_BLOCKS[1]) backend, htod = mocked([circle()], [jet()], FINE_BLOCKS[1])
before = backend._level_boundaries[1].tables["surface_velocity"].copy() before = backend._device_boundaries[1].tables["surface_velocity"].copy()
device = backend._level_boundaries[1].device["surface_velocity"] device = backend._device_boundaries[1].device["surface_velocity"]
with patch.object(cuda, "memcpy_htod", htod): with patch.object(cuda, "memcpy_htod", htod):
htod.reset_mock() htod.reset_mock()
backend.set_jet_amplitude(0, -SPEED, level=1) backend.set_jet_amplitude(0, -SPEED, level=1)
@@ -115,7 +115,7 @@ class JetStepperCPU(unittest.TestCase):
self.assertEqual(htod.call_count, 1) self.assertEqual(htod.call_count, 1)
self.assertIs(ptr, device) self.assertIs(ptr, device)
np.testing.assert_allclose(sent, -before, rtol=1e-6, atol=0) np.testing.assert_allclose(sent, -before, rtol=1e-6, atol=0)
np.testing.assert_array_equal(backend._level_boundaries[1].tables["surface_velocity"], sent) np.testing.assert_array_equal(backend._device_boundaries[1].tables["surface_velocity"], sent)
with self.assertRaises(ValueError): with self.assertRaises(ValueError):
backend.set_jet_amplitude(0, 1.5, level=1) backend.set_jet_amplitude(0, 1.5, level=1)
self.assertEqual(htod.call_count, 1) self.assertEqual(htod.call_count, 1)
@@ -133,7 +133,7 @@ class JetStepperCPU(unittest.TestCase):
with self.assertRaises(IndexError): with self.assertRaises(IndexError):
backend.set_jet_amplitude(1, -SPEED, level=1) backend.set_jet_amplitude(1, -SPEED, level=1)
self.assertEqual(backend._status, status) self.assertEqual(backend._status, status)
before = backend._level_boundaries[1].tables["surface_velocity"].copy() before = backend._device_boundaries[1].tables["surface_velocity"].copy()
boom = Mock(side_effect=RuntimeError("copy failed")) boom = Mock(side_effect=RuntimeError("copy failed"))
with patch.object(cuda, "memcpy_htod", boom), patch.object(cuda, "memcpy_htod_async", boom), \ with patch.object(cuda, "memcpy_htod", boom), patch.object(cuda, "memcpy_htod_async", boom), \
patch.object(cuda, "pagelocked_empty", side_effect=lambda shape, dtype: np.empty(shape, dtype)): patch.object(cuda, "pagelocked_empty", side_effect=lambda shape, dtype: np.empty(shape, dtype)):
@@ -141,7 +141,7 @@ class JetStepperCPU(unittest.TestCase):
backend.set_jet_amplitude(0, -SPEED, level=1) backend.set_jet_amplitude(0, -SPEED, level=1)
self.assertEqual(backend._status, "failed", fine_block) self.assertEqual(backend._status, "failed", fine_block)
self.assertEqual(backend.jets[1][0].amplitude, SPEED) self.assertEqual(backend.jets[1][0].amplitude, SPEED)
np.testing.assert_array_equal(backend._level_boundaries[1].tables["surface_velocity"], before) np.testing.assert_array_equal(backend._device_boundaries[1].tables["surface_velocity"], before)
with self.assertRaisesRegex(RuntimeError, "failed"): with self.assertRaisesRegex(RuntimeError, "failed"):
backend.step(1) backend.step(1)
with self.assertRaisesRegex(RuntimeError, "failed"): with self.assertRaisesRegex(RuntimeError, "failed"):
@@ -102,7 +102,7 @@ class Ledger:
def step(self): def step(self):
run=self.run;self.samples.clear();run.step() run=self.run;self.samples.clear();run.step()
b=run._level_boundaries[self.index];h=run._exchange_spacing[self.index] b=run._device_boundaries[self.index];h=run._exchange_spacing[self.index]
self.test.assertEqual(len(self.samples),2**self.index) self.test.assertEqual(len(self.samples),2**self.index)
expected=np.zeros_like(self.samples[0],dtype=np.float64) expected=np.zeros_like(self.samples[0],dtype=np.float64)
weights=np.array([h*h,h*h,h*h*h,h*h]) weights=np.array([h*h,h*h,h*h*h,h*h])
@@ -140,7 +140,7 @@ class ShiftedJetsGPU(unittest.TestCase):
level_objects=[(),(obj,)],level_jets=[(),(whole,)],observe_raw_exchange=True) as run: level_objects=[(),(obj,)],level_jets=[(),(whole,)],observe_raw_exchange=True) as run:
run._initialize_raw(initial(run,rho,u=.013,v=0.)) run._initialize_raw(initial(run,rho,u=.013,v=0.))
expected=run._download_raw();ledger=Ledger(self,run,1) expected=run._download_raw();ledger=Ledger(self,run,1)
self.assertTrue(run._level_boundaries[1].tables['aperture'].all()) self.assertTrue(run._device_boundaries[1].tables['aperture'].all())
for _ in range(4):ledger.step() for _ in range(4):ledger.step()
error=max(float(abs(a[:,leaf]-b[:,leaf]).max()) for a,b,leaf in error=max(float(abs(a[:,leaf]-b[:,leaf]).max()) for a,b,leaf in
zip(run._download_raw().values,expected.values,run.leaf)) zip(run._download_raw().values,expected.values,run.leaf))
@@ -201,7 +201,7 @@ class ShiftedJetsGPU(unittest.TestCase):
self.exact_raw(before,old._download_raw()) self.exact_raw(before,old._download_raw())
self.exact_raw(_StoredPopulations(expected.g,True),new._download_raw()) self.exact_raw(_StoredPopulations(expected.g,True),new._download_raw())
self.assertEqual(new.jets[2][0].amplitude,-.001) self.assertEqual(new.jets[2][0].amplitude,-.001)
np.testing.assert_array_equal(new._level_boundaries[2].tables['surface_velocity'],old._level_boundaries[2].tables['surface_velocity']) np.testing.assert_array_equal(new._device_boundaries[2].tables['surface_velocity'],old._device_boundaries[2].tables['surface_velocity'])
newledger=Ledger(self,new,2);ledger.step();newledger.step() newledger=Ledger(self,new,2);ledger.step();newledger.step()
self.assertEqual((old.step_count,new.step_count),(count+1,count+1)) self.assertEqual((old.step_count,new.step_count),(count+1,count+1))
ledger.report('L2 source');newledger.report('L2 remapped') ledger.report('L2 source');newledger.report('L2 remapped')
+15 -6
View File
@@ -81,7 +81,7 @@ def public_jets(center, radius, scenario):
) )
def jet_interval(handle, scenario, begin, n, diameter, *, level=None): def jet_interval(handle, scenario, begin, n, diameter, level=None):
"""Apply a reversal only at a completed interval boundary, before advancing.""" """Apply a reversal only at a completed interval boundary, before advancing."""
if scenario != "none": if scenario != "none":
boundaries = [round(t * diameter / SPEED) for t in (250, 260, 270, 280, 290, 300, 305)] boundaries = [round(t * diameter / SPEED) for t in (250, 260, 270, 280, 290, 300, 305)]
@@ -92,9 +92,10 @@ def jet_interval(handle, scenario, begin, n, diameter, *, level=None):
if begin < switch: if begin < switch:
return min(n, switch - begin), None return min(n, switch - begin), None
if begin == switch: if begin == switch:
target = {} if level is None else dict(level=level) # FixedPatchStepper needs an explicit Jet level; Simulation has one level.
handle.set_jet_amplitude(0, -JET_PEAK, **target) kw = {} if level is None else {"level": level}
handle.set_jet_amplitude(1, JET_PEAK, **target) handle.set_jet_amplitude(0, -JET_PEAK, **kw)
handle.set_jet_amplitude(1, JET_PEAK, **kw)
return n, dict(step=begin, tstar=300.0, amplitudes=[-JET_PEAK, JET_PEAK]) return n, dict(step=begin, tstar=300.0, amplitudes=[-JET_PEAK, JET_PEAK])
return n, None return n, None
@@ -257,7 +258,7 @@ def host_geometry(source, h, collision="SRT", jet_scenario="none", placement="wi
validate_channel_config, validate_channel_config,
validate_channel_transfers, validate_channel_transfers,
) )
from CelerisLab.lbm.refinement_tables import blocked_base_tables from CelerisLab.lbm import refinement_tables
from CelerisLab.body.geometry.circle import CircleGeometry from CelerisLab.body.geometry.circle import CircleGeometry
from CelerisLab.body.objects import SimObject from CelerisLab.body.objects import SimObject
@@ -269,7 +270,15 @@ def host_geometry(source, h, collision="SRT", jet_scenario="none", placement="wi
pair = FixedPatch2D.from_config(cfg, PLACEMENTS[placement]) pair = FixedPatch2D.from_config(cfg, PLACEMENTS[placement])
layout = BlockLayout2D.build(cfg.ny, cfg.nx, *BASE_BLOCK) layout = BlockLayout2D.build(cfg.ny, cfg.nx, *BASE_BLOCK)
validate_channel_config(cfg, [BASE_BLOCK, FINE_BLOCK]) validate_channel_config(cfg, [BASE_BLOCK, FINE_BLOCK])
targets = blocked_base_tables(pair, layout, np.argsort(pair.coalesce_coarse))[1] # Select by the imported source's own API, not by scenario: the pinned historical
# source (verify_source) defines blocked_base_tables, the current checkout (host
# tests, public Jet scenarios) blocked_coarse_tables. Exactly one must exist.
names = [n for n in ("blocked_base_tables", "blocked_coarse_tables")
if hasattr(refinement_tables, n)]
if len(names) != 1:
raise RuntimeError(f"Unexpected coarse transfer table API: {names}")
coarse_tables = getattr(refinement_tables, names[0])
targets = coarse_tables(pair, layout, np.argsort(pair.coalesce_coarse))[1]
validate_channel_transfers(pair, layout, targets) validate_channel_transfers(pair, layout, targets)
BlockLayout2D.build(*pair.fine_shape, *FINE_BLOCK, periodic=(False, False)) BlockLayout2D.build(*pair.fine_shape, *FINE_BLOCK, periodic=(False, False))
body = SimObject( body = SimObject(