Compare commits

..
Author SHA1 Message Date
ClaudeandClaude Opus 5.5 6d046f91c1 docs: sync src_coherence with the landed slices
Stage history and verification, historical label on the migration map,
host-preparation/device order, pair accessors, per-level state split into
fixed plan structure / mutable boundary data / runtime _Level, the two
block_collide families, and the deferred items with their evidence needs.
Doc-only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 12:13:06 +00:00
ClaudeandClaude Opus 5.5 3361e97d5f docs: src_coherence migration map for the combined _boundaries record
Final private names (_device_boundaries/_level_jets/_exchange_*/_level_blocks
-> _boundaries record fields / layout.block), public boundaries = host view,
and the plan-versus-device (_Level) lifecycle. Doc-only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 12:09:26 +00:00
ClaudeandClaude Opus 5.5 442979b4db FixedPatchStepper: exchange metadata on the boundary record; _boundaries name
- _exchange_spacing/_exchange_centers (two more per-level lists, only read
  for levels that have a boundary) become _BlockedBoundary.spacing/centers.
- _level_blocks dropped: remap_finest reads each layout's own .block.
- _device_boundaries -> _boundaries: the record holds host geometry, tables,
  Jets and device handles, so 'device' was misleading (the old host-only
  _boundaries list no longer exists). Tests migrated mechanically.
Remaining per-level lists (_layouts, _level_cfgs, _leaf, _state, _boundaries)
are the immutable construction plan, built before the memory preflight and
kept after close; _levels holds device handles only and is cleared on close.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 11:41:13 +00:00
ClaudeandClaude Opus 5.5 c5a3cdb757 test_refinement_mrt: check force reads in executable text, not comments
The block_collision.cuh header comment names the absent d_params.fx/fy
read, which tripped the source-text assertion. Strip // comments first
and also require the real d_params.omega read. Test-only; no kernel source
change.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 11:27:35 +00:00
ClaudeandClaude Opus 5.5 a4fc46d849 block_collide: one unforced uniform-equivalent branch; drop collide_mrt_unforced
The blocked collision had three arms (MRT wrapper, channel SRT, level SRT)
and a one-caller MRT wrapper in collision_mrt.cuh. Now block_collide computes
rho/u once and has two compile-time families:
- uniform-equivalent (MRT any level, SRT channel root): collide_dispatch
  arithmetic with zero force, no LES: one clamped omega, collide_mrt_no_force
  or collide_srt with its zero Guo term (unchanged);
- level SRT (BLOCK_ROOT_BC 0/1): unclamped, collide_srt_no_force, shifted
  momentum feq exception kept visible.
No runtime d_params.fx/fy is read. Uniform kernels unchanged (only the
unused inline wrapper is removed from collision_mrt.cuh). Source-text test
updated to the new single call sites.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-28 10:55:02 +00:00
22 changed files with 165 additions and 130 deletions
+49 -13
View File
@@ -1,8 +1,11 @@
# src coherence stage (from main 914906e)
Owner direction 2026-10-28: unified style over historical compatibility, inside `src`.
No numerical change is intended; kernel text was moved verbatim (PTX/bit parity is
pending a GPU run).
No numerical change is intended. Stages after 5f69eba: 65339ab (one boundary/Jet record,
host preparation helpers, pair accessors), c5a3cdb (`block_collide` families), 442979b
(exchange metadata on the boundary record). Each was checked by parent/candidate
config.h/PTX comparison and host/GPU equality runs; c5a3cdb changed the PTX of the two
channel-root SRT modules only (omega load/clamp moved, outputs exact).
## Conventions (defined once, in `lbm.hierarchy.Hierarchy`)
@@ -28,20 +31,25 @@ pending a GPU run).
## FixedPatchStepper migration map
Historical: left column = names at or after 914906e that no longer exist; right column =
the current name.
| removed / changed | replacement |
|---|---|
| `patch`, `child_patch` | `hierarchy.pairs[0]`, `hierarchy.pairs[1]` |
| `base_layout`, `fine_layout`, `child_layout` (stepper properties) | `layouts[0]`, `layouts[1]`, `layouts[2]` (`FixedPatchSet2D.fine_layout(block)` is unchanged) |
| `fine_boundary`, `child_boundary` | `boundaries[1]`, `boundaries[2]` (host geometry; `None` = no bodies/Jets) |
| `_jets` | `jets[1]` |
| `_boundaries[l]` (dict from L1, level keys) | `_device_boundaries[l]` (list from L0, `None` at L0 unless single-level static circles) |
| `_level_boundaries` (private, b4eacf7) | `_device_boundaries` |
| `_fine_tables`, `_actuated`, `_jet_cache` | `_device_boundaries[1].tables` / `.actuated` / `.jet_cache` |
| `_boundaries[l]` (dict from L1, level keys), `_level_boundaries` (b4eacf7), `_device_boundaries` (65339ab) | `_boundaries[l]`: the one combined per-level `_BlockedBoundary` record (list from L0, `None` at L0 unless single-level static circles): host geometry `.host`, `.tables`, `.jets`, `.jet_cache`, exchange `.spacing`/`.centers`, device handles |
| `_fine_tables`, `_actuated`, `_jet_cache` | `_boundaries[1].tables` / `.actuated` / `.jet_cache` |
| `_level_jets` (private) | `_boundaries[l].jets` (public `jets[l]` is derived) |
| `_exchange_spacing[l]`, `_exchange_centers[l]` (private) | `_boundaries[l].spacing` (h_l, root units) / `.centers` (root edge frame); only levels with a boundary report exchange |
| `_level_blocks` (private) | `_layouts[l].block` |
| `_remap_cfg` (private) | `_root_cfg` |
| `refinement_tables.blocked_base_tables` (also via `refinement_stepper`) | `refinement_tables.blocked_coarse_tables` (import from `refinement_tables`; the stepper no longer re-exports) |
| `_fine_leaf`, `_child_leaf` | `leaf[1]`, `leaf[2]` (copies) |
| `_middle_kind` | `hierarchy.kinds[1]` |
| `_boundary_hosts` (private) | `_boundaries` (private list behind `boundaries`) |
| `_boundary_hosts` (private) | public `boundaries[l]` = `_boundaries[l].host` (host geometry only; the private record is the combined one) |
| `remap_child(rects, child_objects=...)` | `remap_finest(rects, objects=...)` |
| `set_jet_amplitude(i, a)` (implicit level 1) | `set_jet_amplitude(i, a, level=l)` (required) |
| `read_raw_exchange(b)` / `read_exchange_window(b)` (implicit level 1) | `read_raw_exchange(b, level=l)` / `read_exchange_window(b, level=l)` (required) |
@@ -56,13 +64,35 @@ the replaced level relative to its parent).
## Responsibilities
- Construction: `_check_levels` (inputs), `_kernel_configs` (per-level kernel configs),
`Hierarchy.edge_frame` (geometry frames), `_build_level` / `_build_pair` / `_upload_boundary`
(device). Observation reads share `_observed_boundary` + `_summarize`.
- Construction: `_check_levels` (inputs); host preparation, numpy only, before any device
work: `_level_layouts`, `_kernel_configs` (per-level kernel configs), `_host_boundaries`,
`_prepare_blocked_boundary`, `_fluid_leaves`; `Hierarchy.edge_frame` (geometry frames);
then `_planned_bytes()` (memory preflight) and `_build_level` / `_build_pair` /
`_upload_boundary` (device). A construction failure closes the partial instance.
Observation reads share `_observed_boundary` + `_summarize`.
- Pair accessors: `rectangles` (a scalar FixedPatch2D gives `(bounds,)`) and
`fine_level_offset` (`(0, 0)` for a set) on both pair types; the types and their
storage frames stay distinct.
- Per-level state, by lifecycle:
- Fixed plan structure (built once in `__init__`, read by `_planned_bytes()` before any
allocation, readable after `close()` through the public views): `_layouts`,
`_level_cfgs`, `_leaf`, `_state`, and in each `_boundaries[l]` record the host geometry
`.host`, link tables, `.jet_cache`, `.spacing`, `.centers`.
- Mutable boundary data in the same record: `.jets` and `.tables["surface_velocity"]`,
changed only by `set_jet_amplitude` (host prepare -> device write -> commit; an uncertain
write fails the instance); device handles (`.flags`, `.device`, `.kernel`, ledgers) are
set by `_upload_boundary` and cleared by `close()`.
- Runtime: `_levels` (`_Level`) holds compiled modules, kernels, the two population
buffers (swapped per step) and transfer/kind handles. It is created only after the
memory preflight and cleared on `close()`; `.layout`/`.boundary` reference the plan for
`_hop`. Folding the plan into `_Level` would mix both lifecycles in one record.
- Blocked kernels: `operators/block_collision.cuh block_collide` is the one collision extension
point of `BlockStep` and `BlockChannelStep`. MRT shares one path; SRT keeps its two explicit
arithmetic variants (level: unclamped, no-force, shifted direct-momentum feq; channel: uniform
OneStep arithmetic with omega clamp), selected by `BLOCK_ROOT_BC`.
point of `BlockStep` and `BlockChannelStep`. It computes rho/u once and reads no runtime
force (`d_params.fx/fy`). Two compile-time families: uniform-equivalent (MRT on any level,
SRT on the channel root `BLOCK_ROOT_BC == 2`): clamped omega, `collide_mrt_no_force` or
`collide_srt` with its zero Guo term; and level SRT (`BLOCK_ROOT_BC` 0/1): unclamped omega,
`collide_srt_no_force`, with the shifted direct-momentum feq under `USE_DDF_SHIFTING`.
The uniform path keeps `collide_dispatch` (runtime force, LES, TRT).
- `explosion="projected_linear"` stays an optional capability; `"constant"` stays default.
## Not done / open
@@ -72,6 +102,12 @@ the replaced level relative to its parent).
`FixedPatchStepper._kernel_configs` (reset only on level 1 of a slot_reservoir root;
generates each level's config.h). Folding them needs a config.h/PTX diff gate.
- No single-`Simulation` facade over FixedPatchStepper: the two own different state models
(uniform node frame, body sync, forcing, checkpoint vs. immutable blocked hierarchy, edge
(uniform node frame, body sync, forcing, checkpoint vs. fixed blocked hierarchy, edge
frame, remap-as-new-instance). A facade would only forward calls; no concrete need found.
- Deferred, each needing its own evidence first: sharing the channel-root cell/donor selector
with OneStep (the closure functions are already shared; a shared selector changes uniform
signatures and needs a uniform PTX gate); merging FixedPatch2D into FixedPatchSet2D
(different origin, ghost band, layout size and addressing; needs a geometry/parity table and
would need tiles that are not grid-aligned); one `clamp_omega` helper for uniform and
blocked paths (it lives in the hash-pinned LES header `turbulence_smag.cuh`).
- Tests, examples and validation callers are migrated separately (pinned runners stay frozen).
@@ -1,10 +1,12 @@
// CelerisLab – operators/block_collision.cuh
// The one collision extension point of the blocked steps (BlockStep, BlockChannelStep).
// Compile-time selection: MRT shares one path; SRT keeps two explicit arithmetic
// variants, unchanged from their original kernels:
// level (BLOCK_ROOT_BC 0/1): unclamped omega, no forcing term; shifted storage
// builds feq from the retained drho and direct momentum (numerical variant);
// channel (BLOCK_ROOT_BC 2): matches uniform OneStep SRT arithmetic (omega clamp,
// zero Guo forcing term, compute_feq).
// Unforced by construction: FixedPatchStepper supplies no forcing, so no runtime
// d_params.fx/fy branch is read. Two compile-time arithmetic families:
// uniform-equivalent (MRT on every level; SRT on the channel root, BLOCK_ROOT_BC 2):
// collide_dispatch arithmetic with zero force and no LES: clamped omega;
// SRT keeps the zero Guo forcing term of collide_srt;
// level SRT (BLOCK_ROOT_BC 0/1): unclamped omega, no forcing term; shifted storage
// builds feq from the retained drho and direct momentum (numerical variant).
#ifndef CELERIS_OPERATORS_BLOCK_COLLISION_CUH
#define CELERIS_OPERATORS_BLOCK_COLLISION_CUH
@@ -16,17 +18,19 @@
#endif
__device__ __forceinline__ void block_collide(float* f) {
#if COLLISION_MODEL == 2
collide_mrt_unforced(f);
#elif BLOCK_ROOT_BC == 2
float rho, ux, uy, feq[NQ], forcing[NQ] = {};
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
compute_feq(rho, ux, uy, feq);
#if COLLISION_MODEL == 2 || BLOCK_ROOT_BC == 2
const float omega = fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega));
#if COLLISION_MODEL == 2
collide_mrt_no_force(f, rho, ux, uy, omega);
#else
float feq[NQ], forcing[NQ] = {};
compute_feq(rho, ux, uy, feq);
collide_srt(f, feq, forcing, omega);
#endif
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
float feq[NQ];
#if USE_DDF_SHIFTING
// Preserve the density perturbation before adding the unit reference.
const float drho = f[0]+f[1]+f[2]+f[3]+f[4]+f[5]+f[6]+f[7]+f[8];
@@ -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,
+23 -24
View File
@@ -317,26 +317,25 @@ class FixedPatchStepper:
self._window_valid = False
self._window_start = None
self._exchange_valid = False
self._exchange_spacing = tuple(self.hierarchy.edge_frame(l)[1] for l in range(depth))
self._exchange_centers = tuple(tuple(tuple(o.geometry.center) for o in objects)
for objects in level_objects)
root = None
if depth == 1 and level_objects[0]:
tables = static_circle_tables(self._layouts[0], level_objects[0])
root = _BlockedBoundary(None, tables, (), False)
# The one boundary/Jet record per level; host geometry and Jets are views of it.
self._device_boundaries = [root] + [
self._boundaries = [root] + [
_prepare_blocked_boundary(hosts[level], self._layouts[level],
self._level_cfgs[level], level_jets[level],
_shifted_blocked=self._ddf_shifting)
for level in range(1, depth)]
self._leaf = _fluid_leaves(self.hierarchy, self._device_boundaries)
for level, b in enumerate(self._boundaries):
if b is not None: # exchange scaling (dt_l = h_l) and root-edge-frame moment centers
b.spacing = self.hierarchy.edge_frame(level)[1]
b.centers = tuple(tuple(o.geometry.center) for o in level_objects[level])
self._leaf = _fluid_leaves(self.hierarchy, self._boundaries)
self._state = [mask.copy() for mask in self._leaf]
if self._channel:
self._leaf[0] &= channel_masks((cfg.ny, cfg.nx))[0]
# Private constructor inputs for remap_finest.
self._root_cfg = replace(cfg)
self._level_blocks = level_blocks
self._root_cfg = replace(cfg) # private constructor input for remap_finest
self.ctx = context
self._lock = RLock()
self._status = "new"
@@ -391,7 +390,7 @@ class FixedPatchStepper:
for index, layout in enumerate(self._layouts):
total += 2 * 4 * int(np.prod(layout.blocked_shape))
total += layout.halo_dst.nbytes + layout.halo_src.nbytes
b = self._device_boundaries[index]
b = self._boundaries[index]
if b is not None:
total += sum(np.asarray(v).nbytes for k, v in b.tables.items()
if k == "flags" or k in BLOCKED_BOUNDARY_DEVICE_COLUMNS)
@@ -414,7 +413,7 @@ class FixedPatchStepper:
grid = ((config.nx + self._block[0] - 1) // self._block[0], config.ny, layout.nblocks)
level = _Level(shape, module, self._allocate(size), self._allocate(size),
module.get_function(name), grid, layout)
level.boundary = self._device_boundaries[index]
level.boundary = self._boundaries[index]
level.halo_copy = module.get_function("BlockHaloCopy")
if len(layout.halo_dst):
level.halo_dst = self._upload(layout.halo_dst)
@@ -576,12 +575,12 @@ class FixedPatchStepper:
pairs = self.hierarchy.pairs
if (self._wall_aware or self._channel or self.hierarchy.depth != 3
or not isinstance(pairs[1], FixedPatchSet2D)
or any(b is not None for b in self._device_boundaries[:2])):
or any(b is not None for b in self._boundaries[:2])):
raise ValueError("remap_finest requires a periodic three-level hierarchy with a "
"level-2 set and no level-0/1 bodies")
cfg = self._root_cfg
objects = tuple(objects)
actual = self._device_boundaries[2]
actual = self._boundaries[2]
jets = () if actual is None else actual.jets
old = RemapTopology.from_config(cfg, pairs[0], pairs[1].bounds,
child_objects=objects, child_jets=jets,
@@ -603,7 +602,7 @@ class FixedPatchStepper:
or not np.array_equal(old.surface_velocity, new.surface_velocity)):
raise ValueError("remap requires the identical pinned physical boundary and BFL support")
# Check sparse block alignment/launch size before any device read/allocation.
layout = new.child.fine_layout(self._level_blocks[2])
layout = new.child.fine_layout(self._layouts[2].block)
if layout.nblocks > MAX_GRID_DIM or layout.local_shape[0] > MAX_GRID_DIM:
raise ValueError("level 2 blocks exceed launch limits")
with self.ctx.activate():
@@ -621,7 +620,7 @@ class FixedPatchStepper:
phase="postcollision_g", **options)
candidate = type(self)(
cfg, context=self.ctx, refinements=[pairs[0].bounds, new.child.bounds],
level_blocks=list(self._level_blocks),
level_blocks=[layout.block for layout in self._layouts],
level_objects=[(), (), objects], level_jets=[(), (), jets],
observe_raw_exchange=self._observe_raw_exchange,
observe_window_exchange=self._observe_window, explosion=self._explosion)
@@ -651,7 +650,7 @@ class FixedPatchStepper:
if (isinstance(level, bool) or not isinstance(level, Integral)
or not 1 <= level < self.hierarchy.depth):
raise ValueError("Jet level must name an existing refined level (1 .. depth-1)")
boundary = self._device_boundaries[level]
boundary = self._boundaries[level]
# Host preparation: validate and build the new column without touching state.
jets = list(() if boundary is None else boundary.jets)
jets[index] = jets[index].with_amplitude(amplitude) # validates first
@@ -711,7 +710,7 @@ class FixedPatchStepper:
b = level.boundary
if b is not None and b.last_raw is not None:
count = len(b.tables["cell"])
h = self._exchange_spacing[index] # acoustic scaling: dt_l = h_l
h = b.spacing # acoustic scaling: dt_l = h_l
b.accumulate(b.last_raw, b.interval_integral, np.uint32(count),
np.float64(h), np.float64(h),
grid=((count + self._block[0] - 1) // self._block[0], 1), **kw)
@@ -755,7 +754,7 @@ class FixedPatchStepper:
self._exchange_valid = False
for _ in range(n):
if self._observe_raw_exchange: # stream-ordered: holds the last interval only
for b in self._device_boundaries:
for b in self._boundaries:
if b is not None:
cuda.memset_d8_async(b.interval_integral, 0,
len(b.tables["cell"]) * 32, self._stream)
@@ -778,7 +777,7 @@ class FixedPatchStepper:
raise RuntimeError("window exchange observation is disabled")
self._window_valid = False
try:
for b in self._device_boundaries:
for b in self._boundaries:
if b is not None and b.window_integral is not None:
cuda.memset_d8(b.window_integral, 0, len(b.tables["cell"]) * 32)
except BaseException:
@@ -791,7 +790,7 @@ class FixedPatchStepper:
if (isinstance(level, bool) or not isinstance(level, Integral)
or not 0 <= level < len(self._levels)):
raise ValueError("level must name an existing level")
b = self._device_boundaries[level]
b = self._boundaries[level]
if (b is None or isinstance(body_id, bool) or not isinstance(body_id, Integral)
or not 0 <= body_id < len(b.host.objects)):
raise ValueError("body_id must name a body on this level")
@@ -799,9 +798,9 @@ class FixedPatchStepper:
def _summarize(self, b, rows, level, body_id, window, **kw):
return summarize_exchange(rows, b.tables["body_id"], b.tables["aperture"],
body_id=body_id, level=level, h=self._exchange_spacing[level],
body_id=body_id, level=level, h=b.spacing,
time_end=self._step_count, window=window,
reference_center=self._exchange_centers[level][body_id], **kw)
reference_center=b.centers[body_id], **kw)
def read_exchange_window(self, body_id, *, level):
"""Non-destructive read of [window start, current step count] in base-lattice units."""
@@ -888,7 +887,7 @@ class FixedPatchStepper:
self._stream.synchronize()
except Exception as exc:
error = exc
for b in getattr(self, "_device_boundaries", ()):
for b in getattr(self, "_boundaries", ()):
if b is not None:
b.flags = b.device = b.kernel = None
b.last_raw = b.interval_integral = b.accumulate = None
@@ -923,11 +922,11 @@ class FixedPatchStepper:
@property
def boundaries(self):
"""Host boundary geometry per level (None: no bodies/Jets on a refined level)."""
return tuple(None if b is None else b.host for b in self._device_boundaries)
return tuple(None if b is None else b.host for b in self._boundaries)
@property
def jets(self):
return tuple(() if b is None else b.jets for b in self._device_boundaries)
return tuple(() if b is None else b.jets for b in self._boundaries)
@property
def leaf(self):
+2
View File
@@ -130,6 +130,8 @@ class _BlockedBoundary:
interval_integral: object = None
accumulate: object = None
window_integral: object = None # optional observation-window FP64 ledger
spacing: float = None # level spacing h_l in root units (exchange scaling)
centers: tuple = () # body reference centers in the root edge frame
def _prepare_blocked_boundary(boundary, layout, cfg, jets, *, _shifted_blocked=False):
+5 -5
View File
@@ -102,7 +102,7 @@ def test_simulation_before_initialize_commits_host_only(monkeypatch):
def stepper_state(run):
boundary = run._device_boundaries[2]
boundary = run._boundaries[2]
return boundary.tables["surface_velocity"].copy(), boundary.jets, run.jets[2]
@@ -116,8 +116,8 @@ def test_stepper_invalid_input_writes_nothing_and_stays_ready():
with pytest.raises(error):
run.set_jet_amplitude(*args, **kw)
module.cuda.memcpy_htod.assert_not_called()
np.testing.assert_array_equal(run._device_boundaries[2].tables["surface_velocity"], velocity)
assert (run._device_boundaries[2].jets, run.jets[2]) == (jets, level_jets)
np.testing.assert_array_equal(run._boundaries[2].tables["surface_velocity"], velocity)
assert (run._boundaries[2].jets, run.jets[2]) == (jets, level_jets)
assert run._status == "ready"
run.set_jet_amplitude(0, .017, level=2)
assert run.jets[2][0].amplitude == .017
@@ -136,8 +136,8 @@ def test_stepper_device_failure_keeps_host_and_blocks_reuse(where):
run._stream.synchronize.side_effect = RuntimeError("upload failed")
with pytest.raises(RuntimeError, match="upload failed"):
run.set_jet_amplitude(0, .017, level=2)
np.testing.assert_array_equal(run._device_boundaries[2].tables["surface_velocity"], velocity)
assert (run._device_boundaries[2].jets, run.jets[2]) == (jets, level_jets)
np.testing.assert_array_equal(run._boundaries[2].tables["surface_velocity"], velocity)
assert (run._boundaries[2].jets, run.jets[2]) == (jets, level_jets)
assert run._status == "failed"
module.cuda.memcpy_htod.side_effect = None
run._stream.synchronize.side_effect = None
+4 -4
View File
@@ -96,11 +96,11 @@ class MultiplePinnedRemap(unittest.TestCase):
before = run.download()
clone = run.remap_finest((PINS+MOVED+ADDED)[::-1], objects=bodies)
try:
self.assertEqual(tuple(j.amplitude for j in clone._device_boundaries[2].jets), (.017, -.011))
np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
run._device_boundaries[2].tables['surface_velocity'])
self.assertEqual(tuple(j.amplitude for j in clone._boundaries[2].jets), (.017, -.011))
np.testing.assert_array_equal(clone._boundaries[2].tables['surface_velocity'],
run._boundaries[2].tables['surface_velocity'])
clone.set_jet_amplitude(1, .005, level=2)
self.assertEqual(run._device_boundaries[2].jets[1].amplitude, -.011)
self.assertEqual(run._boundaries[2].jets[1].amplitude, -.011)
clone.step()
self.assertEqual((run.step_count, clone.step_count), (7, 8))
for got, want in zip(run.download(), before):
+5 -5
View File
@@ -62,9 +62,9 @@ class PinnedBodyRemap(unittest.TestCase):
run.set_jet_amplitude(0,.017,level=2)
before=run.download()
clone=run.remap_finest(B,objects=objects)
self.assertEqual(clone._device_boundaries[2].jets[0].amplitude,.017)
np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
run._device_boundaries[2].tables['surface_velocity'])
self.assertEqual(clone._boundaries[2].jets[0].amplitude,.017)
np.testing.assert_array_equal(clone._boundaries[2].tables['surface_velocity'],
run._boundaries[2].tables['surface_velocity'])
self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2]))
clone.step()
self.assertEqual((run.step_count,clone.step_count),(7,8))
@@ -74,8 +74,8 @@ class PinnedBodyRemap(unittest.TestCase):
def test_static_body_without_actuation(self):
run,module,allocations,storage,objects=ready(jets=())
clone=run.remap_finest(A,objects=objects)
self.assertFalse(clone._device_boundaries[2].actuated)
self.assertNotIn('surface_velocity',clone._device_boundaries[2].tables)
self.assertFalse(clone._boundaries[2].actuated)
self.assertNotIn('surface_velocity',clone._boundaries[2].tables)
for a,b in zip(run.download(),clone.download()): np.testing.assert_array_equal(a,b)
clone.close();run.close()
+2 -2
View File
@@ -85,7 +85,7 @@ class RawExchange(unittest.TestCase):
def test_four_level_schedule_publication_and_fixed_storage(self):
run, host, calls, arrays = backend()
b = run._device_boundaries[3]
b = run._boundaries[3]
rows = len(b.tables['cell'])
self.assertEqual(b.last_raw.size, rows * 16)
self.assertEqual(b.interval_integral.size, rows * 32)
@@ -128,7 +128,7 @@ class RawExchange(unittest.TestCase):
def test_disabled_has_nulls_and_no_observation_allocation_or_launch(self):
run, host, calls, _ = backend(False)
b = run._device_boundaries[3]
b = run._boundaries[3]
self.assertIs(b.last_raw, b.interval_integral)
self.assertIs(b.interval_integral, None)
run.step()
+14 -10
View File
@@ -1,5 +1,6 @@
"""Host guards and source contract for MRT reuse in refined/blocked stepping (no CUDA)."""
from pathlib import Path
import re
import unittest
from CelerisLab.config import LBMConfig
@@ -32,22 +33,25 @@ class RefinementMRTHost(unittest.TestCase):
"collide_srt_no_force(f, feq, d_params.omega);"):
self.assertIn(text, block_collision)
self.assertIn("collide_srt(f, feq, forcing, omega);", block_collision)
self.assertIn("#elif BLOCK_ROOT_BC == 2", block_collision)
self.assertIn("#if COLLISION_MODEL == 2 || BLOCK_ROOT_BC == 2", block_collision)
self.assertEqual(channel.count("block_collide(f);"), 1)
self.assertNotIn("collide_mrt_unforced(f);", channel)
# Executable text only: the header comment explains the absent force read.
code = re.sub(r"//[^\n]*", "", channel + "\n" + block_collision)
self.assertIn("d_params.omega", code)
for text in ("collide_mrt_unforced", "collide_dispatch(", "d_params.fx", "d_params.fy"):
self.assertNotIn(text, code)
for source in (blocks,):
self.assertIn("(COLLISION_MODEL != 0 && COLLISION_MODEL != 2)", source)
self.assertIn("USE_LES", source)
self.assertIn("STORE_PRECISION != 0", source)
self.assertIn('#if COLLISION_MODEL == 2\n#include "operators/collision_mrt.cuh"\n#endif', source)
# 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 + "\n" + block_collision, channel + "\n" + block_collision):
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.assertNotRegex(source, r"compute_rho_u\(f, rho, ux, uy\);\s*collide_mrt")
# One unforced MRT call, in block_collide only; the steps never collide directly.
self.assertEqual(block_collision.count("collide_mrt_no_force(f, rho, ux, uy, omega);"), 1)
self.assertEqual(block_collision.count("compute_rho_u(f, rho, ux, uy);"), 1)
self.assertEqual(block_collision.count("fmaxf(OMEGA_COLLISION_MIN, d_params.omega)"), 1)
for source in (blocks, channel):
self.assertNotIn("collide_mrt", source)
self.assertNotIn("collide_srt", source)
if __name__ == "__main__": unittest.main()
+3 -3
View File
@@ -70,11 +70,11 @@ class RuntimeChildRemap(unittest.TestCase):
run._status = 'ready'
for bad in ([], [(40,16,48,40)], [(16,16,32,24),(16,27,32,55)]):
with self.assertRaises(ValueError): run.remap_finest(bad)
records, pairs = run._device_boundaries, run.hierarchy.pairs
records, pairs = run._boundaries, run.hierarchy.pairs
for attribute, value in (('_wall_aware', True),
# level-2 record whose host geometry differs from objects=()
('_device_boundaries', [*records[:2], SimpleNamespace(host=object(), jets=())]),
('_device_boundaries', [records[0], object(), records[2]]), # level-1 body
('_boundaries', [*records[:2], SimpleNamespace(host=object(), jets=())]),
('_boundaries', [records[0], object(), records[2]]), # level-1 body
('hierarchy', replace(run.hierarchy, pairs=(pairs[0], None)))):
with patch.object(run, attribute, value):
with self.assertRaises(ValueError): run.remap_finest(B)
+4 -4
View File
@@ -67,7 +67,7 @@ class ShiftedJets(unittest.TestCase):
def test_switch_upload_only_then_raw_remap_keeps_live_velocity(self):
run,module,allocations,storage,objects=ready();clone=None
try:
boundary=run._device_boundaries[2]
boundary=run._boundaries[2]
before=run._download_raw()
original=boundary.tables['surface_velocity'].copy()
module.cuda.memcpy_htod.reset_mock()
@@ -79,7 +79,7 @@ class ShiftedJets(unittest.TestCase):
for a,b in zip(before.values,run._download_raw().values):self.assertEqual(a.tobytes(),b.tobytes())
clone=run.remap_finest(B,objects=objects)
self.assertEqual(clone.jets[2][0].amplitude,-jet().amplitude)
np.testing.assert_array_equal(clone._device_boundaries[2].tables['surface_velocity'],
np.testing.assert_array_equal(clone._boundaries[2].tables['surface_velocity'],
boundary.tables['surface_velocity'])
self.assertEqual(clone.step_count,7)
saved=run._download_raw()
@@ -92,12 +92,12 @@ class ShiftedJets(unittest.TestCase):
def test_invalid_amplitude_does_not_upload_or_commit(self):
run,module,_,_,_=ready()
try:
before=run._device_boundaries[2].tables['surface_velocity'].copy()
before=run._boundaries[2].tables['surface_velocity'].copy()
module.cuda.memcpy_htod.reset_mock()
for amplitude in (float('nan'),float('inf'),1.):
with self.assertRaises(ValueError):run.set_jet_amplitude(0,amplitude,level=2)
module.cuda.memcpy_htod.assert_not_called()
np.testing.assert_array_equal(before,run._device_boundaries[2].tables['surface_velocity'])
np.testing.assert_array_equal(before,run._boundaries[2].tables['surface_velocity'])
finally:run.close()
+1 -1
View File
@@ -73,7 +73,7 @@ class SingleLevel(unittest.TestCase):
def test_circle_solids_are_not_state(self):
run, _, _ = backend((circle(),))
solid = run._device_boundaries[0].tables["solid"]
solid = run._boundaries[0].tables["solid"]
self.assertTrue(solid.any())
np.testing.assert_array_equal(run._state[0], ~solid)
+1 -1
View File
@@ -37,7 +37,7 @@ class UnionBodySupport(unittest.TestCase):
run.set_jet_amplitude(0,.017,level=2)
clone=run.remap_finest(((16,16,24,40),(24,16,32,40),A[1]),objects=objects)
try:
self.assertEqual(clone._device_boundaries[2].jets[0].amplitude,.017)
self.assertEqual(clone._boundaries[2].jets[0].amplitude,.017)
self.assertTrue(same_physical_boundary(run.boundaries[2],clone.boundaries[2]))
for a,b in zip(run.download(),clone.download()):np.testing.assert_array_equal(a,b)
finally:clone.close()
@@ -240,7 +240,7 @@ class SingleLevelGPU(unittest.TestCase):
with CudaContext(0) as context:
for split in SPLITS:
with single(context=context, block=split, objects=[circle()]) as backend:
tables = backend._device_boundaries[0].tables
tables = backend._boundaries[0].tables
for key, values, dtype in (("q", [l.q for l in links], np.float32),
("dir", [l.dir for l in links], np.uint8),
("fallback", [l.fallback for l in links], np.uint8)):
@@ -59,7 +59,7 @@ class ExchangeWindowCPU(unittest.TestCase):
stack.enter_context(patch.object(cuda, "memcpy_dtoh",
side_effect=lambda out, src: out.fill(0)))
backend = build(m)
b = backend._device_boundaries[1]
b = backend._boundaries[1]
self.assertIsNotNone(b.window_integral)
with self.assertRaises(RuntimeError):
backend.start_exchange_window() # not initialized
@@ -96,7 +96,7 @@ class ExchangeWindowCPU(unittest.TestCase):
with mocked_backend() as m, ExitStack() as stack:
stack.enter_context(patch.object(cuda, "memset_d8_async"))
backend = build(m, window=False)
self.assertIsNone(backend._device_boundaries[1].window_integral)
self.assertIsNone(backend._boundaries[1].window_integral)
backend.initialize(*[np.full(l.shape, 1 / 9, np.float32) for l in backend._levels])
backend.step(1)
self.assertEqual(len(accumulations(m)), 2)
@@ -120,7 +120,7 @@ class ExchangeWindowGPU(unittest.TestCase):
windowed.step(2)
plain.step(2)
windowed.start_exchange_window()
b = windowed._device_boundaries[1]
b = windowed._boundaries[1]
count = len(b.tables["cell"])
factors = np.array([.25, .25, .125, .25]) # h*dt, h*h*dt, h*h for level-1 h=.5
scale = np.zeros(4)
@@ -212,7 +212,7 @@ class ExchangeWindowShiftedMRTGPU(unittest.TestCase):
run.step(2)
if run is not plain:
run.start_exchange_window()
b = windowed._device_boundaries[1]
b = windowed._boundaries[1]
count = len(b.tables["cell"])
record = dict(scale=np.zeros(4), substeps=0)
audit_substeps(windowed, b, count, record)
@@ -228,7 +228,7 @@ class ExchangeWindowShiftedMRTGPU(unittest.TestCase):
same_bits(self, a, e)
rows = [np.empty((count, 4), np.float64) for _ in range(2)]
cuda.memcpy_dtoh(rows[0], b.window_integral)
cuda.memcpy_dtoh(rows[1], split._device_boundaries[1].window_integral)
cuda.memcpy_dtoh(rows[1], split._boundaries[1].window_integral)
np.testing.assert_array_equal(rows[0].view(np.uint64), rows[1].view(np.uint64))
out = windowed.read_exchange_window(0, level=1)
@@ -57,7 +57,7 @@ def make(context):
def topology(run, rectangles, objects):
return RemapTopology.from_config(
CFG, run.hierarchy.pairs[0], rectangles, child_objects=objects,
child_jets=run._device_boundaries[2].jets,
child_jets=run._boundaries[2].jets,
)
@@ -74,7 +74,7 @@ class MultiplePinnedRemapGPU(unittest.TestCase):
def change_amplitudes(self, run, amplitudes):
# Each jet selects a different body. Updating one must preserve the
# other body's current per-link velocities, including its driven links.
boundary = run._device_boundaries[2]
boundary = run._boundaries[2]
body_ids = run.boundaries[2].curved_host[3]
for index, amplitude in enumerate(amplitudes):
before = boundary.tables['surface_velocity'].copy()
@@ -91,7 +91,7 @@ class MultiplePinnedRemapGPU(unittest.TestCase):
self.assertTrue(same_physical_boundary(new.boundaries[2], mesh.boundary))
# Compare physical cut-link rows, never layout-dependent blocked addresses.
for run in (old, new):
boundary = run._device_boundaries[2]
boundary = run._boundaries[2]
self.assertEqual(tuple(jet.amplitude for jet in boundary.jets), amplitudes)
exact(self, boundary.tables['surface_velocity'], mesh.surface_velocity,
'latest physical per-link velocities')
@@ -39,9 +39,9 @@ class PinnedBodyRemapGPU(unittest.TestCase):
exact(self,getattr(after.boundaries[2],name),getattr(before.boundaries[2],name),name)
for a,b in zip(after.boundaries[2].curved_host,before.boundaries[2].curved_host):
exact(self,a,b,'physical cut links')
exact(self,after._device_boundaries[2].tables['surface_velocity'],
before._device_boundaries[2].tables['surface_velocity'],'current per-link velocities')
self.assertEqual(after._device_boundaries[2].jets[0].amplitude,.017)
exact(self,after._boundaries[2].tables['surface_velocity'],
before._boundaries[2].tables['surface_velocity'],'current per-link velocities')
self.assertEqual(after._boundaries[2].jets[0].amplitude,.017)
def test_noop_next_step_bitwise(self):
with CudaContext(0) as ctx, ExitStack() as owned:
@@ -63,7 +63,7 @@ class PinnedBodyRemapGPU(unittest.TestCase):
current.set_jet_amplitude(0,.017,level=2)
for rects in (B,A):
previous=current;before=previous.download();count=previous.step_count
jets=previous._device_boundaries[2].jets
jets=previous._boundaries[2].jets
old=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],previous.hierarchy.pairs[1].bounds,
child_objects=objects,child_jets=jets)
new=RemapTopology.from_config(CFG,previous.hierarchy.pairs[0],rects,
@@ -105,7 +105,7 @@ class RawExchangeGPU(unittest.TestCase):
arrays = [state(l.shape) for l in measured._levels]
measured.initialize(*arrays)
plain.initialize(*arrays)
b = measured._device_boundaries[1]
b = measured._boundaries[1]
count = len(b.tables['cell'])
kernel = measured._levels[1].kernel
samples = []
@@ -158,7 +158,7 @@ class RawExchangeGPU(unittest.TestCase):
FixedPatchStepper(cfg, context=ctx, **settings) as plain:
arrays = [state(l.shape, ux=0.) for l in run._levels]
run.initialize(*arrays); plain.initialize(*arrays)
b = run._device_boundaries[3]
b = run._boundaries[3]
count = len(b.tables['cell'])
kernel = run._levels[3].kernel
samples, independent_samples = [], []
@@ -213,7 +213,7 @@ class RawExchangeGPU(unittest.TestCase):
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._device_boundaries[1]
b = amr._boundaries[1]
count = len(b.tables['cell'])
links = obj.geometry.build_cut_links(64, 48)
np.testing.assert_array_equal([v.dir for v in links], b.tables['dir'])
@@ -221,7 +221,7 @@ class RawExchangeGPU(unittest.TestCase):
# Same fine-grid force arms after translation to the same center: the
# single level's static BFL gets the AMR Jet column and force arms; the raw
# output pointer is substituted at launch (no stepper-side accumulation).
root = uniform._device_boundaries[0]
root = uniform._boundaries[0]
velocity = uniform._upload(b.tables['surface_velocity'])
root.device.update(surface_velocity=velocity, rx=uniform._upload(b.tables['rx']),
ry=uniform._upload(b.tables['ry']))
@@ -187,8 +187,8 @@ class BlockedFineCPU(unittest.TestCase):
backend, calls = self.mocked_backend([circle(spec)])
c, f = backend._levels
coarse, cur, spare = c.current, f.current, f.spare
flags, kind = backend._device_boundaries[1].flags, backend._levels[1].kind
velocity = backend._device_boundaries[1].device.get("surface_velocity")
flags, kind = backend._boundaries[1].flags, backend._levels[1].kind
velocity = backend._boundaries[1].device.get("surface_velocity")
self.assertEqual(velocity is None, spec is None)
backend.step(1)
backend.close()
@@ -213,8 +213,8 @@ class BlockedFineCPU(unittest.TestCase):
def test_actuation_without_velocity_buffer_fails_before_bfl(self):
for speed in (0.0, .002): # any spec, even zero amplitude
backend, calls = self.mocked_backend([circle(actuation(speed))])
self.assertIsNotNone(backend._device_boundaries[1].tables["surface_velocity"])
del backend._device_boundaries[1].device["surface_velocity"]
self.assertIsNotNone(backend._boundaries[1].tables["surface_velocity"])
del backend._boundaries[1].device["surface_velocity"]
with self.assertRaisesRegex(RuntimeError, "Surface velocity device buffer is missing"):
backend.step(1)
names = [name for name, _ in calls]
@@ -62,8 +62,8 @@ class JetStepperCPU(unittest.TestCase):
C = (16.0, 12.0)
whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.)
backend, _ = mocked([circle()], [whole], FINE_BLOCKS[1])
_, jet_id = backend._device_boundaries[1].jet_cache
self.assertGreater(len(np.unique(backend._device_boundaries[1].tables["block"][jet_id == 0])), 1)
_, jet_id = backend._boundaries[1].jet_cache
self.assertGreater(len(np.unique(backend._boundaries[1].tables["block"][jet_id == 0])), 1)
backend.close()
def test_fine_ownership_uses_remapped_physical_owner(self):
@@ -76,10 +76,10 @@ class JetStepperCPU(unittest.TestCase):
source = SimObject(99, geometry=CircleGeometry(*C, 2.15), center=C, radius=2.15)
whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.)
backend, _ = mocked([source], [whole], FINE_BLOCKS[1])
_, jet_id = backend._device_boundaries[1].jet_cache
_, jet_id = backend._boundaries[1].jet_cache
owners = backend.boundaries[1].curved_host[3][jet_id == 0]
self.assertEqual(np.unique(owners).tolist(), [0])
self.assertGreater(len(np.unique(backend._device_boundaries[1].tables["block"][jet_id == 0])), 1)
self.assertGreater(len(np.unique(backend._boundaries[1].tables["block"][jet_id == 0])), 1)
self.assertEqual(source.obj_id, 99)
backend.close()
cfg = config("periodic")
@@ -97,8 +97,8 @@ class JetStepperCPU(unittest.TestCase):
def test_blocked_jet_velocity_matches_legacy_rows(self):
legacy, _ = mocked([circle(actuation(SPEED))], (), FINE_BLOCKS[1])
driven, _ = mocked([circle()], [jet()], FINE_BLOCKS[1])
a = legacy._device_boundaries[1].tables["surface_velocity"]
b = driven._device_boundaries[1].tables["surface_velocity"]
a = legacy._boundaries[1].tables["surface_velocity"]
b = driven._boundaries[1].tables["surface_velocity"]
diff = np.abs(a - b).max()
print(f"fine jet vs legacy velocity rows: max |dU| = {diff:.3e} (exact={diff == 0})")
self.assertLessEqual(diff, 4 * np.finfo(np.float32).eps * SPEED)
@@ -106,8 +106,8 @@ class JetStepperCPU(unittest.TestCase):
def test_blocked_amplitude_update_rewrites_only_velocity_column(self):
backend, htod = mocked([circle()], [jet()], FINE_BLOCKS[1])
before = backend._device_boundaries[1].tables["surface_velocity"].copy()
device = backend._device_boundaries[1].device["surface_velocity"]
before = backend._boundaries[1].tables["surface_velocity"].copy()
device = backend._boundaries[1].device["surface_velocity"]
with patch.object(cuda, "memcpy_htod", htod):
htod.reset_mock()
backend.set_jet_amplitude(0, -SPEED, level=1)
@@ -115,7 +115,7 @@ class JetStepperCPU(unittest.TestCase):
self.assertEqual(htod.call_count, 1)
self.assertIs(ptr, device)
np.testing.assert_allclose(sent, -before, rtol=1e-6, atol=0)
np.testing.assert_array_equal(backend._device_boundaries[1].tables["surface_velocity"], sent)
np.testing.assert_array_equal(backend._boundaries[1].tables["surface_velocity"], sent)
with self.assertRaises(ValueError):
backend.set_jet_amplitude(0, 1.5, level=1)
self.assertEqual(htod.call_count, 1)
@@ -133,7 +133,7 @@ class JetStepperCPU(unittest.TestCase):
with self.assertRaises(IndexError):
backend.set_jet_amplitude(1, -SPEED, level=1)
self.assertEqual(backend._status, status)
before = backend._device_boundaries[1].tables["surface_velocity"].copy()
before = backend._boundaries[1].tables["surface_velocity"].copy()
boom = Mock(side_effect=RuntimeError("copy failed"))
with patch.object(cuda, "memcpy_htod", boom), patch.object(cuda, "memcpy_htod_async", boom), \
patch.object(cuda, "pagelocked_empty", side_effect=lambda shape, dtype: np.empty(shape, dtype)):
@@ -141,7 +141,7 @@ class JetStepperCPU(unittest.TestCase):
backend.set_jet_amplitude(0, -SPEED, level=1)
self.assertEqual(backend._status, "failed", fine_block)
self.assertEqual(backend.jets[1][0].amplitude, SPEED)
np.testing.assert_array_equal(backend._device_boundaries[1].tables["surface_velocity"], before)
np.testing.assert_array_equal(backend._boundaries[1].tables["surface_velocity"], before)
with self.assertRaisesRegex(RuntimeError, "failed"):
backend.step(1)
with self.assertRaisesRegex(RuntimeError, "failed"):
@@ -102,7 +102,7 @@ class Ledger:
def step(self):
run=self.run;self.samples.clear();run.step()
b=run._device_boundaries[self.index];h=run._exchange_spacing[self.index]
b=run._boundaries[self.index];h=b.spacing
self.test.assertEqual(len(self.samples),2**self.index)
expected=np.zeros_like(self.samples[0],dtype=np.float64)
weights=np.array([h*h,h*h,h*h*h,h*h])
@@ -140,7 +140,7 @@ class ShiftedJetsGPU(unittest.TestCase):
level_objects=[(),(obj,)],level_jets=[(),(whole,)],observe_raw_exchange=True) as run:
run._initialize_raw(initial(run,rho,u=.013,v=0.))
expected=run._download_raw();ledger=Ledger(self,run,1)
self.assertTrue(run._device_boundaries[1].tables['aperture'].all())
self.assertTrue(run._boundaries[1].tables['aperture'].all())
for _ in range(4):ledger.step()
error=max(float(abs(a[:,leaf]-b[:,leaf]).max()) for a,b,leaf in
zip(run._download_raw().values,expected.values,run.leaf))
@@ -201,7 +201,7 @@ class ShiftedJetsGPU(unittest.TestCase):
self.exact_raw(before,old._download_raw())
self.exact_raw(_StoredPopulations(expected.g,True),new._download_raw())
self.assertEqual(new.jets[2][0].amplitude,-.001)
np.testing.assert_array_equal(new._device_boundaries[2].tables['surface_velocity'],old._device_boundaries[2].tables['surface_velocity'])
np.testing.assert_array_equal(new._boundaries[2].tables['surface_velocity'],old._boundaries[2].tables['surface_velocity'])
newledger=Ledger(self,new,2);ledger.step();newledger.step()
self.assertEqual((old.step_count,new.step_count),(count+1,count+1))
ledger.report('L2 source');newledger.report('L2 remapped')