Compare commits

..
Author SHA1 Message Date
ClaudeandClaude Opus 5.5 18f018befe README: several regions are refinements[0], not bounds
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 02:02:38 +00:00
ClaudeandClaude Opus 5.5 ec9b3d721c P1 (12/12): one FixedPatchStepper construction and initialize API
Removes the legacy entry points now that every caller is migrated:
FixedPatchStepper(cfg, *, context, refinements, coarse_boundary,
level_blocks, level_objects, level_jets, observe_*, explosion) and
initialize(*level_g). Deleted: positional bounds, fine_block, base_block,
fine_objects, jets, child, child_block, child_objects, child_jets,
child_rectangles_in_parent, initialize(coarse_g/fine_g/child_g),
_legacy_levels and _initial_levels. Docstring, README and docs describe
the level lists. Compatibility views (patch, fine_boundary, ...) and the
PeriodicPatchStepper alias are unchanged (separate follow-up).

Error -> rule table (legacy keyword check -> level-list rule now firing):
- several patches without fine_block/base_block -> _check_levels
  'several patches require fine_block and base_block' (unchanged rule)
- base_block without fine_block -> 'a blocked level requires blocked
  refined levels below it'
- nested child without all blocks -> 'three or more levels must all be
  blocked' / 'a blocked level requires ...'
- fine_objects/jets with a nested child ('only on L2') -> 'bodies/Jets are
  supported only on the finest level of a nested hierarchy'
- child_objects/child_jets without a nested child -> 'level_objects/
  level_jets needs one entry per level'
- scalar child with several L1 patches ('scalar child requires one parent
  patch') -> hierarchy 'a scalar child requires a scalar parent'
- child vs child_rectangles_in_parent mutually exclusive -> no level-list
  equivalent (cannot be expressed); covered by the new explicit
  RemovedLegacyAPI test: removed keywords and positional bounds raise
  TypeError before any context use, initialize has only *level_g.
The InitialLevels test of the removed adapter is replaced by that test.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 02:02:29 +00:00
ClaudeandClaude Opus 5.5 84b76f49cc P1 (11/n): multipatch GPU/CPU tests use the level-list constructor
One local levels(bounds, blocks, objects, jets, **rest) builds the
two-level lists; run_steps and every rejection case use it. All five
before-CUDA rejections, the margin and one-jet checks keep their inputs.
Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 02:00:59 +00:00
ClaudeandClaude Opus 5.5 050d12b76c P1 (10/n): blocked-fine GPU/CPU tests use the level-list constructor
Dense root + blocked fine level is level_blocks=[None, block]; bodies go
to level_objects[1]. Both before-CUDA rejections (non-dividing and
malformed fine block) keep their inputs. Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 02:00:15 +00:00
ClaudeandClaude Opus 5.5 ec3e3b928a P1 (9/n): refinement follow-up tests pass refinements=[bounds]
FixedPatchStepper and the PeriodicPatchStepper alias are called with
refinements=[bounds] instead of the positional legacy bounds (same
_legacy_levels result). Assertions, including the before-CUDA and
four-cell/omega/layout rejections, unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:59:56 +00:00
ClaudeandClaude Opus 5.5 2be80324ad P1 (8/n): three-level GPU/CPU tests use the level-list constructor
KW becomes levels() (scalar L1 + scalar L2, per-level blocks, objects,
Jets). All rejection cases are kept per level. Two cases change branch:
L2 bodies/Jets with no nested child (legacy 'child_objects/child_jets
require a nested child') are now an entry for a level that does not
exist ('needs one entry per level'); several L1 patches + scalar child
now meets 'a scalar child requires a scalar parent'. Both still reject
before CUDA. Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:59:37 +00:00
ClaudeandClaude Opus 5.5 0fcd4ac4cb P1 (7/n): MRT, shifted-jet and blocked-base GPU tests migrated
test_refinement_mrt_gpu, test_shifted_jets_gpu and
test_refinement_blocked_base pass refinements/level_* lists. The three
blocked-base rejection cases keep their inputs per level; base-only
blocks now meet the level-list rule (blocked refined levels below) where
the legacy keyword check fired before. Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:58:57 +00:00
ClaudeandClaude Opus 5.5 af0566db98 P1 (6/n): exchange-window, raw-exchange, Jet and pinned-remap GPU tests migrated
test_exchange_window, test_raw_exchange_gpu (one fine(block) helper for
the repeated circle+jet construction), test_refinement_jets,
test_multiple_pinned_remap_gpu and test_pinned_body_remap_gpu pass
refinements/level_* lists. Dense fine levels stay level_blocks=[None,
None|block] as before. Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:58:15 +00:00
ClaudeandClaude Opus 5.5 0b1b3d938f P1 (5/n): six GPU validation tests use the level-list constructor
test_refinement_gpu (incl. PeriodicPatchStepper alias calls),
test_channel_root_gpu, test_fine_surface_gpu, test_shifted_remap_gpu,
test_runtime_child_remap_gpu and test_refinement_three_level_multipatch_gpu
pass refinements/level_* lists; scalar and set children map to the second
refinement entry exactly as _legacy_levels did. Assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:48:09 +00:00
ClaudeandClaude Opus 5.5 5167c79378 P1 (4/n): unit tests and examples use the level-list constructor
test_shifted_jets, test_shifted_storage, test_raw_exchange and
test_channel_root pass refinements/level_* lists; examples
fixed_patch_jet and fixed_periodic_patch pass refinements=[BOUNDS]
instead of the positional bounds. Assertions kept; the shifted-storage
case base_block-without-fine_block now expects the level-list rule
'blocked refined levels below' (same constraint).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:46:52 +00:00
ClaudeandClaude Opus 5.5 aedddc35b1 P1 (3/n): three-level multipatch unit tests use the level-list constructor
KW (legacy base/fine/child blocks + child_rectangles_in_parent) becomes
levels(): refinements=[PARENTS, CHILDREN], level_blocks, level_objects,
level_jets. Every negative case is kept and expressed per level (child
rectangles, per-level blocks); the L1-body case now expects the level-list
rule's message 'only on the finest level' instead of the legacy keyword
message 'only on L2' (same constraint).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:46:06 +00:00
ClaudeandClaude Opus 5.5 7821c6ad71 P1 (2/n): runners and example use the level-list constructor
run_refined_dfg, run_dfg_re20, time_window_dfg, time_refined_dfg,
run_jet_circle_cv and examples/fixed_patch_circle pass refinements /
level_blocks / level_objects / level_jets instead of bounds, base_block,
fine_block, fine_objects and jets. One-to-one with _legacy_levels: same
rectangle (FixedPatch2D frame), same blocks, bodies and Jets.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:35:04 +00:00
ClaudeandClaude Opus 5.5 392314e1ae P1 (1/n): remap_child builds its candidate through the level-list API
The only in-library caller of the legacy constructor keywords now passes
refinements=[L1 bounds, L2 set], level_blocks=[base, fine, child],
level_objects=[(), (), child_objects], level_jets=[(), (), jets] -- the
exact lists _legacy_levels produced for it. Legacy keywords stay until
every caller is migrated (final P1 commit).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:34:33 +00:00
ClaudeandClaude Opus 5.5 f3f1e3d588 PeriodicBlockStepper.step: mid-batch launch failure test; async-error docstring note
Review follow-up (non-blocking items): a launch failure on the second of
three steps fails closed without the batch wait, keeps step_count at its
entry value and blocks download; the docstring now states, as for
FixedPatchStepper, that async faults surface at the single wait.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:31:09 +00:00
ClaudeandClaude Opus 5.5 3648ad70ef PeriodicBlockStepper.step(n): wait on the stream once per call
Same change as the reviewed FixedPatchStepper end-sync (cherry-picked
ad1d82c+512dcc9 above): the n steps are enqueued in order and the host
waits once; step(0) launches nothing; step_count advances after the wait
and keeps its entry value on failure. Launches, buffer roles and bits are
those of n step(1) calls. Mocked-driver host test covers the single
wait, launch/role equality with n step(1), step(0) and failure.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:25:16 +00:00
ClaudeandClaude Opus 5.5 5399c9c0a6 step(n) docstring: raw-exchange ledger holds the final interval only
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:24:36 +00:00
ClaudeandClaude Opus 5.5 54be33e006 FixedPatchStepper.step(n): wait on the stream once per call
Enqueue all n base intervals in order on the one stream and synchronize
once before returning (step(0): no launch, no wait). step_count, the raw
exchange ledger and an open window are published only after that wait;
any failure leaves step_count at its entry value and invalidates both
ledgers. Async device errors are no longer attributed to an interval;
the docstring says to use step(1) loops to localize one. Launch order and
buffers are unchanged. Adds a mocked-driver host test of the contract.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:24:36 +00:00
43 changed files with 642 additions and 386 deletions
+13 -12
View File
@@ -112,8 +112,9 @@ from CelerisLab.lbm.refinement_stepper import FixedPatchStepper
cfg = LBMConfig(nx=32, ny=24, viscosity=.1)
body = SimObject(0, geometry=CircleGeometry(16.0, 12.0, 2.15), center=(16.0, 12.0), radius=2.15,
surface_actuation=SurfaceActuation(0.7, 0.85, 0.002)) # angle, half-width, speed
with CudaContext(0) as ctx, FixedPatchStepper(cfg, (8, 6, 24, 18), context=ctx,
fine_objects=[body], fine_block=(16, 20)) as run:
with CudaContext(0) as ctx, FixedPatchStepper(cfg, context=ctx, refinements=[(8, 6, 24, 18)],
level_blocks=[None, (16, 20)],
level_objects=[(), (body,)]) as run:
rest = np.asarray(compute_feq_d2q9(1.0, 0.0, 0.0), np.float32)[:, None, None]
# (9, ny, nx) postcollision populations per level; non-leaf cells are discarded.
run.initialize(np.broadcast_to(rest, (9, *run.patch.coarse_shape)),
@@ -133,19 +134,19 @@ Scope and limits:
- Bodies must be stationary circles whose whole boundary stencil lies in the fine
active region, clear of the coarse/fine transfer cells; violations are rejected
before any GPU allocation.
- `fine_block=(by, bx)` stores the fine level in fixed-size blocks; `by`/`bx` must
- `level_blocks=[None, (by, bx)]` stores the fine level in fixed-size blocks; `by`/`bx` must
divide the WHOLE fine allocation, i.e. `2*(y1-y0)+8` by `2*(x1-x0)+8` including the
four-cell ghost band. It gives results identical to dense storage and is
infrastructure for future block AMR, not a speedup.
- Several regions: pass a sequence of `(x0, y0, x1, y1)` rectangles as
`bounds`, with both `fine_block` and `base_block`. Transfers are built once
`refinements[0]`, with both levels blocked (`level_blocks`). Transfers are built once
from their union; adjacent, overlapping and concave unions are accepted when
their actual dependencies close. There is no fixed rectangle-gap rule.
The fine level is then one level-wide grid `(2*ny, 2*nx)`:
`run.patch.fine_shape` gives it, `initialize`/`download` use it, and only the
grid-aligned tiles intersecting the union plus its four-cell ghost band are stored.
Block sizes must divide the level-wide grid; rectangle edges need not align
with tiles. Order of rectangles does not matter. `fine_objects`/`jets` work
with tiles. Order of rectangles does not matter. L1 bodies/Jets (`level_objects[1]`/`level_jets[1]`) work
as for one patch: each circle (same coarse cell-edge frame) must lie wholly
inside the union's active region with
its boundary stencil clear of transfer cells; a jet may not select links of two
@@ -156,8 +157,8 @@ Scope and limits:
cfg = LBMConfig(nx=32, ny=16, viscosity=.1)
rest = np.asarray(compute_feq_d2q9(1.0, 0.0, 0.0), np.float32)[:, None, None]
with CudaContext(0) as ctx, FixedPatchStepper(
cfg, ((4, 4, 12, 12), (16, 4, 24, 12)), context=ctx,
fine_block=(4, 4), base_block=(4, 4)) as run:
cfg, context=ctx, refinements=[((4, 4, 12, 12), (16, 4, 24, 12))],
level_blocks=[(4, 4), (4, 4)]) as run:
run.initialize(np.broadcast_to(rest, (9, *run.patch.coarse_shape)),
np.broadcast_to(rest, (9, *run.patch.fine_shape))) # (9, 32, 64)
run.step(32)
@@ -166,12 +167,12 @@ Scope and limits:
```
- Legacy three-level scalar chain: add `child=(x0, y0, x1, y1)` in
**L1 window cells** (the patch's fine grid including its 4-cell ghost band) and
`child_block`, with `fine_block` and `base_block`. L2 must sit inside L1's
all three `level_blocks` entries. L2 must sit inside L1's
ACTIVE cells, away from L1's ghost band and transfer footprint (checked before
any GPU allocation). `initialize` and `download` take/return three arrays; each
level is NaN on its non-ACTIVE (and child-covered) cells. One coarse step = 2 L1
steps = 4 L2 steps. One L1 patch only. Optional stationary circles and Jets
live on L2 via `child_objects` and `child_jets`; `fine_objects`/`jets` retain
live on L2 via `level_objects[2]` and `level_jets[2]`; L1 bodies/Jets retain
their L1 meaning and must be empty in this mode. All object and Jet geometry
uses the **L0 edge frame**, not the child bounds' L1 window frame. L2 has
dx=dt=0.25, so prescribed lattice velocity is unchanged. Each whole body and
@@ -188,9 +189,9 @@ Scope and limits:
jet = Jet(Tube((6.5, 8.), 0., 3., .3), Parallel(0.), amplitude=.002)
rest = np.asarray(compute_feq_d2q9(1.0, 0.0, 0.0), np.float32)[:, None, None]
with CudaContext(0) as ctx, FixedPatchStepper(
cfg, (4, 4, 12, 12), context=ctx, fine_block=(4, 4), base_block=(4, 4),
child=(8, 8, 16, 16), child_block=(4, 4),
child_objects=[body], child_jets=[jet]) as run:
cfg, context=ctx, refinements=[(4, 4, 12, 12), (8, 8, 16, 16)],
level_blocks=[(4, 4)] * 3, level_objects=[(), (), (body,)],
level_jets=[(), (), (jet,)]) as run:
shapes = (run.patch.coarse_shape, run.patch.fine_shape, run.child_patch.fine_shape)
run.initialize(*(np.broadcast_to(rest, (9, *s)) for s in shapes)) # 16², 24², 24²
run.step(4)
+1 -1
View File
@@ -1,6 +1,6 @@
# Blocked AMR root channel increment
Opt in with `FixedPatchStepper(..., coarse_boundary="channel", base_block=(by,bx), fine_block=(fy,fx))`, or the existing level-list form. Every level must be blocked. The supported boundary is the existing uniform D2Q9/SRT/FP32 double-buffer channel: local Zou-He west inlet, pressure-NEQ east outlet, stationary halfway bounce-back y walls. Uniform and parabolic inlet profiles and shifted/unshifted storage are supported; the existing configuration guards still reject MRT/LES and incompatible storage. This does not change the external DFG MRT reference.
Opt in with `FixedPatchStepper(..., coarse_boundary="channel", refinements=[bounds], level_blocks=[(by,bx), (fy,fx)])`. Every level must be blocked. The supported boundary is the existing uniform D2Q9/SRT/FP32 double-buffer channel: local Zou-He west inlet, pressure-NEQ east outlet, stationary halfway bounce-back y walls. Uniform and parabolic inlet profiles and shifted/unshifted storage are supported; the existing configuration guards still reject MRT/LES and incompatible storage. This does not change the external DFG MRT reference.
`step/block_channel.cu` is a small root-only kernel. Fine kernels, transfer kernels and recursive hop scheduling are unchanged. Global coordinates and the global profile height/velocity are separate from block-local storage strides. Corner ownership and collision phases match `OneStep`:
+6 -6
View File
@@ -16,8 +16,8 @@ Keep the existing scalar `child=(x0,y0,x1,y1)` for a single scalar parent patch.
Its coordinates retain their original meaning: half-open integer indices in the
allocated L1 window, including its ghost band.
The new keyword `child_rectangles_in_parent=[...]` is mutually exclusive with
`child`. It accepts a scalar or a set parent and interprets every rectangle in
A rectangle list as the second `refinements` entry (formerly the keyword
`child_rectangles_in_parent=[...]`) accepts a scalar or a set parent and interprets every rectangle in
the actual allocated L1 array, including ghost bands and holes. It never infers
region-local coordinates. The physical root cell edge corresponding to an L1
index `(x,y)` is `parent.fine_origin + (x,y)/2`:
@@ -47,9 +47,9 @@ tile edges.
```python
run = FixedPatchStepper(
cfg, [(4, 4, 20, 28), (24, 4, 44, 28)], context=ctx,
base_block=(4, 4), fine_block=(4, 4), child_block=(4, 4),
child_rectangles_in_parent=[(16, 16, 32, 40), (56, 16, 80, 40)],
cfg, context=ctx,
refinements=[[(4, 4, 20, 28), (24, 4, 44, 28)], [(16, 16, 32, 40), (56, 16, 80, 40)]],
level_blocks=[(4, 4)] * 3,
)
```
@@ -57,7 +57,7 @@ This example uses a 48x32 root. `initialize` and `download` use three Q-major
dense arrays. Public arrays and masks still scale with the full allocated
extent; sparse device buffers do not remove that host-memory limit.
Stationary whole-circle `child_objects` and `child_jets` use root edge
Stationary whole-circle L2 bodies and Jets (`level_objects[2]`, `level_jets[2]`) use root edge
coordinates and a single L2 boundary object. Body IDs are physical object-list
positions, not region IDs or caller `obj_id`. Complete boundary stencils must
avoid the transfer footprint. Change a child Jet with
+2 -2
View File
@@ -90,8 +90,8 @@ def main():
# The caller owns the CUDA context; the stepper must be closed first (the
# nested with-statement does that).
with CudaContext(0) as context, FixedPatchStepper(
cfg, BOUNDS, context=context, coarse_boundary=args.mode,
fine_objects=[body], fine_block=fine_block) as run:
cfg, context=context, coarse_boundary=args.mode, refinements=[BOUNDS],
level_blocks=[None, fine_block], level_objects=[(), (body,)]) as run:
p = run.patch
# Rest equilibrium is already postcollision G, so no initial collision is
# needed. Non-leaf entries (covered, ghost, solid) are discarded by initialize.
+2 -2
View File
@@ -148,8 +148,8 @@ def main():
samples=[], status='running', last_synchronized_coarse_step=0)
start, run = time.monotonic(), None
try:
with CudaContext(0) as context, FixedPatchStepper(cfg, BOUNDS, context=context,
coarse_boundary='slot_reservoir') as run:
with CudaContext(0) as context, FixedPatchStepper(cfg, context=context,
coarse_boundary='slot_reservoir', refinements=[BOUNDS]) as run:
p = run.patch
masks = [np.ones(p.coarse_shape, bool), p.fine_kind == ACTIVE]
masks[0][[0, -1]] = False
+1 -1
View File
@@ -48,7 +48,7 @@ def main():
nx=16, ny=16, viscosity=.1, velocity=1e-3, rho=1.0,
lattice_model="D2Q9", collision="SRT", streaming="double_buffer",
store_precision="FP32", ddf_shifting=False, les_enabled=False)
with CudaContext(0) as context, PeriodicPatchStepper(cfg, (4, 4, 12, 12), context=context) as run:
with CudaContext(0) as context, PeriodicPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)]) as run:
p = run.patch
coarse = shear_equilibrium(p.coarse_shape, (0, 0), 1.0, cfg.ny)
fine = shear_equilibrium(p.fine_shape, p.fine_origin, .5, cfg.ny)
+14 -4
View File
@@ -204,11 +204,21 @@ class PeriodicBlockStepper:
raise
def step(self, n=1):
"""Advance n steps: halo refresh, then core pull/collide, per step."""
"""Advance n steps: halo refresh, then core pull/collide, per step.
All n steps are enqueued in order on the one stream and the host waits
once before returning (same launches and bits as n step(1) calls).
Asynchronous device errors surface at that single wait, so they are not
attributed to the step that caused them; call step(1) repeatedly to
localize one. On any failure the instance is failed and step_count keeps
its entry value.
"""
if isinstance(n, bool) or not isinstance(n, Integral) or n < 0:
raise ValueError("n must be a nonnegative integer")
with self._lock, self.ctx.activate():
self._require("ready")
if n == 0:
return
kw = dict(block=self._block, stream=self._stream)
try:
self._status = "advancing"
@@ -224,9 +234,9 @@ class PeriodicBlockStepper:
self._step(self._current, self._spare, np.uint32(self.layout.halo),
self._tables.get("flags", np.uintp(0)), np.uintp(0), # no kind
grid=self._grid, **kw)
self._stream.synchronize()
self._current, self._spare = self._spare, self._current
self._step_count += 1
self._current, self._spare = self._spare, self._current # host handles only
self._stream.synchronize()
self._step_count += n
self._status = "ready"
except BaseException:
self._status = "failed"
+5 -1
View File
@@ -178,7 +178,11 @@ extern "C" __global__ void BlockStep(const fpxx* in, fpxx* out, unsigned int hal
else if (gy == d_block_wall.ny - 2) reconstruct_top_reservoir(f, RHO);
}
#if COLLISION_MODEL == 2
collide_mrt_unforced(f); // SRT path below unchanged.
// Uniform MRT sequence (collide_dispatch, zero forcing); SRT path below unchanged.
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
@@ -114,16 +114,6 @@ __device__ __forceinline__ void collide_mrt_no_force(float* __restrict__ g,
collide_mrt(g, rho, ux, uy, Fin, omega);
}
// Unforced MRT with the clamped device omega, as in collide_dispatch with zero
// forcing. Shared by the blocked, channel-root and dense-patch collisions.
__device__ __forceinline__ void collide_mrt_unforced(float* __restrict__ g)
{
float rho, ux, uy;
compute_rho_u(g, rho, ux, uy);
collide_mrt_no_force(g, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
}
#elif LATTICE_MODEL == LATTICE_D3Q19
__device__ __forceinline__ void collide_mrt(float* __restrict__ g,
+4 -1
View File
@@ -37,7 +37,10 @@ __device__ unsigned int incoming_counts[NQ * TOTAL_CELLS];
__device__ __forceinline__ void patch_collide(float* f, unsigned long k) {
#if COLLISION_MODEL == 2
collide_mrt_unforced(f);
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
@@ -55,7 +55,11 @@ extern "C" __global__ void BlockChannelStep(
}
if (!wall && !outlet) {
#if COLLISION_MODEL == 2
collide_mrt_unforced(f);
// Match uniform MRT arithmetic (collide_dispatch, zero forcing).
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
#else
// Match uniform SRT arithmetic, including its shifted equilibrium path.
float rho, ux, uy, feq[NQ], forcing[NQ] = {};
+81 -142
View File
@@ -65,52 +65,8 @@ def _level_name(level):
return ("base", "fine", "child")[level] if level < 3 else f"level {level}"
def _legacy_levels(bounds, *, fine_block, base_block, child, child_block, child_rectangles_in_parent,
fine_objects, jets, child_objects, child_jets):
"""Normalise the original two/three-level keywords to per-level lists (original errors)."""
if child is not None and child_rectangles_in_parent is not None:
raise ValueError("child and child_rectangles_in_parent are mutually exclusive")
nested = child is not None or child_rectangles_in_parent is not None
bounds = tuple(bounds)
if not is_rectangle(bounds) and (fine_block is None or base_block is None):
raise ValueError("several patches require fine_block and base_block")
fine_objects, jets = tuple(fine_objects), tuple(jets)
child_objects, child_jets = tuple(child_objects), tuple(child_jets)
if not nested and (child_objects or child_jets):
raise ValueError("child_objects/child_jets require a nested child")
if base_block is not None and fine_block is None:
raise ValueError("base_block requires fine_block")
refinements, blocks = [bounds], [base_block, fine_block]
objects, level_jets = [(), fine_objects], [(), jets]
if nested:
if child is not None and not is_rectangle(bounds):
raise ValueError("scalar child requires one parent patch; use child_rectangles_in_parent")
if fine_block is None or base_block is None or child_block is None:
raise ValueError("a nested child needs fine_block, base_block and child_block")
if fine_objects or jets:
raise ValueError("a nested child supports bodies/Jets only on L2; use child_objects/child_jets")
refinements.append(tuple(child) if child is not None else tuple(child_rectangles_in_parent))
blocks.append(child_block)
objects.append(child_objects)
level_jets.append(child_jets)
return refinements, blocks, objects, level_jets
def _initial_levels(level_g, coarse_g, fine_g, child_g):
"""One tuple of per-level arrays from initialize()'s positional or legacy named form."""
if coarse_g is None and fine_g is None:
if child_g is not None and len(level_g) != 2:
raise TypeError("child_g needs exactly the coarse and fine arrays before it")
return tuple(level_g) + (() if child_g is None else (child_g,))
if level_g:
raise TypeError("pass level arrays positionally or as coarse_g/fine_g/child_g, not both")
if coarse_g is None or fine_g is None:
raise TypeError("coarse_g and fine_g are both required in the named form")
return (coarse_g, fine_g) + (() if child_g is None else (child_g,))
def _check_levels(refinements, blocks, objects, jets):
"""Candidate limits shared by the list and the legacy entry points."""
"""Candidate limits of the level lists, checked before any device work."""
depth = len(refinements) + 1
if depth < 2:
raise ValueError("at least one refinement level is required")
@@ -134,7 +90,14 @@ def _check_levels(refinements, blocks, objects, jets):
class FixedPatchStepper:
"""D2Q9/SRT/FP32, one unforced internal 2:1 patch in an x-periodic box.
"""D2Q9 SRT/MRT FP32 fixed 2:1 refinement hierarchy, built from level lists.
``refinements[l]`` gives level l+1 in level-l integer cells: one rectangle
(x0, y0, x1, y1) (FixedPatch2D, a fine window with its own origin) or a
sequence of rectangles (FixedPatchSet2D, a dependency-validated union on one
level-wide grid, public fine arrays (2*ny, 2*nx), NaN outside ACTIVE cells);
see lbm.hierarchy. ``level_blocks``, ``level_objects`` and ``level_jets``
hold one entry per level, L0 included (L0 bodies/Jets are not supported).
``coarse_boundary`` is ``periodic`` (default) or ``slot_reservoir``:
a bottom on-site slot and top fixed-density reservoir from ``cfg``.
@@ -145,68 +108,53 @@ class FixedPatchStepper:
actual channel closure footprint. Channel remapping is not supported.
Periodic mode does not apply the wall/slot boundary fields in ``cfg``.
Outer walls remain coarse; fine halo donors must be ordinary fluid.
``fine_objects`` accepts fixed circles in coarse physical coordinates;
velocities are in coarse lattice units (dx=dt=1). Independent local copies
are projected at 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 (coarse centers I+.5), whereas ordinary Simulation objects use the integer-node frame.
``initialize(coarse_g, fine_g)`` accepts physical Q-major (9,ny,nx) arrays
already in postcollision phase. It does not collide or build channel flags.
Every ``step()`` publishes only after both levels complete on one stream.
On a device/launch failure the instance becomes unusable: close and recreate
it, rather than pretending the partially advanced fine buffer rolled back.
Bodies are fixed circles in root (L0) physical coordinates; velocities are
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
the actual two-hop transfer footprint. These use the edge/cell-center frame
(coarse centers I+.5), whereas ordinary Simulation objects use the
integer-node frame. Jets (body.jets.Jet) use the same frame and drive cut
links selected at their true hits; ``set_jet_amplitude(index, amplitude,
level=l)`` re-uploads the per-link velocity column only.
``initialize(*level_g)`` accepts one physical Q-major (9, ny, nx) array per
level, already in postcollision phase. It does not collide or build channel
flags. Every ``step()`` publishes only after all levels complete on one
stream. On a device/launch failure the instance becomes unusable: close and
recreate it, rather than pretending a partially advanced buffer rolled back.
A ``level_blocks`` entry (by, bx) stores that level in blocks.cu storage:
the ENTIRE level shape, including the four-cell patch ghost band, is
partitioned by a non-periodic BlockLayout2D with one-cell same-level halos.
Patch kinds stay per cell (OUTER invalid, INNER stream only, ACTIVE
collide); halos are gathered before each hop. Public arrays keep the dense
shapes; only the allocation differs. A blocked level requires blocked
refined levels below it; several L1 rectangles and depth >= 3 require all
levels blocked. Blocked bodies are stationary circles, optionally
surface-actuated: prepare_fine_boundary decides legality on the dense grid,
then cut links are staged receiver-local per hop: halo -> BFL -> step, same
buffer. A blocked parent gathers halos before receiver-local coalescence
staging and pull/collision; covered cells are poisoned on every step.
``fine_block=(by, bx)`` (opt-in) stores the fine level in
blocks.cu storage: the ENTIRE fine_shape, including the four-cell patch
ghost band, is partitioned by a non-periodic BlockLayout2D with one-cell
same-level halos. Patch kinds stay per cell (OUTER invalid, INNER stream
only, ACTIVE collide); halos are gathered before each fine hop. Public
arrays keep the dense (9, ny, nx) fine shape; only the allocation differs.
Blocked fine_objects are stationary circles, optionally surface-actuated:
prepare_fine_boundary decides legality on the dense fine grid, then cut
links are staged receiver-local per hop: halo -> BFL -> step, same buffer.
``base_block=(by, bx)`` also blocks the parent; it requires ``fine_block``.
Parent halos precede receiver-local coalescence staging and pull/collision;
covered cells are explicitly poisoned on every step. Public shapes stay dense.
``observe_raw_exchange=True`` records blocked physical cut-link exchange
without changing populations. read_raw_exchange reports the last complete
base interval or last substep in base-lattice units per unit depth. This
is boundary momentum reaction, not a validated nozzle/stress load. Dense
body layouts are not supported by this new diagnostic in this candidate.
body levels are not supported by this diagnostic.
``observe_window_exchange=True`` (requires observe_raw_exchange) adds a second
per-link FP64 ledger accumulated by the same kernel between
start_exchange_window() and read_exchange_window(); disabled by default.
``jets`` (body.jets.Jet) use the same edge frame as ``fine_objects`` and
drive fine cut links selected at their true hits; ``set_jet_amplitude``
re-uploads the per-link velocity column only.
Several patches: pass ``bounds`` as a sequence of (x0, y0, x1, y1)
rectangles (FixedPatchSet2D, dependency-validated union). This requires
fine_block and base_block; fine_objects/jets are projected once onto the
level-wide fine grid (each circle inside the union, body ids = object
positions, jet and writer checks level-wide), and
uses one level-wide fine grid (2*ny, 2*nx): public fine arrays have that
shape, NaN outside the patches' ACTIVE cells. The step sequence is the
blocked one above, with unique union transfers in one launch each.
Three-level separated unions: ``child_rectangles_in_parent=[...]`` is
mutually exclusive with scalar ``child`` and accepts either parent form.
Coordinates index the allocated L1 array, including ghost bands and holes:
root edge = patch.fine_origin + L1 index / 2. A parent set has origin zero;
a scalar parent retains its shifted window origin. The child set's public
extent is twice the allocated L1 shape, even though storage is sparse.
Each pair requires gap >= 4 in its own parent-cell units. All regions share
one buffer pair per level and the fixed 1:2:4 schedule; L2-only bodies/Jets
below apply unchanged, with object IDs assigned over the whole level.
Three-level single chain: ``child`` bounds use L1 window cells; all three
levels must be blocked. ``child_objects`` and ``child_jets`` use the L0
edge frame, with L2 dx=dt=.25 and unchanged lattice velocity. L1 bodies
and Jets remain unsupported in a chain. Every whole circle and its BFL
stencil must fit strictly away from the L2 transfer footprint. L2 runs
halo -> BFL -> step four times per base step. Runtime child Jet changes
require ``set_jet_amplitude(index, amplitude, level=2)``.
Nested levels: a child set's coordinates index the allocated parent array,
including ghost bands and holes: root edge = parent fine_origin + index / 2
(zero origin for a parent set; a scalar parent keeps its shifted window
origin). A scalar child requires a scalar parent. Each pair requires gap >= 4
in its own parent-cell units. All regions share one buffer pair per level
and the fixed schedule of 2**l hops per base step at level l. Bodies/Jets are supported only on the
finest level of a nested hierarchy (L2 dx=dt=.25 in root units, unchanged
lattice velocity); every whole circle and its BFL stencil must fit strictly
away from that level's transfer footprint.
``explosion="projected_linear"`` selects experimental transverse linear
reconstruction on each adjacent pair. Default "constant" is unchanged.
@@ -214,43 +162,21 @@ class FixedPatchStepper:
bodies/Jets or LES; this backend supplies no forcing. Geometry is generic,
but union/multilevel GPU and matched-physics qualification remain pending.
Level lists (any depth): ``refinements[l]`` gives level l+1 in level-l
integer cells (scalar window or rectangle list, see lbm.hierarchy);
``level_blocks``/``level_objects``/``level_jets`` hold one entry per level,
L0 included. They are mutually exclusive with the legacy keywords above,
which are normalised to the same lists. One recursive hop advances every
level: a child explodes from its parent's current hop source, takes two
hops, then the parent gathers halos, stages restriction, applies BFL and
steps. Depth >= 3 must be fully blocked with bodies/Jets on the finest level.
One recursive hop advances every level: a child explodes from its parent's
current hop source, takes two hops, then the parent gathers halos, stages
restriction, applies BFL and steps.
"""
def __init__(self, cfg, bounds=None, *, context, coarse_boundary="periodic", fine_objects=(),
fine_block=None, base_block=None, jets=(), child=None, child_block=None,
child_objects=(), child_jets=(), child_rectangles_in_parent=None,
refinements=None, level_blocks=None, level_objects=None, level_jets=None,
def __init__(self, cfg, *, context, refinements, coarse_boundary="periodic",
level_blocks=None, level_objects=None, level_jets=None,
observe_raw_exchange=False, observe_window_exchange=False, explosion="constant"):
if HALO != 4:
raise RuntimeError("refinement.cu requires a four-cell fine halo")
if coarse_boundary not in ("periodic", "slot_reservoir", "channel"):
raise ValueError("coarse_boundary must be periodic, slot_reservoir or channel")
if refinements is None:
if bounds is None:
raise ValueError("bounds or refinements is required")
if level_blocks is not None or level_objects is not None or level_jets is not None:
raise ValueError("level_blocks/level_objects/level_jets require refinements")
refinements, level_blocks, level_objects, level_jets = _legacy_levels(
bounds, fine_block=fine_block, base_block=base_block, child=child,
child_block=child_block, child_rectangles_in_parent=child_rectangles_in_parent,
fine_objects=fine_objects, jets=jets, child_objects=child_objects,
child_jets=child_jets)
else:
if (bounds is not None or child is not None or child_rectangles_in_parent is not None
or fine_block is not None or base_block is not None or child_block is not None
or tuple(fine_objects) or tuple(jets) or tuple(child_objects) or tuple(child_jets)):
raise ValueError("refinements and the legacy patch keywords are mutually exclusive")
depth = len(refinements) + 1
level_blocks = [None] * depth if level_blocks is None else list(level_blocks)
level_objects = [()] * depth if level_objects is None else list(level_objects)
level_jets = [()] * depth if level_jets is None else list(level_jets)
depth = len(refinements) + 1
level_blocks = [None] * depth if level_blocks is None else list(level_blocks)
level_objects = [()] * depth if level_objects is None else list(level_objects)
level_jets = [()] * depth if level_jets is None else list(level_jets)
level_objects = [tuple(o) for o in level_objects]
level_jets = [tuple(j) for j in level_jets]
level_blocks = [None if b is None else tuple(b) for b in level_blocks]
@@ -540,13 +466,13 @@ class FixedPatchStepper:
if self._ddf_shifting and any(layout is None for layout in self._layouts):
raise ValueError("shifted populations require all-blocked levels")
def initialize(self, *level_g, coarse_g=None, fine_g=None, child_g=None):
def initialize(self, *level_g):
"""Upload physical postcollision G, encoding before any FP32 rounding.
Shifted input should be physical FP64 to retain small perturbations.
Public initialization resets the step count; non-leaves are discarded.
"""
given = _initial_levels(level_g, coarse_g, fine_g, child_g)
given = tuple(level_g)
with self._lock, self.ctx.activate():
self._require("new", "ready")
self._check_storage()
@@ -670,10 +596,9 @@ class FixedPatchStepper:
old, new, source, base_count=count, level_steps=(count, 2*count, 4*count),
phase="postcollision_g", **options)
candidate = type(self)(
cfg, pairs[0].bounds, context=self.ctx, base_block=base_block,
fine_block=fine_block, child_block=child_block,
child_rectangles_in_parent=new.child.bounds,
child_objects=child_objects, child_jets=jets,
cfg, context=self.ctx, refinements=[pairs[0].bounds, new.child.bounds],
level_blocks=[base_block, fine_block, child_block],
level_objects=[(), (), child_objects], level_jets=[(), (), jets],
observe_raw_exchange=self._observe_raw_exchange,
observe_window_exchange=self._observe_window, explosion=self._explosion)
try:
@@ -822,27 +747,41 @@ class FixedPatchStepper:
grid=((pairs + self._block[0] - 1) // self._block[0], 1), **kw)
def step(self, n=1):
"""Advance G→G by n base intervals: level l takes 2**l hops per interval."""
"""Advance G→G by n base intervals: level l takes 2**l hops per interval.
All n intervals are enqueued in order on the one stream, and the host waits
once, before returning; step(0) launches and waits for nothing. On return
the device work is complete, step_count and an open window cover all n
intervals, and the raw-exchange ledger holds the final interval only.
Asynchronous device errors surface at that single wait, so they are not
attributed to the interval that caused them; call step(1) repeatedly to
localize one (same launches, same bits). On any failure the instance is
failed, the exchange ledger and window are invalid, and step_count keeps
its value from entry to the failing call.
"""
if isinstance(n, bool) or not isinstance(n, Integral) or n < 0:
raise ValueError("n must be a nonnegative integer")
with self._lock, self.ctx.activate():
self._require("ready")
if n == 0:
return
root = self._levels[0]
kw = dict(block=self._block, stream=self._stream)
try:
self._status = "advancing"
self._exchange_valid = False
for _ in range(n):
if self._observe_raw_exchange:
self._exchange_valid = False
if self._observe_raw_exchange: # stream-ordered: holds the last interval only
for b in self._level_boundaries:
if b is not None:
cuda.memset_d8_async(b.interval_integral, 0,
len(b.tables["cell"]) * 32, self._stream)
self._hop(0, root.current, root.spare, **kw)
self._stream.synchronize()
root.current, root.spare = root.spare, root.current
self._step_count += 1
self._exchange_valid = self._observe_raw_exchange
root.current, root.spare = root.spare, root.current # host handles only
self._stream.synchronize()
self._step_count += n
self._exchange_valid = self._observe_raw_exchange
self._status = "ready"
except BaseException:
self._window_valid = False
+102
View File
@@ -0,0 +1,102 @@
"""Host contract of PeriodicBlockStepper.step(n): one stream wait per call (mocked driver)."""
import ast
from contextlib import nullcontext
from pathlib import Path
import types
from unittest.mock import Mock
import unittest
import numpy as np
import CelerisLab.lbm.block_stepper as b
from CelerisLab.config import LBMConfig
from CelerisLab.lbm.blocks import BlockLayout2D
def host_module():
"""block_stepper with pycuda/compiler replaced by mocks; no GPU context."""
path = Path(b.__file__)
tree = ast.parse(path.read_text())
tree.body = [node for node in tree.body if not (
isinstance(node, ast.Import) and any(a.name == "pycuda.driver" for a in node.names)
or isinstance(node, ast.ImportFrom) and node.module == "cuda")]
module = types.ModuleType("CelerisLab.lbm._test_host_block_stepper")
module.__package__ = "CelerisLab.lbm"
module.__dict__.update(cuda=Mock(), compiler=Mock())
exec(compile(tree, str(path), "exec"), module.__dict__)
return module
def backend():
host, log = host_module(), []
host.cuda.mem_get_info.return_value = (1 << 40, 1 << 40)
host.cuda.mem_alloc.side_effect = lambda size: Mock(size=size) # distinct buffers
module = Mock()
module.get_global.return_value = (Mock(), 20)
module.get_function.side_effect = lambda name: Mock(
side_effect=lambda *a, **k: log.append((name, a[:2])))
host.compiler.load_module.return_value = module
cfg = LBMConfig(nx=16, ny=16, viscosity=.1)
run = host.PeriodicBlockStepper(cfg, BlockLayout2D.build(16, 16, 8, 8),
context=Mock(activate=nullcontext, sm_arch="sm_70"))
run._stream.synchronize.side_effect = lambda: log.append(("sync", ()))
run.initialize(np.full((9, 16, 16), 1 / 9))
log.clear()
return run, log
class BlockStepSync(unittest.TestCase):
def test_step_n_waits_once_and_matches_n_single_steps(self):
many, log_many = backend()
many.step(3)
single, log_single = backend()
for _ in range(3):
single.step(1)
names = [name for name, _ in log_many]
self.assertEqual(names.count("sync"), 1)
self.assertEqual(names[-1], "sync")
self.assertEqual(names.count("BlockStep"), 3)
# Same launch order and buffer roles as three step(1) calls.
strip = lambda log: [(n, tuple(id(x) for x in a)) for n, a in log if n != "sync"]
roles = lambda run, log: [(n, tuple({id(run._current): "c", id(run._spare): "s"}.get(x, "-")
for x in ids)) for n, ids in strip(log)]
self.assertEqual([n for n, _ in strip(log_many)], [n for n, _ in strip(log_single)])
self.assertEqual(roles(many, log_many), roles(single, log_single))
self.assertIsNot(many._current, many._spare)
self.assertEqual((many.step_count, single.step_count), (3, 3))
def test_step_zero_launches_and_waits_for_nothing(self):
run, log = backend()
run.step(0)
self.assertEqual((log, run.step_count), ([], 0))
def test_failure_keeps_entry_step_count(self):
run, log = backend()
run.step(1)
run._stream.synchronize.side_effect = RuntimeError("device fault")
with self.assertRaises(RuntimeError):
run.step(2)
self.assertEqual((run.step_count, run._status), (1, "failed"))
def test_launch_failure_mid_batch_fails_closed_without_waiting(self):
run, log = backend()
calls = []
step = run._step
def failing(*args, **kwargs):
calls.append(1)
if len(calls) == 2:
raise RuntimeError("launch failure")
return step(*args, **kwargs)
run._step = failing
with self.assertRaises(RuntimeError):
run.step(3)
self.assertNotIn("sync", [name for name, _ in log])
self.assertEqual((run.step_count, run._status), (0, "failed"))
with self.assertRaises(RuntimeError):
run.download()
if __name__ == "__main__":
unittest.main()
+2 -2
View File
@@ -80,8 +80,8 @@ def test_evolved_ghost_snapshot_and_physical_download():
host.cuda.memcpy_dtoh.side_effect = lambda a, ptr: np.copyto(a, storage[id(ptr)])
context = Mock(activate=nullcontext, sm_arch="sm_70")
cfg = LBMConfig(nx=24, ny=20, viscosity=.1, velocity=.025)
backend = host.FixedPatchStepper(cfg, (6, 6, 18, 14), context=context,
coarse_boundary="channel", base_block=(4, 6), fine_block=(4, 4))
backend = host.FixedPatchStepper(cfg, context=context, coarse_boundary="channel",
refinements=[(6, 6, 18, 14)], level_blocks=[(4, 6), (4, 4)])
try:
weights = np.array([4/9]+[1/9]*4+[1/36]*4)[:, None, None]
arrays = [np.broadcast_to(weights, level.shape).copy() for level in backend._levels]
+3 -2
View File
@@ -157,8 +157,9 @@ class RawExchange(unittest.TestCase):
host = host_module()
ctx = Mock()
with self.assertRaisesRegex(ValueError, 'blocked body'):
host.FixedPatchStepper(LBMConfig(nx=32, ny=32, viscosity=.1), (4, 4, 28, 28),
context=ctx, fine_objects=[body(radius=2.15)], observe_raw_exchange=True)
host.FixedPatchStepper(LBMConfig(nx=32, ny=32, viscosity=.1), context=ctx,
refinements=[(4, 4, 28, 28)],
level_objects=[(), (body(radius=2.15),)], observe_raw_exchange=True)
ctx.activate.assert_not_called()
def test_all_depth_weights(self):
+1 -6
View File
@@ -37,13 +37,8 @@ class RefinementMRTHost(unittest.TestCase):
self.assertIn("STORE_PRECISION != 0", source)
self.assertIn('#if COLLISION_MODEL == 2\n#include "operators/collision_mrt.cuh"\n#endif', source)
self.assertIn("USE_DDF_SHIFTING ||", refinement) # dense/mixed shifted stays rejected
# One shared unforced MRT body; each blocked/patch path calls it exactly once.
mrt = (KERNELS / "operators/collision_mrt.cuh").read_text()
self.assertEqual(mrt.count("void collide_mrt_unforced("), 1)
for source in (blocks, channel, refinement):
self.assertEqual(source.count("collide_mrt_unforced(f);"), 1)
self.assertNotIn("collide_mrt_no_force(", source)
self.assertNotIn("compute_rho_u(f, rho, ux, uy);\n collide_mrt", source)
self.assertEqual(source.count("collide_mrt_no_force("), 1)
if __name__ == "__main__": unittest.main()
@@ -173,8 +173,10 @@ def bodies_and_jets(speed=.002):
return bodies,jets
KW=dict(base_block=(4,4),fine_block=(4,4),child_block=(4,4),
child_rectangles_in_parent=CHILDREN)
def levels(parents=PARENTS,children=CHILDREN,blocks=((4,4),)*3,objects=((),(),()),jets=((),(),())):
"""Level-list constructor keywords for the L0 set, L1 set and L2 set used here."""
return dict(refinements=[parents,children],level_blocks=list(blocks),
level_objects=list(objects),level_jets=list(jets))
def mock_backend(*, children=CHILDREN, bodies=False, fail_compile=False, fail_upload=False):
@@ -201,9 +203,8 @@ def mock_backend(*, children=CHILDREN, bodies=False, fail_compile=False, fail_up
module.cuda.memcpy_htod.side_effect=upload
objects,jets=bodies_and_jets() if bodies else ((),())
try:
run=module.FixedPatchStepper(CFG,PARENTS,context=context,
child_objects=objects,child_jets=jets,
**dict(KW,child_rectangles_in_parent=children))
run=module.FixedPatchStepper(CFG,context=context,
**levels(children=children,objects=((),(),objects),jets=((),(),jets)))
except RuntimeError:
for allocation in allocations:
allocation.free.assert_called_once_with()
@@ -214,28 +215,29 @@ def mock_backend(*, children=CHILDREN, bodies=False, fail_compile=False, fail_up
class ConstructorAndSchedule(unittest.TestCase):
def test_invalid_inputs_fail_before_context_or_allocations(self):
context=Mock()
cases=[dict(child=(8,8,16,16)),dict(child_rectangles_in_parent=[]),
dict(child_rectangles_in_parent=[(16,16,32)]),
dict(child_rectangles_in_parent=[(16.5,16,32,40)]),
dict(child_rectangles_in_parent=[(16,16,32,24),(16,27,32,55)]),
dict(child_rectangles_in_parent=[(40,16,48,40)]),
dict(child_rectangles_in_parent=[(8,16,32,40)]),
dict(child_block=None),dict(fine_block=None),dict(base_block=None),
dict(child_block=(5,4)),dict(fine_block=(5,4)),
dict(child_block=(4,)),dict(child_rectangles_in_parent=(16,16,32,40))]
b=(4,4)
cases=[dict(children=(8,8,16,16)),dict(children=[]),
dict(children=[(16,16,32)]),
dict(children=[(16.5,16,32,40)]),
dict(children=[(16,16,32,24),(16,27,32,55)]),
dict(children=[(40,16,48,40)]),
dict(children=[(8,16,32,40)]),
dict(blocks=(b,b,None)),dict(blocks=(b,None,b)),dict(blocks=(None,b,b)),
dict(blocks=(b,b,(5,4))),dict(blocks=(b,(5,4),b)),
dict(blocks=(b,b,(4,))),dict(children=(16,16,32,40))]
for changes in cases:
with self.subTest(changes=changes),self.assertRaises((ValueError,TypeError)):
host.FixedPatchStepper(CFG,PARENTS,context=context,**dict(KW,**changes))
host.FixedPatchStepper(CFG,context=context,**levels(**changes))
with self.assertRaisesRegex(ValueError,'margin'):
host.FixedPatchStepper(CFG,[(1,4,20,28),(23,4,44,28)],context=context,**KW)
host.FixedPatchStepper(CFG,context=context,**levels(parents=[(1,4,20,28),(23,4,44,28)]))
objects,jets=bodies_and_jets()
with self.assertRaisesRegex(ValueError,'only on L2'):
host.FixedPatchStepper(CFG,PARENTS,context=context,fine_objects=objects,**KW)
with self.assertRaisesRegex(ValueError,'only on the finest level'):
host.FixedPatchStepper(CFG,context=context,**levels(objects=((),objects,())))
from CelerisLab.body.jets import Jet,Tube,Parallel
spanning=Jet(Tube((10.5,14.),0.,25.,.3),Parallel(0.),amplitude=0.)
with self.assertRaisesRegex(ValueError,'one jet'):
host.FixedPatchStepper(CFG,PARENTS,context=context,child_objects=objects,
child_jets=[spanning],**KW)
host.FixedPatchStepper(CFG,context=context,
**levels(objects=((),(),objects),jets=((),(),[spanning])))
context.activate.assert_not_called()
host.cuda.mem_alloc.assert_not_called()
@@ -313,21 +315,23 @@ class ConstructorAndSchedule(unittest.TestCase):
class InitialLevels(unittest.TestCase):
"""initialize() entry adapter: positional per-level arrays or the legacy names, never both."""
class RemovedLegacyAPI(unittest.TestCase):
"""The pre-level-list constructor/initialize forms are gone and fail loudly."""
def test_forms(self):
a, b, c = object(), object(), object()
f = host._initial_levels
self.assertEqual(f((a, b), None, None, None), (a, b))
self.assertEqual(f((a, b, c), None, None, None), (a, b, c))
self.assertEqual(f((a, b), None, None, c), (a, b, c)) # old initialize(c, f, child_g=g)
self.assertEqual(f((), a, b, None), (a, b))
self.assertEqual(f((), a, b, c), (a, b, c))
for args in (((a,), None, None, c), ((a, b, c), None, None, c), ((a,), a, b, None),
((), a, None, None), ((), None, b, c)):
with self.subTest(args=args), self.assertRaises(TypeError):
f(*args)
def test_removed_constructor_forms_raise_type_error(self):
import inspect
context = Mock()
with self.assertRaises(TypeError): # positional bounds
host.FixedPatchStepper(CFG, PARENTS, context=context)
removed = dict(bounds=PARENTS, fine_block=(4, 4), base_block=(4, 4), fine_objects=(),
jets=(), child=(8, 8, 16, 16), child_block=(4, 4), child_objects=(),
child_jets=(), child_rectangles_in_parent=CHILDREN)
for name, value in removed.items():
with self.subTest(keyword=name), self.assertRaises(TypeError):
host.FixedPatchStepper(CFG, context=context, **levels(), **{name: value})
context.activate.assert_not_called()
self.assertEqual(list(inspect.signature(host.FixedPatchStepper.initialize).parameters),
["self", "level_g"]) # no coarse_g/fine_g/child_g
if __name__=='__main__':
+5 -6
View File
@@ -54,14 +54,13 @@ class ShiftedJets(unittest.TestCase):
(Jet(Tube((10.,14.),0.,26.,.3),Parallel(0.)),))]
for objects,jets in cases:
with self.subTest(objects=len(objects),jets=len(jets)),self.assertRaises(ValueError):
module.FixedPatchStepper(CFG,backend.PARENTS,context=None,
base_block=(4,4),fine_block=(4,4),child_block=(4,4),
child_rectangles_in_parent=A,child_objects=objects,child_jets=jets)
module.FixedPatchStepper(CFG,context=None,refinements=[backend.PARENTS,A],
level_blocks=[(4,4)]*3,level_objects=[(),(),objects],level_jets=[(),(),jets])
for blocks in ((None,None,None),((4,4),None,(4,4))):
with self.assertRaises(ValueError):
module.FixedPatchStepper(CFG,backend.PARENTS,context=None,
base_block=blocks[0],fine_block=blocks[1],child_block=blocks[2],
child_rectangles_in_parent=A,child_objects=(body(),),child_jets=(jet(),))
module.FixedPatchStepper(CFG,context=None,refinements=[backend.PARENTS,A],
level_blocks=list(blocks),level_objects=[(),(),(body(),)],
level_jets=[(),(),(jet(),)])
module.cuda.mem_alloc.assert_not_called()
module.compiler.compile_kernel.assert_not_called()
+4 -3
View File
@@ -129,9 +129,10 @@ class ShiftedStorage(unittest.TestCase):
def test_dense_or_mixed_shifted_rejected_before_compilation(self):
module=backend.host_module()
cfg=replace(backend.CFG,ddf_shifting=True)
for kwargs in ({},dict(base_block=(4,4)),dict(fine_block=(4,4))):
with self.assertRaisesRegex(ValueError,'all-blocked|base_block requires fine_block'):
module.FixedPatchStepper(cfg,(4,4,20,28),context=None,**kwargs)
for blocks in ([None,None],[(4,4),None],[None,(4,4)]):
with self.assertRaisesRegex(ValueError,'all-blocked|blocked refined levels below'):
module.FixedPatchStepper(cfg,context=None,refinements=[(4,4,20,28)],
level_blocks=blocks)
module.compiler.compile_kernel.assert_not_called()
def test_uniform_block_codec_and_raw_snapshot(self):
+170
View File
@@ -0,0 +1,170 @@
"""Host contract of FixedPatchStepper.step(n): one stream wait per call (mocked driver).
All n base intervals are enqueued in order on the backend stream and the host
waits once before step returns. Launch order and buffer arguments must equal n
step(1) calls; completion metadata, raw-exchange and window ledgers, and failure
invalidation are published only after that wait.
"""
from contextlib import nullcontext
from unittest.mock import Mock
import unittest
import numpy as np
from CelerisLab.config import LBMConfig
from tests.unit.test_raw_exchange import body
from tests.unit.test_refinement_three_level_multipatch import host_module
BODY_LEVEL, HOPS = 3, 8 # body on L3: 2**3 hops per base interval
def backend(window=False):
"""Four-level host backend; every launch, memset and stream wait goes to one log."""
host = host_module()
log, allocations = [], []
def allocate(size):
allocations.append(Mock(size=size))
return allocations[-1]
module = Mock()
module.get_global.return_value = (Mock(), 20)
module.get_function.side_effect = lambda name: Mock(
side_effect=lambda *a, **k: log.append((name, a, k)))
host.compiler.load_module.return_value = module
host.cuda.mem_alloc.side_effect = allocate
host.cuda.memcpy_dtoh.side_effect = lambda array, pointer: array.fill(0)
host.cuda.memset_d8_async.side_effect = lambda *a: log.append(("memset", a, {}))
run = host.FixedPatchStepper(
LBMConfig(nx=48, ny=48, viscosity=.02), context=Mock(activate=nullcontext, sm_arch="sm_70"),
refinements=[[(6, 6, 26, 26)], [(20, 20, 44, 44)], [(56, 56, 72, 72)]],
level_blocks=[(4, 4)] * 4, level_objects=[(), (), (), (body(),)],
observe_raw_exchange=True, observe_window_exchange=window)
run._stream.synchronize.side_effect = lambda: log.append(("sync", (), {}))
run.initialize(*[np.full(level.shape, 1 / 9, np.float32) for level in run._levels])
log.clear()
return run, log, allocations
def signature(log, allocations):
"""Event names with device-pointer arguments replaced by allocation index."""
index = {id(a): i for i, a in enumerate(allocations)}
return [(name, tuple(index.get(id(a), "-") for a in args)) for name, args, _ in log]
def names(log):
return [name for name, _, _ in log]
class StepSync(unittest.TestCase):
def test_step1_waits_once_after_all_launches(self):
run, log, _ = backend()
run.step(1)
self.assertEqual(names(log).count("sync"), 1)
self.assertEqual(names(log)[-1], "sync")
self.assertEqual(names(log)[0], "memset")
self.assertEqual(names(log).count("BlockStep"), 1 + 2 + 4 + 8)
self.assertEqual(run.step_count, 1)
run.close()
def test_step_n_equals_n_step1_launches_with_one_wait(self):
single, single_log, single_alloc = backend()
for _ in range(3):
single.step(1)
batched, batched_log, batched_alloc = backend()
batched.step(3)
expected = [e for e in signature(single_log, single_alloc) if e[0] != "sync"]
got = signature(batched_log, batched_alloc)
self.assertEqual(got[-1][0], "sync")
self.assertEqual(names(batched_log).count("sync"), 1)
self.assertEqual(got[:-1], expected) # same order, same ping-pong buffers
self.assertEqual(names(batched_log).count("memset"), 3)
self.assertEqual(batched.step_count, single.step_count)
index = [{id(a): i for i, a in enumerate(alloc)} for alloc in (single_alloc, batched_alloc)]
self.assertEqual(index[1][id(batched._levels[0].current)],
index[0][id(single._levels[0].current)])
single.close()
batched.close()
def test_step0_launches_and_waits_for_nothing(self):
run, log, _ = backend()
run.step(2)
log.clear()
run.step(0)
self.assertEqual(log, [])
self.assertEqual(run.step_count, 2)
value = run.read_raw_exchange(0, level=BODY_LEVEL) # step(0) keeps the ledger valid
self.assertEqual((value["time_start"], value["time_end"]), (1, 2))
run.close()
def test_completion_metadata_covers_all_intervals(self):
run, _, _ = backend()
with self.assertRaisesRegex(RuntimeError, "completed"):
run.read_raw_exchange(0, level=BODY_LEVEL)
run.step(3)
value = run.read_raw_exchange(0, level=BODY_LEVEL)
self.assertEqual((value["time_start"], value["time_end"], value["sample_count"]),
(2, 3, HOPS)) # step(n) retains only the final base interval
run.step(2)
self.assertEqual(run.read_raw_exchange(0, level=BODY_LEVEL)["time_end"], 5)
run.close()
def test_window_accounting_across_batched_calls(self):
windows = []
for calls in ((1, 3), (4,), (1, 1, 1, 1)):
run, log, _ = backend(window=True)
run.step(1)
run.start_exchange_window()
log.clear()
for n in calls:
run.step(n)
self.assertEqual(names(log).count("BlockAccumulateExchange"), 2 * 4 * HOPS)
self.assertEqual(names(log).count("sync"), len(calls))
out = run.read_exchange_window(0, level=BODY_LEVEL)
windows.append((out["time_start"], out["time_end"], out["sample_count"]))
run.close()
self.assertEqual(windows, [(1., 5, 4 * HOPS)] * 3)
def test_failed_wait_publishes_nothing(self):
run, _, _ = backend(window=True)
run.step(2)
run.start_exchange_window()
run._stream.synchronize.side_effect = RuntimeError("device failure")
with self.assertRaisesRegex(RuntimeError, "device failure"):
run.step(3)
self.assertEqual(run.step_count, 2) # value at entry to the failing call
self.assertFalse(run._exchange_valid)
self.assertFalse(run._window_valid)
for read in (lambda: run.read_raw_exchange(0, level=BODY_LEVEL),
lambda: run.read_exchange_window(0, level=BODY_LEVEL), run.download,
lambda: run.step(1)):
with self.assertRaisesRegex(RuntimeError, "failed"):
read()
run._stream.synchronize.side_effect = None
run.close()
def test_launch_error_mid_batch_skips_wait_and_fails(self):
run, log, _ = backend()
run.step(1)
log.clear()
hop, calls = run._hop, []
def failing(index, *args, **kw):
if index == 0:
calls.append(1)
if len(calls) == 2:
raise RuntimeError("launch failure")
return hop(index, *args, **kw)
run._hop = failing
with self.assertRaisesRegex(RuntimeError, "launch failure"):
run.step(3)
self.assertNotIn("sync", names(log))
self.assertEqual(run.step_count, 1)
self.assertFalse(run._exchange_valid)
with self.assertRaisesRegex(RuntimeError, "failed"):
run.read_raw_exchange(0, level=BODY_LEVEL)
run.close()
if __name__ == "__main__":
unittest.main()
@@ -106,9 +106,10 @@ class ChannelRootGPU(unittest.TestCase):
cfg = LBMConfig(nx=48, ny=32, viscosity=.1, velocity=.015, ddf_shifting=shifted)
reference = None
for block in ((32,48), (8,8), (4,1)):
with FixedPatchStepper(cfg, (8,6,40,26), context=context,
coarse_boundary="channel", base_block=block, fine_block=(8,8),
fine_objects=(body,) if bodies else (), jets=(jet,) if bodies else ()) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary="channel",
refinements=[(8,6,40,26)], level_blocks=[block, (8,8)],
level_objects=[(), (body,) if bodies else ()],
level_jets=[(), (jet,) if bodies else ()]) as backend:
arrays = [state(*level.shape[1:]) for level in backend._levels]
backend.initialize(*arrays)
backend.step(4)
@@ -25,8 +25,8 @@ EPS = np.finfo(np.float64).eps
def build(m, window=True, raw=True):
return FixedPatchStepper(config("periodic"), BOUNDS, context=m["context"], base_block=(3, 4),
fine_block=(16, 20), fine_objects=[circle()],
return FixedPatchStepper(config("periodic"), context=m["context"], refinements=[BOUNDS],
level_blocks=[(3, 4), (16, 20)], level_objects=[(), (circle(),)],
observe_raw_exchange=raw, observe_window_exchange=window)
@@ -109,10 +109,11 @@ class ExchangeWindowGPU(unittest.TestCase):
from test_raw_exchange_gpu import BOUNDS as RB, CFG, circle as rc, jet, state
from CelerisLab.cuda.context import CudaContext
with CudaContext(0) as ctx:
kw = dict(context=ctx, fine_block=(16, 20), fine_objects=[rc()], jets=[jet()],
kw = dict(context=ctx, refinements=[RB], level_blocks=[None, (16, 20)],
level_objects=[(), (rc(),)], level_jets=[(), (jet(),)],
observe_raw_exchange=True)
with FixedPatchStepper(CFG, RB, observe_window_exchange=True, **kw) as windowed, \
FixedPatchStepper(CFG, RB, **kw) as plain:
with FixedPatchStepper(CFG, observe_window_exchange=True, **kw) as windowed, \
FixedPatchStepper(CFG, **kw) as plain:
arrays = [state(level.shape) for level in windowed._levels]
windowed.initialize(*arrays)
plain.initialize(*arrays)
@@ -200,11 +201,11 @@ class ExchangeWindowShiftedMRTGPU(unittest.TestCase):
from CelerisLab.cuda.context import CudaContext
cfg = LBMConfig(nx=32, ny=24, viscosity=.1, collision="MRT", ddf_shifting=True)
with CudaContext(0) as ctx:
kw = dict(context=ctx, base_block=(12, 8), fine_block=(16, 20), fine_objects=[rc()],
observe_raw_exchange=True)
with FixedPatchStepper(cfg, RB, observe_window_exchange=True, **kw) as windowed, \
FixedPatchStepper(cfg, RB, observe_window_exchange=True, **kw) as split, \
FixedPatchStepper(cfg, RB, **kw) as plain:
kw = dict(context=ctx, refinements=[RB], level_blocks=[(12, 8), (16, 20)],
level_objects=[(), (rc(),)], observe_raw_exchange=True)
with FixedPatchStepper(cfg, observe_window_exchange=True, **kw) as windowed, \
FixedPatchStepper(cfg, observe_window_exchange=True, **kw) as split, \
FixedPatchStepper(cfg, **kw) as plain:
arrays = [state(level.shape) for level in windowed._levels]
for run in (windowed, split, plain):
run.initialize(*arrays)
@@ -140,7 +140,8 @@ def run_patch(speed, *, audit=True):
build = (patch.object(compiler, "compile_kernel", audit_compiler(compiler.compile_kernel))
if audit else nullcontext())
with build:
backend = implementation.FixedPatchStepper(cfg, BOUNDS, context=context, fine_objects=[source])
backend = implementation.FixedPatchStepper(cfg, context=context, refinements=[BOUNDS],
level_objects=[(), (source,)])
with backend, context.activate():
c, f = backend._levels
cl = backend._fine_curved
@@ -325,6 +326,7 @@ def test_cross_interface_donor_rejects_before_cuda():
context.activate.side_effect = AssertionError("must reject before CUDA activation")
with patch.object(cuda, "mem_alloc") as allocate:
with pytest.raises(ValueError, match="donor crosses"):
implementation.FixedPatchStepper(cfg, BOUNDS, context=context, fine_objects=[source])
implementation.FixedPatchStepper(cfg, context=context, refinements=[BOUNDS],
level_objects=[(), (source,)])
context.activate.assert_not_called()
allocate.assert_not_called()
@@ -48,9 +48,8 @@ def make(context):
for (x, y), amplitude in zip(centers, (.003, -.002))
)
run = FixedPatchStepper(
CFG, PARENTS, context=context, base_block=(4, 4), fine_block=(4, 4),
child_block=(4, 4), child_rectangles_in_parent=START,
child_objects=objects, child_jets=jets,
CFG, context=context, refinements=[PARENTS, START], level_blocks=[(4, 4)] * 3,
level_objects=[(), (), objects], level_jets=[(), (), jets],
)
return run, objects
@@ -25,9 +25,8 @@ B = ((16,16,32,40), (60,16,76,40))
def make(ctx):
body = SimObject(71, geometry=CircleGeometry(12.,14.,.8),center=(12.,14.),radius=.8)
jet = Jet(Tube((10.5,14.),0.,3.,.3),Parallel(0.),'parabolic',.003)
run = FixedPatchStepper(CFG,PARENTS,context=ctx,base_block=(4,4),fine_block=(4,4),
child_block=(4,4),child_rectangles_in_parent=A,
child_objects=(body,),child_jets=(jet,))
run = FixedPatchStepper(CFG,context=ctx,refinements=[PARENTS,A],level_blocks=[(4,4)]*3,
level_objects=[(),(),(body,)],level_jets=[(),(),(jet,)])
return run,(body,)
@@ -90,14 +90,18 @@ def census(run, index, buffer):
return totals(dense, leaf), external, independent
def fine(block):
"""Two-level circle+jet construction on BOUNDS; block None keeps the dense fine level."""
return dict(refinements=[BOUNDS], level_blocks=[None, block],
level_objects=[(), (circle(),)], level_jets=[(), (jet(),)])
class RawExchangeGPU(unittest.TestCase):
def test_row_ledger_interval_jet_switch_and_no_population_effect(self):
with CudaContext(0) as ctx:
for block in ((32, 40), (16, 20), (8, 10)):
with FixedPatchStepper(CFG, BOUNDS, context=ctx, fine_block=block,
fine_objects=[circle()], jets=[jet()], observe_raw_exchange=True) as measured, \
FixedPatchStepper(CFG, BOUNDS, context=ctx, fine_block=block,
fine_objects=[circle()], jets=[jet()]) as plain:
with FixedPatchStepper(CFG, context=ctx, **fine(block), observe_raw_exchange=True) as measured, \
FixedPatchStepper(CFG, context=ctx, **fine(block)) as plain:
arrays = [state(l.shape) for l in measured._levels]
measured.initialize(*arrays)
plain.initialize(*arrays)
@@ -146,9 +150,8 @@ class RawExchangeGPU(unittest.TestCase):
"""Same fine resolution/initial data with dense storage: a numerical
reference for cut-link rows, not a uniform-whole-domain convergence study."""
with CudaContext(0) as ctx:
with FixedPatchStepper(CFG, BOUNDS, context=ctx, fine_objects=[circle()], jets=[jet()]) as dense, \
FixedPatchStepper(CFG, BOUNDS, context=ctx, fine_block=(16, 20),
fine_objects=[circle()], jets=[jet()], observe_raw_exchange=True) as blocked:
with FixedPatchStepper(CFG, context=ctx, **fine(None)) as dense, \
FixedPatchStepper(CFG, context=ctx, **fine((16, 20)), observe_raw_exchange=True) as blocked:
arrays = [state(l.shape) for l in blocked._levels]
dense.initialize(*arrays)
blocked.initialize(*arrays)
@@ -229,8 +232,7 @@ class RawExchangeGPU(unittest.TestCase):
layout = BlockLayout2D.build(48, 64, 8, 8)
with CudaContext(0) as ctx:
with PeriodicBlockStepper(fine_cfg, layout, context=ctx, objects=[obj]) as uniform, \
FixedPatchStepper(CFG, BOUNDS, context=ctx, fine_block=(16, 20),
fine_objects=[circle()], jets=[jet()], 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))
uniform.initialize(state((9, 48, 64)))
b = amr._level_boundaries[1]
@@ -83,11 +83,10 @@ def mocked_backend(mode="periodic", fail_alloc=None):
class BlockedBaseCPU(unittest.TestCase):
def test_reject_before_cuda(self):
context = Mock(activate=Mock(side_effect=AssertionError("CUDA reached")))
for kwargs in ({"base_block": (12, 8)},
{"base_block": (7, 8), "fine_block": FINE_BLOCK},
{"base_block": (12,), "fine_block": FINE_BLOCK}):
with self.subTest(kwargs=kwargs), self.assertRaises(ValueError):
FixedPatchStepper(config("periodic"), BOUNDS, context=context, **kwargs)
for blocks in ([(12, 8), None], [(7, 8), FINE_BLOCK], [(12,), FINE_BLOCK]):
with self.subTest(blocks=blocks), self.assertRaises(ValueError):
FixedPatchStepper(config("periodic"), context=context, refinements=[BOUNDS],
level_blocks=blocks)
context.activate.assert_not_called()
def test_all_arrival_rows_donors_and_kind_masks(self):
@@ -193,9 +192,10 @@ class BlockedBaseCPU(unittest.TestCase):
def test_launch_order_stream_wall_constants_and_buffer_swaps(self):
with mocked_backend("slot_reservoir") as m:
backend = FixedPatchStepper(config("slot_reservoir"), BOUNDS, context=m["context"],
coarse_boundary="slot_reservoir", fine_block=FINE_BLOCK,
base_block=(3, 4), fine_objects=[circle(actuation(.002))])
backend = FixedPatchStepper(config("slot_reservoir"), context=m["context"],
coarse_boundary="slot_reservoir", refinements=[BOUNDS],
level_blocks=[(3, 4), FINE_BLOCK],
level_objects=[(), (circle(actuation(.002)),)])
try:
backend.initialize(state(NY, NX, 1), state(*backend.patch.fine_shape, 2))
coarse, fine = backend._levels
@@ -236,16 +236,16 @@ class BlockedBaseCPU(unittest.TestCase):
for fail_at in (1, 6, 10):
with self.subTest(fail_at=fail_at), mocked_backend(fail_alloc=fail_at) as m:
with self.assertRaisesRegex(MemoryError, "injected allocation"):
FixedPatchStepper(config("periodic"), BOUNDS, context=m["context"],
fine_block=FINE_BLOCK, base_block=(3, 4))
FixedPatchStepper(config("periodic"), context=m["context"], refinements=[BOUNDS],
level_blocks=[(3, 4), FINE_BLOCK])
self.assertEqual(len(m["allocations"]), fail_at)
for allocation in m["allocations"]:
allocation.free.assert_called_once_with()
def test_failed_parent_stage_cannot_publish_or_resume(self):
with mocked_backend() as m:
backend = FixedPatchStepper(config("periodic"), BOUNDS, context=m["context"],
fine_block=FINE_BLOCK, base_block=(3, 4))
backend = FixedPatchStepper(config("periodic"), context=m["context"], refinements=[BOUNDS],
level_blocks=[(3, 4), FINE_BLOCK])
backend.initialize(state(NY, NX, 1), state(*backend.patch.fine_shape, 2))
m["functions"]["BlockCoalesce"].side_effect = RuntimeError("injected parent stage failure")
current = backend._levels[0].current
@@ -289,9 +289,10 @@ class BlockedBaseGPU(unittest.TestCase):
for base_block in (None, *layouts):
label = f"{mode} width={width} speed={speed} base={base_block}"
snapshots = []
with FixedPatchStepper(cfg, BOUNDS, context=context, coarse_boundary=mode,
fine_block=FINE_BLOCK, base_block=base_block,
fine_objects=() if speed is None else [circle(actuation(speed))]) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[BOUNDS], level_blocks=[base_block, FINE_BLOCK],
level_objects=[(), () if speed is None
else (circle(actuation(speed)),)]) as backend:
backend.initialize(*arrays)
for step in (0, 1, 2, 8):
if step == 1 and base_block is not None:
@@ -105,9 +105,10 @@ class BlockedFineCPU(unittest.TestCase):
return Mock(activate=Mock(side_effect=AssertionError("CUDA reached")), sm_arch="sm_70")
def test_rejections_before_cuda(self):
for kw in ({"fine_block": (7, 7)}, {"fine_block": (16,)}):
for kw in ({"level_blocks": [None, (7, 7)]}, {"level_blocks": [None, (16,)]}):
with self.subTest(kw=kw), self.assertRaises(ValueError):
FixedPatchStepper(config("periodic"), BOUNDS, context=self.forbidden_context(), **kw)
FixedPatchStepper(config("periodic"), context=self.forbidden_context(),
refinements=[BOUNDS], **kw)
def test_single_block_layout_allocates_no_empty_halo_tables(self):
# A 1x1 non-periodic fine layout has zero halo pairs; CUDA rejects mem_alloc(0).
@@ -125,8 +126,8 @@ class BlockedFineCPU(unittest.TestCase):
patch.object(compiler, "compile_kernel"), patch.object(compiler, "cleanup_build_dir"), \
patch.object(compiler, "load_module", return_value=module), \
patch.object(cuda, "mem_alloc", side_effect=alloc), patch.object(cuda, "memcpy_htod"):
backend = FixedPatchStepper(config("periodic"), BOUNDS, context=context,
fine_block=FINE_BLOCKS[0])
backend = FixedPatchStepper(config("periodic"), context=context, refinements=[BOUNDS],
level_blocks=[None, FINE_BLOCKS[0]])
self.assertEqual(len(backend.fine_layout.halo_dst), 0)
self.assertIsNone(backend._levels[1].halo_dst)
backend._gather_halo(backend._levels[1], Mock()) # no launch without pairs
@@ -177,8 +178,9 @@ class BlockedFineCPU(unittest.TestCase):
patch.object(cuda, "mem_alloc", side_effect=lambda size: Mock()), \
patch.object(cuda, "memcpy_htod"), \
patch.object(curved_links, "CurvedLinkSoA", side_effect=AssertionError("dense SoA built")):
backend = FixedPatchStepper(config("periodic"), BOUNDS, context=context,
fine_objects=objects, fine_block=FINE_BLOCKS[1])
backend = FixedPatchStepper(config("periodic"), context=context, refinements=[BOUNDS],
level_blocks=[None, FINE_BLOCKS[1]],
level_objects=[(), tuple(objects)])
backend.initialize(state(NY, NX, 1), state(*backend.patch.fine_shape, 2))
return backend, calls
@@ -289,8 +291,9 @@ class BlockedFineGPU(unittest.TestCase):
results = []
with CudaContext(0) as context:
for fine_block in (None, *FINE_BLOCKS):
with FixedPatchStepper(cfg, BOUNDS, context=context, coarse_boundary=mode,
fine_block=fine_block) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[BOUNDS],
level_blocks=[None, fine_block]) as backend:
backend.initialize(g_coarse, g_fine)
backend.step(STEPS)
self.assertEqual(backend.step_count, STEPS)
@@ -318,8 +321,9 @@ class BlockedFineGPU(unittest.TestCase):
results = []
with CudaContext(0) as context:
for fine_block in (None, *FINE_BLOCKS):
with FixedPatchStepper(cfg, BOUNDS, context=context, coarse_boundary=mode,
fine_objects=[circle()], fine_block=fine_block) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[BOUNDS], level_blocks=[None, fine_block],
level_objects=[(), (circle(),)]) as backend:
backend.initialize(g_coarse, g_fine)
backend.step(STEPS)
results.append((backend.download(), raw(backend)))
@@ -346,9 +350,9 @@ class BlockedFineGPU(unittest.TestCase):
results = []
with CudaContext(0) as context:
for fine_block in (None, *FINE_BLOCKS):
with FixedPatchStepper(cfg, BOUNDS, context=context, coarse_boundary=mode,
fine_objects=[circle(actuation(speed))],
fine_block=fine_block) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[BOUNDS], level_blocks=[None, fine_block],
level_objects=[(), (circle(actuation(speed)),)]) as backend:
backend.initialize(g_coarse, g_fine)
backend.step(STEPS)
results.append((backend.download(), raw(backend)))
@@ -94,7 +94,8 @@ class FaultGuardsCPU(unittest.TestCase):
("slot_reservoir", (12, 4, 20, 21))):
with self.subTest(mode=mode, bounds=bounds):
with self.assertRaises(ValueError):
implementation.FixedPatchStepper(cfg, bounds, context=context, coarse_boundary=mode)
implementation.FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[bounds])
context.activate.assert_not_called()
# Both exact channel limits and the original periodic two-cell margin
# must reach context activation, with no allocation or device required.
@@ -102,7 +103,8 @@ class FaultGuardsCPU(unittest.TestCase):
("periodic", (2, 2, 6, 5))):
with self.subTest(mode=mode, bounds=bounds):
with self.assertRaisesRegex(RuntimeError, "validated before CUDA"):
implementation.FixedPatchStepper(cfg, bounds, context=context, coarse_boundary=mode)
implementation.FixedPatchStepper(cfg, context=context, coarse_boundary=mode,
refinements=[bounds])
def fake_backend(self):
backend = PeriodicPatchStepper.__new__(PeriodicPatchStepper)
@@ -154,12 +156,12 @@ class FaultGuardsCPU(unittest.TestCase):
context = Mock(activate=nullcontext, sm_arch="sm_70")
with patch.object(implementation, "HALO", 2):
with self.assertRaisesRegex(RuntimeError, "four-cell"):
PeriodicPatchStepper(cfg, (4, 4, 12, 12), context=context)
PeriodicPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)])
p = implementation.FixedPatch2D.from_config(cfg, (4, 4, 12, 12))
with patch.object(cuda, "Stream"), patch.object(implementation.FixedPatch2D, "from_config",
return_value=replace(p, tau_coarse=.9)):
with self.assertRaisesRegex(RuntimeError, "omega/tau"):
PeriodicPatchStepper(cfg, p.bounds, context=context)
PeriodicPatchStepper(cfg, context=context, refinements=[p.bounds])
module = Mock()
module.get_global.return_value = (Mock(), 16)
with patch.object(cuda, "Stream"), patch.object(compiler, "generate_config"), \
@@ -167,7 +169,7 @@ class FaultGuardsCPU(unittest.TestCase):
patch.object(compiler, "load_module", return_value=module), \
patch.object(cuda, "memcpy_htod") as upload:
with self.assertRaisesRegex(RuntimeError, "layout mismatch"):
PeriodicPatchStepper(cfg, p.bounds, context=context)
PeriodicPatchStepper(cfg, context=context, refinements=[p.bounds])
upload.assert_not_called()
@@ -241,7 +243,7 @@ class NormalBuildGPU(unittest.TestCase):
max_error = max_bc_error = 0.0
with CudaContext(0) as context:
with patch.object(compiler, "generate_config", side_effect=record_config):
backend = implementation.FixedPatchStepper(cfg, bounds, context=context,
backend = implementation.FixedPatchStepper(cfg, context=context, refinements=[bounds],
coarse_boundary="slot_reservoir")
with backend:
self.assertEqual(len(generated), 2)
@@ -311,7 +313,8 @@ class NormalBuildGPU(unittest.TestCase):
for a in (initial.coarse_f, initial.fine_f)))
max_pop = max_drift = max_moment = 0.0
with CudaContext(0) as context, PeriodicPatchStepper(
LBMConfig(nx=16, ny=16, viscosity=.1), (4, 4, 12, 12), context=context) as backend:
LBMConfig(nx=16, ny=16, viscosity=.1), context=context,
refinements=[(4, 4, 12, 12)]) as backend:
backend.initialize(*postcollision(state, spec, maps))
before = totals(backend.download(), maps)
for step in range(1, 9):
@@ -342,7 +345,7 @@ class NormalBuildGPU(unittest.TestCase):
fine_mask[4:-4, 4:-4] = True
initial[1][:, ~fine_mask] = np.nan
expected = scalar_interval(*[a.astype(np.float64) for a in initial], bounds)
with CudaContext(0) as context, PeriodicPatchStepper(cfg, bounds, context=context) as backend:
with CudaContext(0) as context, PeriodicPatchStepper(cfg, context=context, refinements=[bounds]) as backend:
self.assertEqual(backend.patch.fine_shape, (14, 16))
p = backend.patch
allocated = 2 * 9 * (17*13+16*14) * 4 + sum(a.nbytes for a in
@@ -84,13 +84,13 @@ class CoupledGPU(unittest.TestCase):
cases = (("dense_unshifted", False, {}),
("blocked_shifted", True,
dict(base_block=(4, 4), fine_block=(6, 6))))
dict(level_blocks=[(4, 4), (6, 6)])))
with CudaContext(0) as context:
for label, shifted, storage in cases:
with self.subTest(storage=label):
cfg = LBMConfig(nx=16, ny=16, viscosity=.0192,
ddf_shifting=shifted)
with FixedPatchStepper(cfg, (4, 4, 12, 12), context=context,
with FixedPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)],
**storage) as backend:
self.assertAlmostEqual(backend.patch.tau_coarse, spec.tau_coarse)
self.assertAlmostEqual(backend.patch.tau_fine, spec.tau_fine)
@@ -148,7 +148,7 @@ class CoupledGPU(unittest.TestCase):
maps = ref.build_maps(spec)
with CudaContext(0) as context:
with PeriodicPatchStepper(LBMConfig(nx=16, ny=16, viscosity=.1),
(4, 4, 12, 12), context=context) as backend:
context=context, refinements=[(4, 4, 12, 12)]) as backend:
backend.initialize(*postcollision(ref.make_state(spec), spec, maps))
coarse, fine = backend._levels
@@ -187,7 +187,7 @@ class CoupledGPU(unittest.TestCase):
with CudaContext(0) as context:
with patch.object(compiler, "compile_kernel", audit_compile):
backend = PeriodicPatchStepper(cfg, (4, 4, 12, 12), context=context)
backend = PeriodicPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)])
with backend:
# Production allocations: four DDF arrays and four compact tables.
device_bytes = 2 * 9 * (16**2 + 24**2) * 4
@@ -270,7 +270,7 @@ class CoupledGPU(unittest.TestCase):
backend.download()
# Compile and execute the actual non-audited entry as well: one step.
with PeriodicPatchStepper(cfg, (4, 4, 12, 12), context=context) as plain:
with PeriodicPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)]) as plain:
rest = postcollision(ref.make_state(spec), spec, maps)
plain.initialize(*rest)
plain.step()
@@ -41,8 +41,9 @@ def mocked(objects, jets, fine_block):
patch.object(compiler, "load_module", return_value=module), \
patch.object(cuda, "mem_alloc", side_effect=lambda size: Mock()), \
patch.object(cuda, "memcpy_htod", htod):
backend = FixedPatchStepper(config("periodic"), BOUNDS, context=context,
fine_objects=objects, fine_block=fine_block, jets=jets)
backend = FixedPatchStepper(config("periodic"), context=context, refinements=[BOUNDS],
level_blocks=[None, fine_block], level_objects=[(), tuple(objects)],
level_jets=[(), tuple(jets)])
return backend, htod
@@ -149,8 +150,9 @@ class JetStepperCPU(unittest.TestCase):
def test_jet_overlapping_legacy_aperture_rejected_before_cuda(self):
with patch.object(cuda, "Stream", side_effect=AssertionError("CUDA touched")):
with self.assertRaisesRegex(ValueError, "SurfaceActuation"):
FixedPatchStepper(config("periodic"), BOUNDS, context=Mock(),
fine_objects=[circle(actuation(SPEED))], jets=[jet()])
FixedPatchStepper(config("periodic"), context=Mock(), refinements=[BOUNDS],
level_objects=[(), (circle(actuation(SPEED)),)],
level_jets=[(), (jet(),)])
def run(mode, fine_block, *, objects, jets, switch=None):
@@ -158,8 +160,9 @@ def run(mode, fine_block, *, objects, jets, switch=None):
cfg = config(mode)
patch_ = FixedPatch2D.from_config(cfg, BOUNDS)
with CudaContext(0) as context, FixedPatchStepper(
cfg, BOUNDS, context=context, coarse_boundary=mode, fine_objects=objects,
fine_block=fine_block, jets=jets) as backend:
cfg, context=context, coarse_boundary=mode, refinements=[BOUNDS],
level_blocks=[None, fine_block], level_objects=[(), tuple(objects)],
level_jets=[(), tuple(jets)]) as backend:
backend.initialize(state(NY, NX, 1), state(*patch_.fine_shape, 2))
if switch is None:
backend.step(STEPS)
@@ -42,8 +42,9 @@ class RefinementMRTGPU(unittest.TestCase):
collision="MRT", ddf_shifting=shifted)
reference = None
for block in ((32, 48), (8, 8), (4, 1)):
with FixedPatchStepper(cfg, (8, 6, 40, 26), context=context, coarse_boundary=boundary,
base_block=block, fine_block=(8, 8), fine_objects=(body,)) as backend:
with FixedPatchStepper(cfg, context=context, coarse_boundary=boundary,
refinements=[(8, 6, 40, 26)], level_blocks=[block, (8, 8)],
level_objects=[(), (body,)]) as backend:
backend.initialize(*[state(*level.shape[1:]) for level in backend._levels])
backend.step(8)
raw = backend._download_raw()
@@ -63,8 +64,8 @@ class RefinementMRTGPU(unittest.TestCase):
results = {}
layouts = ((None, None), *((base, FINE_BLOCK) for base in (None, *BASE_BLOCKS)))
for base_block, fine_block in layouts:
with FixedPatchStepper(cfg, BOUNDS, context=context, fine_block=fine_block,
base_block=base_block) as backend:
with FixedPatchStepper(cfg, context=context, refinements=[BOUNDS],
level_blocks=[base_block, fine_block]) as backend:
backend.initialize(*[fixture.state(*level.shape[1:], seed) for seed, level in enumerate(backend._levels)])
backend.step(4)
results[base_block, fine_block] = raw(backend)
@@ -75,9 +76,9 @@ class RefinementMRTGPU(unittest.TestCase):
def test_shifted_all_blocked_three_level_mrt_smoke(self):
"""Per-level shear rate omega_l = 1/(3*2**l*nu + 1/2); finite ACTIVE leaves."""
cfg = LBMConfig(nx=16, ny=16, viscosity=.1, collision="MRT", ddf_shifting=True)
kw = dict(fine_block=(4, 4), base_block=(4, 4), child=(8, 8, 16, 16), child_block=(4, 4))
kw = dict(refinements=[(4, 4, 12, 12), (8, 8, 16, 16)], level_blocks=[(4, 4)] * 3)
with CudaContext(0) as context:
with FixedPatchStepper(cfg, (4, 4, 12, 12), context=context, **kw) as backend:
with FixedPatchStepper(cfg, context=context, **kw) as backend:
self.assertEqual(len(backend._level_cfgs), 3)
for level, level_cfg in enumerate(backend._level_cfgs):
checked = replace(level_cfg)
@@ -26,6 +26,12 @@ FACE = ((4, 4, 8, 12), (12, 4, 16, 12))
S1 = ((4, 4, 12, 12), (20, 4, 28, 12))
def levels(bounds, blocks=(None, None), objects=(), jets=(), **rest):
"""Two-level constructor keywords: one L1 rectangle or set, per-level blocks, L1 bodies/Jets."""
return dict(refinements=[bounds], level_blocks=list(blocks),
level_objects=[(), tuple(objects)], level_jets=[(), tuple(jets)], **rest)
def level_fine(patch_set, windows):
"""Level-wide fine array with each patch window filled from `windows`."""
out = np.full((9, *patch_set.fine_shape), np.nan, np.float32)
@@ -45,21 +51,21 @@ class MultipatchCPU(unittest.TestCase):
def test_rejections_before_cuda(self):
cfg = config("periodic")
context = object() # never activated: every case must fail first
for kw in (dict(), dict(fine_block=(8, 4)), dict(base_block=(4, 4)),
dict(fine_block=(8, 4), base_block=(4, 4), fine_objects=[circle()]),
dict(fine_block=(5, 4), base_block=(4, 4))):
for kw in (dict(), dict(blocks=(None, (8, 4))), dict(blocks=((4, 4), None)),
dict(blocks=((4, 4), (8, 4)), objects=[circle()]),
dict(blocks=((4, 4), (5, 4)))):
with self.subTest(kw=kw), self.assertRaises(ValueError):
FixedPatchStepper(cfg, ((4, 4, 8, 12), (14, 4, 18, 12)), context=context, **kw)
FixedPatchStepper(cfg, context=context, **levels(((4, 4, 8, 12), (14, 4, 18, 12)), **kw))
with self.assertRaisesRegex(ValueError, "margin"):
FixedPatchStepper(cfg, ((1, 4, 8, 12), (10, 4, 14, 12)), context=context,
fine_block=(8, 4), base_block=(4, 4))
FixedPatchStepper(cfg, context=context,
**levels(((1, 4, 8, 12), (10, 4, 14, 12)), blocks=((4, 4), (8, 4))))
from CelerisLab.body.jets import Jet, Parallel, Tube
cfg32 = LBMConfig(nx=32, ny=16, viscosity=.1)
span = Jet(Tube((5., 8.), 0., 22., .6), Parallel(0.), amplitude=0.)
with self.assertRaisesRegex(ValueError, "one jet"): # same obj_id cannot bypass it
FixedPatchStepper(cfg32, ((2, 2, 14, 14), (18, 2, 30, 14)), context=context,
fine_block=(4, 4), base_block=(4, 4),
fine_objects=[body(8., 8.), body(24., 8.)], jets=[span])
FixedPatchStepper(cfg32, context=context,
**levels(((2, 2, 14, 14), (18, 2, 30, 14)), blocks=((4, 4), (4, 4)),
objects=[body(8., 8.), body(24., 8.)], jets=[span]))
class MultipatchGPU(unittest.TestCase):
@@ -67,7 +73,7 @@ class MultipatchGPU(unittest.TestCase):
def run_steps(self, context, cfg, bounds, coarse, fine, **kw):
out = []
with FixedPatchStepper(cfg, bounds, context=context, **kw) as backend:
with FixedPatchStepper(cfg, context=context, **levels(bounds, **kw)) as backend:
backend.initialize(coarse, fine)
for step in self.STEPS:
backend.step(step - backend.step_count)
@@ -82,7 +88,7 @@ class MultipatchGPU(unittest.TestCase):
p = FixedPatch2D.from_config(cfg, BOUNDS)
s = FixedPatchSet2D.from_config(cfg, (BOUNDS,))
coarse, fine = state(NY, NX, 1), state(*p.fine_shape, 2)
kw = dict(coarse_boundary=mode, fine_block=(8, 4), base_block=(3, 4))
kw = dict(coarse_boundary=mode, blocks=((3, 4), (8, 4)))
old = self.run_steps(context, cfg, BOUNDS, coarse, fine, **kw)
new = self.run_steps(context, cfg, (BOUNDS,), coarse, level_fine(s, [fine]), **kw)
for step, a, b in zip(self.STEPS, old, new):
@@ -99,7 +105,7 @@ class MultipatchGPU(unittest.TestCase):
one = FixedPatchSet2D.from_config(cfg16, (S1[0],))
two = FixedPatchSet2D.from_config(cfg32, S1)
coarse, fine = state(16, 16, 3), state(*one.patches[0].fine_shape, 4)
kw = dict(fine_block=(4, 4), base_block=(4, 4))
kw = dict(blocks=((4, 4), (4, 4)))
with CudaContext(0) as context:
ref = self.run_steps(context, cfg16, (S1[0],), coarse, level_fine(one, [fine]), **kw)
tiled = np.concatenate([coarse, coarse], axis=2)
@@ -118,7 +124,7 @@ class MultipatchGPU(unittest.TestCase):
s = FixedPatchSet2D.from_config(cfg, FACE)
coarse = state(16, 24, 5)
fine = level_fine(s, [state(*p.fine_shape, 6 + i) for i, p in enumerate(s.patches)])
kw = dict(fine_block=(4, 4), base_block=(4, 4))
kw = dict(blocks=((4, 4), (4, 4)))
leaf = s.fine_kind == ACTIVE
support = s.fine_layout((4, 4)).owned
outer = support & (s.fine_kind == OUTER_GHOST)
@@ -164,7 +170,7 @@ class MultipatchBodiesGPU(unittest.TestCase):
def run_steps(self, context, cfg, bounds, coarse, fine, switch=None, **kw):
out = []
with FixedPatchStepper(cfg, bounds, context=context, **kw) as backend:
with FixedPatchStepper(cfg, context=context, **levels(bounds, **kw)) as backend:
backend.initialize(coarse, fine)
for step in self.STEPS:
if switch is not None and step == 8:
@@ -182,8 +188,7 @@ class MultipatchBodiesGPU(unittest.TestCase):
p = FixedPatch2D.from_config(cfg, BOUNDS)
s = FixedPatchSet2D.from_config(cfg, (BOUNDS,))
coarse, fine = state(NY, NX, 1), state(*p.fine_shape, 2)
kw = dict(coarse_boundary=mode, fine_block=(8, 4), base_block=(3, 4),
fine_objects=[circle()])
kw = dict(coarse_boundary=mode, blocks=((3, 4), (8, 4)), objects=[circle()])
old = self.run_steps(context, cfg, BOUNDS, coarse, fine, **kw)
new = self.run_steps(context, cfg, (BOUNDS,), coarse, level_fine(s, [fine]), **kw)
for step, a, b in zip(self.STEPS, old, new):
@@ -202,17 +207,17 @@ class MultipatchBodiesGPU(unittest.TestCase):
two = FixedPatchSet2D.from_config(cfg32, self.PATCHES)
coarse, fine = state(16, 16, 3), state(*one.patches[0].fine_shape, 4)
jet = lambda x: Jet(Tube((x, 8.), 0., 6., .6), Parallel(0.), amplitude=.002)
kw = dict(fine_block=(4, 4), base_block=(4, 4))
kw = dict(blocks=((4, 4), (4, 4)))
with CudaContext(0) as context:
for jets in ((), (jet(5.),)):
ref = self.run_steps(context, cfg16, self.PATCHES[:1], coarse, level_fine(one, [fine]),
switch=[(0, -.001)] if jets else None,
fine_objects=[body(8., 8.)], jets=jets, **kw)
objects=[body(8., 8.)], jets=jets, **kw)
both = jets and (jet(5.), jet(21.))
got = self.run_steps(context, cfg32, self.PATCHES, np.concatenate([coarse, coarse], 2),
level_fine(two, [fine, fine]),
switch=[(0, -.001), (1, -.001)] if jets else None,
fine_objects=[body(8., 8., 5), body(24., 8., 5)],
objects=[body(8., 8., 5), body(24., 8., 5)],
jets=both or (), **kw)
for step, a, b in zip(self.STEPS, ref, got):
label = f"jets={len(jets)} step {step}"
@@ -27,7 +27,16 @@ from test_refinement_blocked_base import raw
CFG = LBMConfig(nx=16, ny=16, viscosity=.1)
L1, L2 = (4, 4, 12, 12), (8, 8, 16, 16)
KW = dict(fine_block=(4, 4), base_block=(4, 4), child=L2, child_block=(4, 4))
B = (4, 4)
def levels(l1=L1, l2=L2, blocks=(B, B, B), objects=((), (), ()), jets=((), (), ())):
"""Level-list keywords for the scalar L1 patch with one scalar L2 child."""
return dict(refinements=[l1, l2], level_blocks=list(blocks),
level_objects=list(objects), level_jets=list(jets))
KW = levels()
def child_body(spec=None):
@@ -52,26 +61,29 @@ def mocked_child():
patch.object(compiler, "load_module", return_value=module), \
patch.object(cuda, "mem_alloc", side_effect=lambda size: Mock()), \
patch.object(cuda, "memcpy_htod"):
return FixedPatchStepper(CFG, L1, context=context, child_objects=[child_body()],
child_jets=[child_jet()], **KW)
return FixedPatchStepper(CFG, context=context,
**levels(objects=((), (), (child_body(),)), jets=((), (), (child_jet(),))))
class ThreeLevelCPU(unittest.TestCase):
def test_rejections_before_cuda(self):
context = object() # never activated: every case must fail first
bad = (dict(KW, child_block=None), dict(KW, fine_block=None),
dict(KW, fine_objects=[circle()]), dict(KW, child=(4, 8, 16, 16)),
dict(KW, child_block=(5, 4)))
bad = (levels(blocks=(B, B, None)), levels(blocks=(B, None, B)),
levels(objects=((), (circle(),), ())), levels(l2=(4, 8, 16, 16)),
levels(blocks=(B, B, (5, 4))))
for kw in bad:
with self.subTest(kw=kw), self.assertRaises(ValueError):
FixedPatchStepper(CFG, L1, context=context, **kw)
for kw in (dict(child_objects=[child_body()]), dict(child_jets=[child_jet()]),
dict(KW, child_jets=[child_jet()])):
FixedPatchStepper(CFG, context=context, **kw)
# L2 bodies/Jets without an L2 level (one entry too many), and an L2 Jet without a body.
two_level = dict(refinements=[L1], level_blocks=[None, None])
for kw in (dict(two_level, level_objects=[(), (), (child_body(),)]),
dict(two_level, level_jets=[(), (), (child_jet(),)]),
levels(jets=((), (), (child_jet(),)))):
with self.assertRaises(ValueError):
FixedPatchStepper(CFG, L1, context=context, **kw)
with self.assertRaises(ValueError): # several L1 patches + child: not yet
FixedPatchStepper(LBMConfig(nx=32, ny=16, viscosity=.1), ((4, 4, 12, 12), (20, 4, 28, 12)),
context=context, **KW)
FixedPatchStepper(CFG, context=context, **kw)
with self.assertRaises(ValueError): # several L1 patches + scalar child: not yet
FixedPatchStepper(LBMConfig(nx=32, ny=16, viscosity=.1), context=context,
**levels(l1=((4, 4, 12, 12), (20, 4, 28, 12))))
def test_four_halo_bfl_steps_and_fresh_middle_feeds(self):
@@ -156,7 +168,7 @@ class ThreeLevelGPU(unittest.TestCase):
def test_equilibrium_is_kept_on_every_level(self):
from CelerisLab.cuda.context import CudaContext
feq = np.asarray(compute_feq_d2q9(1.0, .03, -.02), np.float32)[:, None, None]
with CudaContext(0) as ctx, FixedPatchStepper(CFG, L1, context=ctx, **KW) as run:
with CudaContext(0) as ctx, FixedPatchStepper(CFG, context=ctx, **KW) as run:
shapes = ((16, 16), run.patch.fine_shape, run.child_patch.fine_shape)
run.initialize(*(np.broadcast_to(feq, (9, *s)).copy() for s in shapes))
run.step(8)
@@ -169,7 +181,7 @@ class ThreeLevelGPU(unittest.TestCase):
def test_nonequilibrium_masks_and_fresh_middle_feeds(self):
from CelerisLab.cuda.context import CudaContext
with CudaContext(0) as ctx, FixedPatchStepper(CFG, L1, context=ctx, **KW) as run:
with CudaContext(0) as ctx, FixedPatchStepper(CFG, context=ctx, **KW) as run:
p1, p2 = run.patch, run.child_patch
run.initialize(state(16, 16, 1), state(*p1.fine_shape, 2), state(*p2.fine_shape, 3))
f = run._levels[1]
@@ -205,8 +217,9 @@ class ThreeLevelBodyGPU(unittest.TestCase):
CFG, velocity=0., wall_normal_velocity=.001, wall_slot_width=4,
y_wall_bc="on_site_velocity", top_wall_bc="reservoir")
with CudaContext(0) as ctx, FixedPatchStepper(
cfg, L1, context=ctx, coarse_boundary=mode, child_objects=[child_body(spec)],
child_jets=jets, **dict(KW, child_block=block)) as run:
cfg, context=ctx, coarse_boundary=mode,
**levels(blocks=(B, B, block), objects=((), (), (child_body(spec),)),
jets=((), (), tuple(jets)))) as run:
run.initialize(*(state(*level.shape[1:], i + 1) for i, level in enumerate(run._levels)))
snapshots = []
for stop in (0, 1, 2, 4, 8):
@@ -78,10 +78,10 @@ class ThreeLevelMultipatchGPU(unittest.TestCase):
def run_case(self, ctx, cfg=CFG, parents=PARENTS, children=CHILDREN,
mode='periodic', blocks=((4, 4), (4, 4)), objects=(), drive=(), switch=False,
scalar_child=None):
kwargs = dict(child_rectangles_in_parent=children) if scalar_child is None else dict(child=scalar_child)
with FixedPatchStepper(config(cfg, mode), parents, context=ctx, coarse_boundary=mode,
base_block=(4, 4), fine_block=blocks[0], child_block=blocks[1],
child_objects=objects, child_jets=drive, **kwargs) as run:
child = children if scalar_child is None else scalar_child
with FixedPatchStepper(config(cfg, mode), context=ctx, coarse_boundary=mode,
refinements=[parents, child], level_blocks=[(4, 4), blocks[0], blocks[1]],
level_objects=[(), (), objects], level_jets=[(), (), drive]) as run:
frame = frames(run)
run.initialize(*(physical_state(level.shape[1:], *f)
for level, f in zip(run._levels, frame)))
@@ -26,9 +26,8 @@ FP64_ATOL = 3e-6
def make(ctx):
return FixedPatchStepper(CFG, PARENTS, context=ctx, base_block=(4, 4),
fine_block=(4, 4), child_block=(4, 4),
child_rectangles_in_parent=A)
return FixedPatchStepper(CFG, context=ctx, refinements=[PARENTS, A],
level_blocks=[(4, 4)] * 3)
def initial(run, uniform):
@@ -136,8 +136,8 @@ class ShiftedJetsGPU(unittest.TestCase):
Parallel(0.),amplitude=.013)
for rho in (.99,1.01):
with self.subTest(rho=rho),CudaContext(0) as ctx:
with FixedPatchStepper(cfg,(8,6,24,18),context=ctx,base_block=(4,4),fine_block=(4,4),
fine_objects=(obj,),jets=(whole,),observe_raw_exchange=True) as run:
with FixedPatchStepper(cfg,context=ctx,refinements=[(8,6,24,18)],level_blocks=[(4,4)]*2,
level_objects=[(),(obj,)],level_jets=[(),(whole,)],observe_raw_exchange=True) as run:
run._initialize_raw(initial(run,rho,u=.013,v=0.))
expected=run._download_raw();ledger=Ledger(self,run,1)
self.assertTrue(run._level_boundaries[1].tables['aperture'].all())
@@ -153,7 +153,8 @@ class ShiftedJetsGPU(unittest.TestCase):
obj=circle((16.,12.),2.15);j=drive((16.,12.),2.15)
for rho in (.99,1.01):
with self.subTest(rho=rho),CudaContext(0) as ctx,ExitStack() as owned:
settings=dict(bounds=(8,6,24,18),context=ctx,base_block=(4,4),fine_block=(4,4),fine_objects=(obj,),jets=(j,))
settings=dict(context=ctx,refinements=[(8,6,24,18)],level_blocks=[(4,4)]*2,
level_objects=[(),(obj,)],level_jets=[(),(j,)])
measured=owned.enter_context(FixedPatchStepper(cfg,observe_raw_exchange=True,**settings))
plain=owned.enter_context(FixedPatchStepper(cfg,**settings))
measured._initialize_raw(initial(measured,rho));plain._initialize_raw(initial(plain,rho))
@@ -186,8 +187,9 @@ class ShiftedJetsGPU(unittest.TestCase):
parents=((4,4,20,28),(24,4,44,28));a=((16,16,32,40),(56,16,72,40));b=(a[0],(60,16,76,40))
obj=circle((12.,14.),.8);j=drive((12.,14.),.8)
with CudaContext(0) as ctx,ExitStack() as owned:
old=owned.enter_context(FixedPatchStepper(cfg,parents,context=ctx,base_block=(4,4),fine_block=(4,4),
child_block=(4,4),child_rectangles_in_parent=a,child_objects=(obj,),child_jets=(j,),observe_raw_exchange=True))
old=owned.enter_context(FixedPatchStepper(cfg,context=ctx,refinements=[parents,a],
level_blocks=[(4,4)]*3,level_objects=[(),(),(obj,)],level_jets=[(),(),(j,)],
observe_raw_exchange=True))
old._initialize_raw(initial(old,.99));ledger=Ledger(self,old,2)
for speed in (.002,-.001):old.set_jet_amplitude(0,speed,level=2);ledger.step()
before=old._download_raw();count=old.step_count
@@ -22,9 +22,8 @@ B=((20,16,36,40),)
def make(ctx):
return FixedPatchStepper(CFG,PARENTS,context=ctx,base_block=(4,4),
fine_block=(4,4),child_block=(4,4),
child_rectangles_in_parent=A)
return FixedPatchStepper(CFG,context=ctx,refinements=[PARENTS,A],
level_blocks=[(4,4)]*3)
def initial(run):
+2 -2
View File
@@ -126,8 +126,8 @@ def run_case(a, h, case):
cfg, pair, body, host = rr.host_geometry(a.source, h, a.collision)
context = stack.enter_context(CudaContext(0))
backend = stack.enter_context(FixedPatchStepper(
cfg, rr.PLACEMENTS["wide"], context=context, coarse_boundary="channel",
base_block=rr.BASE_BLOCK, fine_block=rr.FINE_BLOCK, fine_objects=(body,),
cfg, context=context, coarse_boundary="channel", refinements=[rr.PLACEMENTS["wide"]],
level_blocks=[rr.BASE_BLOCK, rr.FINE_BLOCK], level_objects=[(), (body,)],
observe_raw_exchange=True, observe_window_exchange=True))
from CelerisLab.lbm.channel import channel_masks
from CelerisLab.lbm.refinement import ACTIVE
+4 -5
View File
@@ -279,12 +279,11 @@ def run_case(c, ctx, output, provenance, continuation=False):
PeriodicBlockStepper(fine_cfg, layout, context=ctx, objects=[fobj]) as uniform,
FixedPatchStepper(
cfg,
c["bounds"],
context=ctx,
base_block=(8, 8),
fine_block=(8, 8),
fine_objects=[obj],
jets=jets,
refinements=[c["bounds"]],
level_blocks=[(8, 8), (8, 8)],
level_objects=[(), (obj,)],
level_jets=[(), tuple(jets)],
observe_raw_exchange=True,
) as amr,
):
+4 -5
View File
@@ -383,13 +383,12 @@ def run_case(a, h, periodic, case):
backend = stack.enter_context(
FixedPatchStepper(
cfg,
PLACEMENTS[a.placement],
context=context,
coarse_boundary="channel",
base_block=BASE_BLOCK,
fine_block=FINE_BLOCK,
fine_objects=(body,),
jets=public_jets((80.5, 81.0), 20.0, a.jet_scenario),
refinements=[PLACEMENTS[a.placement]],
level_blocks=[BASE_BLOCK, FINE_BLOCK],
level_objects=[(), (body,)],
level_jets=[(), public_jets((80.5, 81.0), 20.0, a.jet_scenario)],
observe_raw_exchange=True,
**({"observe_window_exchange": True} if a.observation_mode == "window" else {}),
)
+3 -2
View File
@@ -173,8 +173,9 @@ def time_refined(a, h, sched):
cfg, pair, body, host = r.host_geometry(a.source, h, "MRT")
context = stack.enter_context(CudaContext(0))
backend = stack.enter_context(FixedPatchStepper(
cfg, r.BOUNDS, context=context, coarse_boundary="channel", base_block=r.BASE_BLOCK,
fine_block=r.FINE_BLOCK, fine_objects=(body,), observe_raw_exchange=True))
cfg, context=context, coarse_boundary="channel", refinements=[r.BOUNDS],
level_blocks=[r.BASE_BLOCK, r.FINE_BLOCK], level_objects=[(), (body,)],
observe_raw_exchange=True))
masks = [channel_masks(pair.coarse_shape)[0] & ~pair.covered, # as run_refined_dfg
(pair.fine_kind == ACTIVE) & ((host.flags.reshape(pair.fine_shape) & SOLID) == 0)]
fresh = (r.initial_state(pair.coarse_shape, (0.0, 0.0), 1.0, root=True),
+2 -2
View File
@@ -130,8 +130,8 @@ def main():
cfg, pair, body, host = r.host_geometry(a.source, h, "MRT")
context = stack.enter_context(CudaContext(0))
backend = stack.enter_context(FixedPatchStepper(
cfg, r.BOUNDS, context=context, coarse_boundary="channel", base_block=r.BASE_BLOCK,
fine_block=r.FINE_BLOCK, fine_objects=(body,), observe_raw_exchange=True,
cfg, context=context, coarse_boundary="channel", refinements=[r.BOUNDS],
level_blocks=[r.BASE_BLOCK, r.FINE_BLOCK], level_objects=[(), (body,)], observe_raw_exchange=True,
observe_window_exchange=True))
masks = [channel_masks(pair.coarse_shape)[0] & ~pair.covered,
(pair.fine_kind == ACTIVE) & ((host.flags.reshape(pair.fine_shape) & SOLID) == 0)]