docs: Dütsch protocol — adopt Grok's cycle-index labels (K = first scored cycle)

Same cycles and pairs as 8505a723, relabelled: proof pairs 3→4..5→6 give K_case = 7, window K..K+2 with remain-settled pairs, last scored cycle 11. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
docs: Dütsch protocol — worked cycle-index example, sparse sampling, driver unit tests
2026-10-30 10:50:22 +00:00 · 2026-10-30 10:48:30 +00:00 · 2026-10-30 10:44:38 +00:00 · 2026-10-30 10:43:48 +00:00 · 2026-10-30 10:43:10 +00:00 · 2026-10-30 10:30:38 +00:00
36 changed files with 2070 additions and 198 deletions
@@ -1,8 +1,14 @@
 # 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). All three received parent/candidate
+generated-config and host/GPU equality checks. Only the collision slice ran a new PTX
+matrix (a4fc46d, whose src is identical to c5a3cdb); 65339ab and 442979b reused prior
+PTX evidence after unchanged CUDA/compiler and generated-input checks. In the collision
+matrix, only the two channel-root SRT modules changed PTX (omega load/clamp moved);
+their parent/candidate outputs matched exactly in the exercised cases.

 ## Conventions (defined once, in `lbm.hierarchy.Hierarchy`)

@@ -28,20 +34,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 +67,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 +105,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).
@@ -0,0 +1,403 @@
+# Moving circles on the uniform dense grid — S0 DRAFT contract
+
+Status: **DRAFT, not frozen** for the force RHS, kernels and physical
+acceptance. Host-only S1 (state, validation, guards; no kernels) is on branch
+`claude/moving-bodies-s1`; positive-step runs with a moving body are rejected
+until S2. When S2 lifts that guard, `run(n > 1)` with a pending target must
+raise `ValueError` before mutation (4a), and `run(0)` must not upload a motion
+record or consume the pending target. Base: main
+`b2034f87e03656c5f67ab88562a70d6c338b9918`. Scope approved by Frank (via TRON,
+thread C0C78UFKGRM/1793335557.460509): uniform-grid prescribed translation and
+single-cylinder free transverse VIV. Out of scope: refinement/blocked path,
+pinball assembly, EsoPull, shifted storage, D3Q19, DynamisLab, manuscript.
+
+Items marked **OPEN** are undecided; items marked **CANDIDATE** are proposals.
+
+## 1. Working standards
+
+- Lean, readable kernels in the existing idiom (`kernels/step/*.cu`, shared
+  `apply_bouzidi_link`); Python manages state and validation. No framework,
+  no compatibility layer, no new dependency.
+- The existing stationary/rotating path (`CurvedBoundaryKernel`, host link
+  packing, `sync_bodies`) stays byte-identical. Moving circles are a separate,
+  explicitly enabled path with its own numerical contract (this file).
+- One authority per state (section 4). No per-step full-domain download,
+  rebuild or recompile; local swept-band device updates are in scope.
+- Small slices (S1 host spec, S2 kernels without motion, S3 motion + refill,
+  S4 oscillator); each reviewed SHA is frozen before targeted runs.
+- Tolerances, benchmark parameters and the force composition used in the
+  structural RHS are declared before the runs that judge them. Fresh vs carried
+  evidence and existing baseline failures are recorded; smoke, convergence and
+  physical validation are reported separately. No undeclared force filtering.
+- External code (waLBerla GPLv3, FluidX3D source-available, VIVSIM) is a
+  concept reference only; nothing is copied.
+
+## 2. Current update order (main b2034f87, dense double buffer)
+
+`Simulation.run(n)`: upload host action (if `upload_act`) → zero obs (if
+`zero_obs`) → `DoubleBufferStepper.step(n)` → download obs → synchronize.
+
+Per step (`lbm/stepper.py` `step`):
+
+1. `CurvedBoundaryKernel` over the packed cut links: reads post-collision
+   `G^{n-1}` at `k_f` (and `k_f - c` for q < 1/2), adds the moving-wall term
+   (`Uw = -omega*ry, Vw = omega*rx` + optional per-link `surface_velocity`),
+   writes the reflected value into slot `opp` of the solid node `k_s`, and
+   accumulates `c(f_toward + f_reflected)` and torque into `obs`.
+2. `OneStep` pull/collide `ddf → temp` (collides FLUID cells; obstacle
+   interiors flagged `SOLID|OBSTACLE|BC_CURVED`), then buffer swap.
+3. `ForceRegionKernel`, `SensorKernel`.
+
+Raw `obs` force is the force **on the body** (link derivation in section 6;
+DFG docs report `Cd > 0` from it). Translation `vx, vy` in the action slots is
+currently stored but not read by any curved kernel.
+
+## 3. Proposed update order (moving circles)
+
+Per step n, before the existing step 1:
+
+- **M1 `MovingCircleCells`** (one thread per swept-band cell): classify with the
+  analytic predicates `in_old = inside(x, c^{n-1})`, `in_new = inside(x, c^n)`;
+  write `flag = in_new ? SOLID|OBSTACLE|BC_CURVED : FLUID`; covered cells
+  (`!in_old && in_new`) accumulate `node_cover`; fresh cells
+  (`in_old && !in_new`) are refilled (section 5) and accumulate `node_uncover`.
+- **M2 `MovingCircleLinks`** (one thread per band cell × direction): compute
+  q and hit point at `c^n` (FP64, same formula as `find_circle_intersection`,
+  same clamp [0.01, 1] and fallback classes as `CircleGeometry.build_cut_links`),
+  skip non-cut slots, call `apply_bouzidi_link` with
+  `u_w = (vx^n, vy^n) + omega × (hit - c^n)`, accumulate the link terms
+  (section 6).
+
+Then the unchanged steps 1–3. M2 must follow M1 because a q < 1/2 link may
+read `k_f - c` at a fresh cell.
+
+Exclusivity: moving bodies are never packed into the `CurvedBoundaryKernel`
+link arrays (`n_curved` covers stationary/rotating bodies only), so M2 is the
+sole writer of their solid slots and the sole source of their link force; no
+double write and no double `obs`. The stationary path's packing, launch and
+arithmetic are unchanged.
+
+Completeness: after M1, every (fluid cell, direction) pair in the swept band
+whose neighbour is `in_new` is visited by M2. This holds because the band is
+the circle AABB at the declared extreme offsets plus 2 cells, and the host
+rejects per-step moves with `|c_new - c_old|_2 > 0.5` and bands touching the domain edge,
+other bodies, sensors or force regions. Each M2 write lands at
+`index_f(k_s, opp)`, the exact slot the adjacent fluid cell pulls, so every
+fluid pull from a solid slot of a moving body is BFL-written.
+
+## 4. State, masks and authority
+
+- Host `ObjectState` is the single authority for one displacement
+  `offset = (dx, dy)` relative to the immutable base centre, plus
+  `velocity = (vx, vy)` and `omega`. Centres are always derived
+  (`c = base + offset`), never stored as independent state. Velocity is
+  supplied, never differenced from offsets.
+- `c^{n-1}` is derived from the last *committed* offset, `c^n` from the
+  offset being applied. Per body a small device **motion record**
+  `{c^{n-1}, c^n, u^n}` is uploaded with the action buffer; action-slot meaning
+  is not overloaded. (CANDIDATE layout; OPEN.)
+- Fail-closed commit: an enqueued upload is not a completed one. The new
+  offset becomes the committed offset only after the stream synchronization
+  that follows the upload and the step group succeeds. Any upload or later
+  synchronization failure leaves host and device motion state unknown: the
+  simulation is marked unusable (existing `_require_usable` lifecycle), and
+  no API may report a host/device match afterwards.
+
+### 4a. Applied geometry vs pending target (CANDIDATE, narrow API)
+
+Two distinct facts per moving body, one authority each:
+`applied` = offset of the geometry last applied on the device, changed only
+by a successful commit; `pending` = an optional requested target
+`(offset, velocity)`, changed only by `set_body`. A setter returning proves
+nothing about the device. The uploaded motion record is derived each run:
+`{c_old = base + applied, c_new = base + (pending.offset if pending is not
+None else applied), u = pending.velocity if pending is not None else 0}`
+(explicit `is not None`, so a zero offset or zero velocity is a valid target). Translational wall velocity exists only on a step
+that advances to a pending target; a hold step (`c_old == c_new`) has
+translational `u_w = 0` (spin unchanged), so a non-advancing circle never
+acts as a sliding wall. One moving step consumes one target.
+
+| situation | contract |
+|---|---|
+| `set_body` target | validated against `applied` (box, `|Δc|_2 ≤ 0.5` cell, speed bound) before mutation; on success replaces `pending` |
+| target with `offset == applied` and nonzero translational velocity | `ValueError` before mutation (would be a hold step with sliding wall velocity) |
+| two setters before a run | last wins; each validated against `applied`, not against the earlier pending target |
+| `run(1)` with pending | upload record, step, sync; on success `applied := pending.offset` (applied stores only the offset), `pending := None` |
+| `run(1)` without pending | hold: `c_old == c_new`, no cover/uncover, links at `applied`, translational `u_w = 0` |
+| `run(n > 1)` with pending | `ValueError` before mutation (single-step-only moving API) |
+| `run(n > 1)`, no pending | allowed: hold for n steps, translational `u_w = 0`, spin only |
+| `run(0)` | no motion upload, no commit; `pending` kept |
+| validation failure | `ValueError`, nothing mutated, retry allowed |
+| upload/launch/sync failure after enqueue | simulation unusable; `applied` and `pending` are not reported as device truth; no re-entry |
+
+Displacement norm: Euclidean per-step `|c_new - c_old|_2 ≤ 0.5` lattice
+cells (declared bound; CANDIDATE value).
+- Membership is never read from a flag that M1 is overwriting. The host
+  verifies at `add_body` that every swept-band cell has base flag FLUID, so the
+  band's base state is known without reading `flag_gpu`.
+- `flag_gpu` is the device authority for band cells while a moving body
+  exists; the host flag mirror is stale for those cells and is not used.
+
+## 5. Refill of fresh cells (CANDIDATE)
+
+- Stored time level: everything in `ddf` at M1 is post-collision `G^{n-1}`; the
+  refill writes all 9 slots at that same level, so every slot of a fresh cell
+  that the following pull or a q < 1/2 link may read is defined. Donor values
+  are read at the same level `G^{n-1}` (old-fluid cells, untouched by M1).
+- `G = f_eq(rho_bar, u_w(x)) + G_neq(x_d)` where `rho_bar` = mean density of
+  valid neighbours and `x_d` = the valid neighbour best aligned with the outward
+  normal `(x - c^n)/|x - c^n|`; `G_neq(x_d) = G(x_d) - f_eq(rho(x_d), u(x_d))`.
+- Valid donor: inside the domain, base FLUID, `!in_old && !in_new`.
+- Sets: `V(x)` = the 8 lattice neighbours that are in the domain, base FLUID
+  and `!in_old && !in_new`. `rho_bar` = mean density over `V(x)`. NEQ donor =
+  the member of `V(x)` with the largest `c·n` among those with `c·n > 0`
+  (`n` = outward normal at `c^n`); none → NEQ part 0 (counted as `eq_only`).
+- `V(x)` empty → no defined density: the cell is refilled with
+  `f_eq(RHO, u_w(x))` so all 9 slots stay defined, the `no_donor` counter is
+  incremented, and the host treats `no_donor > 0` after the run as a hard
+  failure (simulation unusable). It is not a silent fallback.
+- Non-emptiness is NOT assumed from the displacement bound. Geometry
+  assumptions: circle radius `R ≥ R_min` (CANDIDATE 4 cells), `|Δc|_2 ≤ 0.5`.
+  S1 includes a host enumeration over sub-cell centre positions and
+  displacement directions on a fine grid within these bounds that checks every
+  fresh cell has nonempty `V(x)`; it is evidence for those bounds, not a
+  proof, and runtime keeps the counter.
+- Fresh cells only write, old-fluid cells only read: no race.
+
+## 6. Discrete momentum ledger (to review; not a selection)
+
+Notation per step n: `F^n` = fluid set after M1; links L = (k_f, dir) at
+`c^n`; `f_t = G_dir(k_f)`; `f_r` = value M2 writes into `(k_s, opp)`.
+Fluid momentum before streaming:
+
+    P_pre = P(G^{n-1} on F^{n-1}) - sum_cover sum_q c_q G_q + sum_fresh sum_q c_q G_q^refill
+
+During the pull, every population entering a fluid cell comes from a fluid
+cell except the BFL-written slots; `f_t` leaves into solid and is never pulled
+by a fluid cell. Collision conserves momentum for SRT/MRT without forcing.
+Hence, in exact arithmetic,
+
+    P_fluid^n - P_fluid^{n-1} = - sum_L c (f_t + f_r) - I_cover + I_uncover + B^n + S^n
+
+with `B^n` external-boundary momentum (inlet/outlet/walls) and `S^n` forcing.
+So **`F_CME + node_cover - node_uncover`** is, by construction, the body
+momentum that balances the fluid ledger of *this* storage/update order.
+Branch details: halfway fallback links (`f_r = f_t` plus the moving-wall
+addend) and q < 1/2 links enter the link accumulators exactly like Bouzidi
+links, with `c(f_t + f_r)`. In the q < 1/2 branch `G_dir(k_f - c)` is an
+interpolation input only, not a fluid-momentum sink. Pull-slot completeness
+is the section 3 contract. The scope excludes `ddf_shifting` (rejected), so
+the unshifted force form applies. The balance is a fluid-side identity; it
+does not close body+fluid globally in a domain with inlet/outlet/forcing
+(`B^n`, `S^n`), and the solid-slot writes are boundary plumbing, not physical
+streaming packets.
+
+The Wen 2014 relative-velocity link sum adds `-u_w (f_t - f_r) = u_w dm` per
+link (`dm = f_r - f_t`, the BFL mass defect). That term has no counterpart in
+the fluid ledger above. Peng et al. (PRE 2016 submitted, §II.C, Eqs. 20–21)
+show for a straight wall that the bounce-back error accumulated before a node
+switch cancels with the switching momentum. Whether `sum u_w dm` approximates
+`I_cover - I_uncover` here (implying double counting if both are used) is
+**NOT demonstrated** for this implementation: it needs a bounded experiment
+(prescribed slow translation; cumulative `sum u_w dm` vs cumulative node
+impulse, per crossing cycle). Momentum conservation is also distinct from
+Galilean invariance, which G4 tests.
+
+External reference note: waLBerla's covered/refill momentum switches
+(`conserveMomentumWhenMapping`, `conserveMomentumWhenRestoring`) default to
+false (Astra), so their existence does not establish which combination its
+validation cases use. A bookkeeping identity alone does not select the
+physical force either.
+
+Mass ledger (same proposition status). With `M^n` the active-fluid mass
+summed over `F^n` at the sample point,
+
+    M^n - M^{n-1} = sum_L dm - m_cover + m_uncover + B_m^n
+
+(`dm = f_r - f_t`, `m_cover`/`m_uncover` = summed `sum_q G` of covered and
+refilled cells, `B_m` external-boundary mass flux; collision conserves mass).
+
+Sample point (fixed): post-collision storage after the swap of step n, i.e.
+the same level `G` that M1 reads at step n+1; both ledgers are evaluated
+there. Ledger residual = left side minus every listed term. Zero mass flux
+does not imply zero momentum exchange: closed walls (bounce-back, free-slip)
+still exchange momentum, so "no inlet/outlet, no forcing" does NOT give
+`B = S = 0`. The rounding-level momentum check therefore uses either a
+momentum-isolated case (periodic in both directions, no forcing) or counts
+the wall impulse explicitly as part of `B`; the mass check may use a closed
+domain. Otherwise `B`, `S` are not separately measured, the residual is
+reported as their estimate and is not a pass criterion.
+Both ledger identities are **proposed, unreviewed**: they are not called
+proven until review and a host replay agree.
+
+Accumulators (separate, per body, per step window; all reported raw):
+`link_cme`, `link_wen_correction` (= `sum u_w dm`), `node_cover`,
+`node_uncover` (momentum), `link_dm`, `mass_cover`, `mass_uncover` (mass),
+plus counters of crossings, `eq_only` and `no_donor`, so force spikes can be
+correlated with crossing events.
+
+**OPEN — structural RHS candidates (not selected):**
+(a) `link_cme + node_cover - node_uncover`; (b) `link_cme + link_wen_correction`
+(no node impulse). Selection follows the ledger review and the bounded
+experiment above.
+
+## 7. Rejection semantics
+
+- Before any mutation (`ValueError`, no device write, host state unchanged):
+  moving body with EsoPull, `ddf_shifting`, refinement/blocked stepper, D3Q19,
+  Jets or surface actuation on that body, `sync_bodies`/checkpoint while a
+  moving body exists; swept band touching domain boundaries, other bodies,
+  sensors or force regions; `set_body` offset outside the declared box,
+  per-step move `|c_new - c_old|_2 > 0.5` cell (the single bound of 4a), speed
+above the declared bound.
+- Host motion state is committed only after the upload *and* the following
+  synchronization of the step group succeed (section 4); any failure in
+  between makes the simulation unusable and no host/device match is reported.
+- The only soft device diagnostic is `eq_only` (NEQ part omitted, counted).
+  `no_donor > 0` is a hard failure: it is checked after the step group and
+  before the step is committed or any results are returned; on failure the
+  simulation is unusable (section 5).
+
+## 8. Coupling and units (CANDIDATE)
+
+Lattice units, rho0 = 1, per unit span. Step n uses `c^n`, `u^n` and yields
+`F^n`. Weak staggered update
+`v^{n+1} = v^n + (F^n - c v^n - k y^n)/m`, `y^{n+1} = y^n + v^{n+1}`
+(Newmark is an alternative; OPEN). Mass, natural frequency, damping and `U_r`
+follow the chosen benchmark's own definitions. Leontini, Thompson & Hourigan,
+JFS 22:857–864 (2006), DOI 10.1016/j.jfluidstructs.2006.04.003, as read by
+Grok from the full text: Re = 200, `m* = m/m_f = 10` (m_f displaced fluid
+mass), `zeta = c/c_crit = 0.01`, `U* = U/(f_N D)` with `f_N` the natural
+frequency in fluid assuming an inviscid added-mass term; the paper also
+reports frequencies normalised by the in-vacuo value. Our runner derives `k`
+from that in-fluid `f_N` (added-mass coefficient as stated in the paper) and
+never from the vacuum formula; `c_crit` convention and the original
+comparison data are still OPEN. Stability/convergence is gated by fluid
+time/space refinement of the coupled run, not assumed from the mass ratio;
+structural sub-stepping does not substitute for fluid-coupling time
+refinement.
+
+## 9. Gates (map; tolerances OPEN, declared before runs)
+
+G0 zero-motion regression, split: G0a populations of the existing stationary
+and spin-only path bitwise unchanged vs main; G0b moving path at rest (with
+and without spin) reproduces the stationary link set and populations; force
+estimator outputs are compared separately from populations, since the moving
+path adds estimator terms · G1 raw force sign (`Cd > 0`) ·
+G2/G3 prescribed translation incl. ledger accumulators · G4 matched
+relative-flow (moving body in quiescent fluid vs fixed body in uniform
+stream), optional boost, co-moving sanity · G5 Dütsch et al. 1998 Re 100,
+KC 5 in-line (needs 2-component offset) · G6 Leontini free transverse VIV
+(m* = 10, ζ = 0.01, Re = 200, definitions pending) · G7 time-level contract
+frozen before the oscillator.
+
+Candidate pre-run engineering tolerances (Astra; for discussion, not
+literature limits or measured outcomes; norms, averaging window, convergence
+protocol and physical reference uncertainty must be fixed before adoption, and
+none may be weakened after seeing runs):
+normalised momentum-ledger residual FP32 ≤ 5e-6 (FP64 oracle ≤ 2e-13);
+10-cycle relative mass drift ≤ 1e-4; prescribed-motion grid study D = 40/80/160
+with change in the force-response amplitude (observed quantity and norm to be
+named before adoption) ≤ 2 % and phase ≤ 2°; enlarged domain
+amplitude change ≤ 1 %; matched-frame force RMS ≤ 2 % and improving with
+refinement; VIV refinement amplitude ≤ 2 %, frequency ≤ 1 %, five consecutive
+cycle amplitudes within 1 %.
+
+Conditions stated before any moving kernel (S2/S3), not only at the
+oscillator stage: (i) fluid refinement protocol for the prescribed-motion
+gates — diffusive scaling (dx halved, dt quartered, physical Re, KC and
+amplitude fixed), the observed quantities and norms named in advance;
+(ii) Galilean constant-boost G4-boost requires boundaries compatible with a
+uniform boost (periodic in the boost direction or boundaries that move with
+it); if the current dense configuration cannot represent that, G4-boost is
+recorded as not runnable rather than replaced by another test.
+
+## 10. Proposed minimal S1 scope (host only; held until cross-review)
+
+- `add_body` moving-circle declaration: base centre, declared offset box
+  (x/y), speed bound; swept-band computation and band-placement validation.
+- `set_body(offset=..., velocity=..., omega=...)` validation before mutation
+  (section 7), applied/pending state and the run-transition table of 4a,
+  and the fail-closed commit of section 4.
+- Host enumeration check of refill-donor non-emptiness (section 5) within
+  the declared `R_min` and displacement bound.
+- Motion-record host layout and upload path; exclusion of moving bodies from
+  the packed `CurvedBoundaryKernel` links (stationary packing unchanged).
+- Explicit rejections for unsupported modes.
+- S1 execution guard: until the S2 moving kernels exist, any `run(n)` with
+  `n > 0` (and `step(n)`) while a moving body is declared raises before any
+  device mutation (no action/motion upload, no launch), including the
+  nominal at-rest case, because the moving body is excluded from the static
+  packed links and nothing else would apply its boundary. `run(0)` keeps its
+  existing behaviour. Simulations without moving bodies are unaffected.
+  Preparation/commit helpers are exercised by host tests with mocks; no
+  operational unsupported stepping path is exposed.
+- Host tests (Astra runs them): rejection-before-mutation, band geometry,
+  derived centres, fail-closed behaviour under an injected upload/sync
+  failure, stationary packing byte-identical, execution guard (positive-step
+  run with a moving body raises before upload/launch; stationary unaffected). No kernels, no GPU runs, no
+  y-only VIV oscillator (S4), no force RHS.
+
+## 10a. S2a — fixed-centre moving links (implemented on `claude/moving-bodies-s2a`)
+
+- `MovingCircleLinks` (`kernels/step/moving_circle.cu`, dense double buffer)
+  runs before the static curved kernel when moving tables are attached: one
+  thread per (band cell, direction) slot, FP64 `{cx, cy, r}` at the applied
+  offset, q/hit/clamp/fallback in the order of `build_cut_links` with explicit
+  `_rn` operations (host/device bit parity is a test target, not a claim),
+  then the unchanged `apply_bouzidi_link` with spin wall velocity only.
+- `obs` keeps its contract (CME and torque from `apply_bouzidi_link`). A
+  separate per-body diagnostics buffer receives `cme`, the Wen term
+  `u_w * dm`, `link_dm` and the link count from the same link outputs; it is
+  zeroed with obs and only exists when moving circles do. Optional
+  `slot_debug` (q, fallback, rx, ry per slot) is null in production.
+- Displacement is disabled: `set_body(offset=..., vx=..., vy=...)` on a
+  translating circle raises before mutation, before and after initialize, so
+  no pending target can exist and `run`/`step`/direct `LBMStepper.step` all
+  hold at the applied geometry. `initialize()` never reads or consumes
+  `pending`. `stepper.blocked` is cleared only when the moving tables are
+  attached (initialize, stepper recreation).
+- S3 rule (not implemented): once targets are re-enabled, any pending target
+  must raise in `run`/`step` and in `LBMStepper.step` alike before any launch,
+  swap or step-count change; zero-step calls never consume it.
+- Not validation: a pass of the zero-motion gates says nothing about moving
+  geometry. GPU parity (G0a/G0b) is unverified until run on a free device.
+
+## 10b. S3a — translation with refill (implemented on `claude/moving-bodies-s3a`)
+
+- Targets are enabled again (4a). `set_body(offset, vx, vy)` validates against
+  `applied`; a pending target blocks direct `LBMStepper.step(n > 0)` before any
+  launch; `run(n > 1)` with a target raises before mutation; `run(0)` never
+  consumes it; `run(1)` advances it.
+- Per advancing step, on one stream: upload the immutable FP64 records
+  `{c_old, c_new, r, u}` (only when changed) → M1 `MovingCircleCells` → M2
+  `MovingCircleLinks` (`u_w = u + omega x r` at the hit point) → static curved →
+  `OneStep` → swap → force region → sensors → synchronize → read `no_donor` →
+  commit (`applied := pending.offset`) → upload the hold record
+  (`c_old = c_new = applied`, `u = 0`) → synchronize. Hold steps launch no M1
+  and are the S2a path. Commit point: the advancing step-group sync plus a
+  zero `no_donor` count; a failure before it commits nothing. The hold-record
+  upload stages the next step after the commit: if it (or its sync) fails,
+  `applied` already matches the moved device geometry, but the Simulation is
+  still marked failed and refuses reuse. Any failure from the first upload on
+  marks the circles and the Simulation failed.
+- M1 classifies with the old/new predicates only (no flag reads). Covered
+  cells add raw `sum c G` / `sum G` of `G^{n-1}` to `cover_*` and become
+  `SOLID|OBSTACLE|BC_CURVED`. Uncovered cells: donors are fluid in both
+  geometries; `rho_bar` over all of them; NEQ from the donor with the largest
+  unnormalized `c·n > 0` (first direction wins ties, `n` from `c_new`);
+  `u_w` (translation + spin) at the cell position; all nine slots written;
+  raw `sum c G`, `sum G` of the refill go to `uncover_*`; no donor →
+  `f_eq(RHO, u_w)` plus `no_donor` (hard failure before commit); no NEQ donor
+  → `eq_only`.
+- `get_flags()` reads the device flags back (synchronized) when moving tables
+  exist; otherwise it returns the host array as before.
+- Not selected: the structural force. CME+node vs CME+Wen stay OPEN until X1
+  and G4 have run.
+
+## 11. Remaining decisions
+
+Structural RHS (section 6); motion-record layout; refill details; coupling
+scheme; benchmark definitions and digitized reference data; tolerances;
+execution/resource feasibility (Astra, pending); whether a quiescent domain
+(zero inlet velocity) runs cleanly in the current config.
@@ -0,0 +1,209 @@
+# Prescribed in-line oscillation: Dütsch Re 100 / KC 5 — PROTOCOL (review draft)
+
+Production under test: `claude/moving-bodies-s3a` at `653fd997` (unchanged).
+Primary source: Dütsch, Durst, Becker & Lienhart, J. Fluid Mech. 360, 249–271
+(1998), doi:10.1017/S002211209800860X. Source facts below were verified from
+the paper (Slack C0C78UFKGRM/1793335557.460509, TRON 1793356317–1793356724,
+page numbers given). Everything marked ENGINEERING is our a-priori choice,
+fixed before any benchmark result; it is not a paper uncertainty.
+
+## 0. Source facts (with provenance)
+
+- Motion (p.251): `x_c = -A sin(2 pi f t)`; paper phase phi = 0 is the mean
+  position moving left at `-U_max`. KC = 5 gives `A/D = 5/(2 pi) = .7957747`.
+- Start (p.254): the computation starts a quarter cycle before phi = 0, at
+  `x = +A` with zero velocity (phi = -90°).
+- Experiment (pp.252, 254–256): D = 10 mm, L = .65 m, tank 1.2 × .7 × .7 m,
+  rho = 998.2 kg/m³, mu = 1.00818 mPa·s, Re = 100, KC = 5; LDA with 1° phase
+  bins and 200 samples per node and degree; no quantified experimental error
+  recovered.
+- Profiles (Fig. 8, p.260): LAB-FIXED stations x1 = -6, 0, 6, 12 mm, i.e.
+  x/D = -.6, 0, .6, 1.2, at phi = 180°, 210°, 330°; u1 and u2 versus x2.
+- Original numerics (§2): cylinder-fixed accelerated frame (`x_i = x~_i +
+  x_i,s`, outer velocity `-u_c`, inertial source `-rho du_c/dt`), O-grid outer
+  diameter 120 D, illustrated 192 circumferential nodes, 720 steps/cycle; tank
+  influence neglected.
+- Force fit (verified): force ON THE BODY per unit span
+  `F1 = -0.5 rho D c_d xdot|xdot| - (pi/4) rho D² c_i xddot` with `x = -A sin phi`;
+  normalized by `0.5 rho D U_max²`: `c_d cos phi |cos phi| - (pi² c_i / KC) sin phi`.
+  `c_i` is the ADDED-MASS coefficient; total inertia `c_m = c_i + 1` is a
+  different convention and is not used.
+- Table 1 computed Morison coefficients (c_d, c_i) for 96 / 192 / 384 perimeter
+  nodes at dt/T = 1/720: (2.15, 1.42), (2.10, 1.45), (2.09, 1.45). These are
+  fitted computed coefficients, not a measured force and not an error bar. The
+  paper states the coarse direct force is closer to the finest direct force
+  than the Morison fit, but gives no waveform-error norm.
+
+## 1. Kinematics, time levels and phase bookkeeping
+
+- Our run is in the LAB frame: the circle moves, the fluid is initially at rest
+  (rho = 1, u = 0); no fictitious inertial source is added. Lab-fixed profile
+  stations are compared directly (no frame transform needed).
+- Geometry reference: the circle's base centre is the INITIAL position, the
+  right turning point `x_mean + A`; the declared offset box is x in [-2A, 0],
+  y = 0. The body is never teleported: offset 0 at s = 0.
+- Elapsed time s = n dt from the start. Analytic position relative to the mean
+  `x(s) = A cos(2 pi s / T)`, velocity `u(s) = -U_max sin(2 pi s / T)`; paper
+  phase `phi = (360 s / T - 90) mod 360`.
+- Step n moves from s_n to s_{n+1}: target offset `x(s_{n+1}) - A`; wall
+  velocity = the single-step interval average `(x(s_{n+1}) - x(s_n)) / dt`
+  (analytic, supplied, not differenced from stored offsets). The per-step
+  force of `run(1)` / `read_force(normalize=True)` is attributed to that same
+  interval.
+- Paper cycle k = 1, 2, ... spans `s in [T/4 + (k-1) T, T/4 + k T)`; the
+  initial quarter cycle (s < T/4) is cycle 0 and never scored. Every logged
+  sample carries (k, phi). Example: phi = 330° of cycle 1 is at s = 7T/6.
+- T is chosen as a multiple of 360 native steps, so phi = 180°, 210°, 330° fall
+  exactly on step boundaries and 1° bins hold whole steps.
+
+## 2. Resolution ladder (uniform grid, fixed Mach), ENGINEERING
+
+| D (cells) | T (steps) | U_max = 5 D / T | nu_lat = U_max D / 100 | tau (3 nu + .5) |
+|---|---|---|---|---|
+| 20 | 2,520 | .039683 | .0079365 | .5238 |
+| 40 | 5,040 | .039683 | .015873 | .5476 |
+| 80 | 10,080 | .039683 | .031746 | .5952 |
+
+- MRT collision. Mach fixed (Ma ≈ .0687); a separate Mach check at D = 40
+  with T = 10,080 (U_max halved, tau .5238) is labelled MACH and not part of the
+  ladder.
+- Per-step displacement ≤ U_max ≈ .04 cells (bound .5).
+
+## 3. Domain and boundaries (finite-domain reproduction), ENGINEERING
+
+The paper's 120 D O-grid is not reproduced; our rectangular box is a
+reproduction choice. The native dense channel (zero-velocity inlet,
+pressure-NEQ outlet, bounce-back y walls) is a finite-domain approximation,
+not infinite still fluid. The G4 x-periodic / y-slip override is a separately
+labelled sensitivity family and never mixed into the primary ladder. The
+runner owner confirms the zero-velocity inlet is supported as configured
+before any run.
+
+Order (closed):
+1. Smoke + ledger at D = 40, native BC, 40 D × 40 D (finite box).
+2. Domain sensitivity before any literature statement: same BC, D = 40,
+   L ∈ {20 D, 40 D, 60 D}; record the change in Fx waveform, c_d, c_i and the
+   scored profiles.
+3. Optional labelled BC sensitivity: one x-periodic / y-slip case at D = 40,
+   40 D; reported as a delta against native only.
+4. Only after step 2 passes: the D = 20 / 40 / 80 ladder on the chosen BC and
+   domain (default 40 D if, against 60 D, `RMS(dFx)/F_ref` and |dc_d|, |dc_i|
+   stay inside the section 6 bands).
+- Memory/runtime: measured by a smoke run before each new size; measurements
+  inform scheduling only, never tolerances.
+
+## 4. Observables (raw, per step, preserved)
+
+- `obs` Fx, Fy (CME) and all `MovingTables.DIAG` columns per step: CME,
+  CME + Wen and CME + cover − uncover reported separately; crossings, eq_only,
+  no_donor (no_donor must stay 0).
+- Velocity u1, u2 on lab-fixed lines x/D = -.6, 0, .6, 1.2 (relative to the
+  mean cylinder position) versus x2, at the scored phases. Lab frame
+  throughout; if any body-frame intermediate is used, position and velocity
+  are both transformed.
+- Sampling cost: raw force and diagnostics counters are read every step
+  (small). Flags and velocity fields are read only at the steps inside the
+  three prescribed 1° phase bins (φ ∈ [180, 181), [210, 211), [330, 331),
+  T/360 steps each) of the cycles needed for settling and scoring, plus
+  sparse declared audits; never every step.
+- Point set: fixed x2/D sample points per station, identical on every grid
+  and declared before running; per point the interpolation weights are
+  recorded; the scored set is the common subset valid on all compared grids
+  at that phase, with its coverage recorded; point sets never change silently.
+- Interpolation: bilinear on cell centres using only cells that are fluid at
+  that step (flags read back via `get_flags()`); weights renormalized over the
+  fluid corners; a sample with no fluid corner, or inside the circle at that
+  step, is MISSING and logged, never set to zero. One-degree bins keep valid
+  fluid weights. The numerical average need not mimic the 200 LDA samples.
+
+## 5. Periodic-regime gate (every case, before any statistic), ENGINEERING
+
+Ledger rows, same schema for every case: dx, dt (in D and T units); U_max, nu,
+Re, KC (lattice and nondimensional); T (native steps, nondimensional);
+duration (native steps, s/T, paper cycles); discard; accepted window; cycle
+ledger; verdict PASS or FAIL-CLOSED.
+
+Scales: `F_ref = 0.5 rho U_max² D`, `U_ref = U_max`; never a near-zero mean
+force as denominator.
+
+- Settling metrics, cycle k−1 → k, phase-aligned by step within the cycle:
+  - waveform `RMS(dFx) / F_ref ≤ .02`;
+  - `|d mean|`, `|d RMS|`, `|d peak|` of Fx, each `/ F_ref ≤ .01`;
+  - scored profiles (fluid-only) `|du| / U_max`, `|dv| / U_max ≤ .03` at the
+    three phases × four stations;
+  - Morison `|dc_d| ≤ .02`, `|dc_i| ≤ .02`.
+- Rule: all metrics for N = 3 consecutive whole paper cycles → settled.
+
+Cycle indexing (worked example, binding). Paper cycle k spans
+`s/T in [0.25 + (k-1), 0.25 + k)`, i.e. φ0 → φ0; cycle 0 is the startup
+quarter `s/T < 0.25`. Pair p_k compares cycle k−1 with cycle k.
+
+- Discarded: cycles 0, 1, 2 (s/T < 2.25). Pairs are formed only between
+  cycles ≥ 3, so the first pair is p_4 (3 → 4).
+- Settling proof for a case: three consecutive pairs all pass; K_case = the
+  first cycle AFTER the last of those pairs. Earliest possible: pairs 3→4,
+  4→5, 5→6 pass ⇒ K_case = 7.
+- Common window: K = max over compared cases of K_case; every case scores the
+  same three complete cycles K, K+1, K+2 and must stay settled throughout:
+  pairs (K−1→K), (K→K+1), (K+1→K+2) must pass on every case, otherwise that
+  case is FAIL-CLOSED (an earlier local settle does not waive this).
+- Cap: cycle k ends at s/T = k + 0.25, so the last admissible scored cycle is
+  11 (ends 11.25 T) and K ≤ 9. Examples: K = 7 → cycles 7–9 (ends 9.25 T);
+  K = 9 → cycles 9–11 (ends 11.25 T); any K_case > 9 → FAIL-CLOSED.
+- Every case is run to the end of cycle K + 2 (proportionally more native
+  steps on finer grids).
+- Discard: the initial quarter cycle plus the first 2 full paper cycles
+  (s/T < 2.25). Settling is searched from cycle 3; cap 12 periods from s = 0;
+  unsettled → FAIL-CLOSED (no literature, cross-resolution or VIV statistic).
+- Common window: K = the largest first-settled cycle index over all compared
+  cases; every case is scored on the same 3 cycles K .. K + 2 and must reach
+  them within the cap (proportionally more native steps on finer grids), or
+  the comparison is FAIL-CLOSED.
+- The earlier G4 timing check (100 coarse vs 400 diffusive-fine steps) was a
+  matched short transient, not periodic evidence.
+
+## 6. Metrics and per-observable criteria (ENGINEERING, not paper uncertainty)
+
+Three criteria stay separate; the G4 same-grid 2% target is not reused.
+
+- Morison fit: least squares of the normalized per-step force over the
+  accepted window against the two basis functions, each averaged over its
+  single step interval, fixed paper phase, no free phase shift; c_d and
+  added-mass c_i per force candidate.
+- Cross-resolution: every scored observable's error against its reference must
+  decrease monotonically D20 → D40 → D80 (necessary, not sufficient).
+- Against Table 1 (computed fit, not measured force): at D ≥ 40,
+  `|c_d - 2.09| ≤ .15` and `|c_i - 1.45| ≤ .15`, and monotone.
+- Force waveform against its own Morison reconstruction: residual reported
+  only, no pass/fail (the paper reports fit error above its coarse numerical
+  error and gives no percentage).
+- Fig. 8 LDA profiles: QUALITATIVE ONLY. The original figure could not be
+  accessed or rendered; no digitized ordinates and no measurement uncertainty
+  exist. After a digitization with figure number, tool and uncertainty note:
+  `RMS(|du|, |dv|) / U_max ≤ .15` at D ≥ 40 and monotone; LDA error stays
+  UNQUANTIFIED and is stated with every comparison.
+- Fy sanity (not literature): `|mean Fy| / F_ref ≤ .05` in the accepted window.
+- No force candidate is selected for an RHS by these results; each is
+  reported.
+
+Superseded (earlier drafts): the velocity-jump start, guessed stations and
+phases, "C_M" for c_i, and the uncertainty wording; cycles 4–6, the 30-cycle
+cap and the earlier tolerance numbers.
+
+## 6a. Driver unit tests (CPU, before any run)
+
+- Phase/cycle arithmetic: s ↔ (k, φ), cycle boundaries, bin membership,
+  window/cap selection for K = 7, 8, 9, 10.
+- Solid-aware interpolation: all-fluid, partial-solid (renormalized weights),
+  no-fluid (MISSING) and inside-circle cases.
+- Interval-averaged Morison bases and the two-parameter fit on a synthetic
+  waveform with known c_d, c_i.
+
+## 7. Order and scope
+
+Protocol review → one-shot Dütsch driver (analytic `set_body` + `run(1)`,
+cycle ledger, lab-fixed fluid-only station sampler, two-basis Morison fit; no
+reusable runner until the first smoke shows shared helpers pay off; no
+production change) → smoke (memory/runtime) → domain sensitivity →
+ladder → physical-validation verdict → only then single-cylinder VIV. No AMR,
+no pinball, no manuscript work.
@@ -405,7 +405,9 @@ class ObjectManager:

        for obj in objects:
            obj.validate_surface_actuation(self.cfg)
-            if (not obj.is_sensor) and (not obj.is_force_region) and hasattr(obj, 'get_curved_list'):
+            # Moving circles never enter the packed static links (docs/moving_bodies.md, section 3).
+            if (not obj.is_sensor) and (not obj.is_force_region) and obj.motion is None \
+                    and hasattr(obj, 'get_curved_list'):
                (
                    fluid_idx, dirs, q_vals, body_ids,
                    rx_vals, ry_vals, rz_vals, fallback_vals
@@ -0,0 +1,259 @@
+# CelerisLab/body/motion.py
+"""Host state of a translating circle (dense double buffer, D2Q9 only).
+
+Contract: docs/moving_bodies.md, sections 4a, 5 and 7. The applied geometry
+(last successfully committed offset) and an optional pending target are two
+separate facts; centres are always ``base + offset``. Nothing here touches
+the device: the caller supplies upload/step/synchronize callables.
+"""
+from __future__ import annotations
+
+import math
+from dataclasses import dataclass
+from typing import Callable, Optional, Tuple
+
+MAX_STEP = 0.5     # Euclidean centre displacement per step (lattice cells)
+BAND_MARGIN = 2    # cells added around the swept circle AABB
+R_MIN = 4.0        # smallest radius covered by the donor enumeration
+
+Vec = Tuple[float, float]
+
+# D2Q9 non-rest directions, paired order of lbm/descriptors.
+_CX = (1, -1, 0, 0, 1, -1, 1, -1)
+_CY = (0, 0, 1, -1, 1, -1, -1, 1)
+
+
+def _vec(value, name: str) -> Vec:
+    x, y = (float(v) for v in value)
+    if not (math.isfinite(x) and math.isfinite(y)):
+        raise ValueError(f"{name} must be finite, got {value!r}")
+    return x, y
+
+
+@dataclass(frozen=True)
+class Translation:
+    """Declared motion envelope: offset box (x_min, x_max, y_min, y_max) and speed bound."""
+    offset_box: Tuple[float, float, float, float]
+    max_speed: float
+
+    def __post_init__(self):
+        box = tuple(float(v) for v in self.offset_box)
+        speed = float(self.max_speed)
+        if len(box) != 4 or not all(math.isfinite(v) for v in box + (speed,)):
+            raise ValueError("offset_box needs 4 finite values and max_speed must be finite")
+        if not (box[0] <= 0.0 <= box[1] and box[2] <= 0.0 <= box[3]):
+            raise ValueError("offset_box must contain the zero offset")
+        if not 0.0 < speed <= MAX_STEP:
+            raise ValueError(f"max_speed must satisfy 0 < max_speed <= {MAX_STEP}")
+        object.__setattr__(self, "offset_box", box)
+        object.__setattr__(self, "max_speed", speed)
+
+
+class MovingCircle:
+    """Applied offset, pending target and the fail-closed commit of one circle."""
+
+    def __init__(self, base, radius: float, translation: Translation,
+                 on_failure: Optional[Callable[[BaseException], None]] = None):
+        self.base = _vec(base, "base")
+        self.radius = float(radius)
+        if not (math.isfinite(self.radius) and self.radius >= R_MIN):
+            raise ValueError(f"moving circle radius must be >= {R_MIN}")
+        if not isinstance(translation, Translation):
+            raise TypeError("translation must be a Translation")
+        self.translation = translation
+        self.applied: Vec = (0.0, 0.0)
+        self.pending: Optional[Tuple[Vec, Vec]] = None
+        self.failed = False
+        self._on_failure = on_failure  # owner lifecycle hook (Simulation marks itself failed)
+
+    def centre(self, offset: Vec) -> Vec:
+        return self.base[0] + offset[0], self.base[1] + offset[1]
+
+    def check_target(self, offset, velocity) -> Tuple[Vec, Vec]:
+        """Validate a target against ``applied``; returns it, mutates nothing."""
+        self._require_usable()
+        offset, velocity = _vec(offset, "offset"), _vec(velocity, "velocity")
+        x0, x1, y0, y1 = self.translation.offset_box
+        if not (x0 <= offset[0] <= x1 and y0 <= offset[1] <= y1):
+            raise ValueError(f"offset {offset} outside declared box {self.translation.offset_box}")
+        step = math.hypot(offset[0] - self.applied[0], offset[1] - self.applied[1])
+        if step > MAX_STEP:
+            raise ValueError(f"per-step displacement {step} exceeds {MAX_STEP}")
+        if math.hypot(*velocity) > self.translation.max_speed:
+            raise ValueError(f"velocity {velocity} exceeds max_speed {self.translation.max_speed}")
+        if offset == self.applied and velocity != (0.0, 0.0):
+            raise ValueError("a hold target (offset == applied) needs zero translational velocity")
+        return offset, velocity
+
+    def set_target(self, offset, velocity) -> None:
+        """Replace the pending target (last wins); invalid input keeps the old one."""
+        self.pending = self.check_target(offset, velocity)
+
+    def record(self) -> Tuple[Vec, Vec, Vec]:
+        """Motion record (c_old, c_new, u) for the next single moving step."""
+        self._require_usable()
+        if self.pending is not None:
+            offset, velocity = self.pending
+        else:
+            offset, velocity = self.applied, (0.0, 0.0)
+        return self.centre(self.applied), self.centre(offset), velocity
+
+    def advance(self, upload: Callable, step: Callable, synchronize: Callable) -> None:
+        """Upload the record, run one step and commit only after a successful sync.
+
+        Any exception after the upload is attempted leaves device geometry
+        unknown: the circle is marked failed, never reports applied state again,
+        and the owner's ``on_failure`` hook is called.
+        """
+        record = self.record()
+        try:
+            upload(record)
+            step()
+            synchronize()
+        except BaseException as exc:
+            self.fail(exc)
+            raise
+        self.commit()
+
+    def commit(self) -> None:
+        """Applied := pending offset (caller has synchronized a successful step)."""
+        if self.pending is not None:
+            self.applied = self.pending[0]
+            self.pending = None
+
+    def fail(self, exc: BaseException) -> None:
+        self.failed = True
+        if self._on_failure is not None:
+            self._on_failure(exc)
+
+    def band(self) -> Tuple[int, int, int, int]:
+        """Inclusive cell AABB (x0, x1, y0, y1) swept over the whole offset box."""
+        x0, x1, y0, y1 = self.translation.offset_box
+        r, m = self.radius, BAND_MARGIN
+        return (math.floor(self.base[0] + x0 - r) - m, math.ceil(self.base[0] + x1 + r) + m,
+                math.floor(self.base[1] + y0 - r) - m, math.ceil(self.base[1] + y1 + r) + m)
+
+    def _require_usable(self) -> None:
+        if self.failed:
+            raise RuntimeError("moving circle state is unknown after a failed step")
+
+
+def _inside(x: int, y: int, c: Vec, r: float) -> bool:
+    return (x - c[0]) ** 2 + (y - c[1]) ** 2 < r * r
+
+
+def empty_donor_cells(c_old: Vec, c_new: Vec, radius: float):
+    """Fresh cells (inside c_old, outside c_new) whose donor set V(x) is empty."""
+    r = radius
+    xs = range(math.floor(min(c_old[0], c_new[0]) - r) - 1, math.ceil(max(c_old[0], c_new[0]) + r) + 2)
+    ys = range(math.floor(min(c_old[1], c_new[1]) - r) - 1, math.ceil(max(c_old[1], c_new[1]) + r) + 2)
+    empty = []
+    for x in xs:
+        for y in ys:
+            if not _inside(x, y, c_old, r) or _inside(x, y, c_new, r):
+                continue
+            if not any(not _inside(x + cx, y + cy, c_old, r) and not _inside(x + cx, y + cy, c_new, r)
+                       for cx, cy in zip(_CX, _CY)):
+                empty.append((x, y))
+    return empty
+
+
+def enumerate_donor_coverage(radius: float, positions: int = 8, directions: int = 16) -> int:
+    """Count empty-donor fresh cells over sub-cell centres and moves of MAX_STEP.
+
+    Evidence for the declared bounds (radius >= R_MIN, |dc|_2 <= MAX_STEP),
+    not a proof; the device keeps the ``no_donor`` counter.
+    """
+    count = 0
+    for i in range(positions):
+        for j in range(positions):
+            c_old = (i / positions, j / positions)
+            for k in range(directions):
+                a = 2.0 * math.pi * k / directions
+                c_new = (c_old[0] + MAX_STEP * math.cos(a), c_old[1] + MAX_STEP * math.sin(a))
+                count += len(empty_donor_cells(c_old, c_new, radius))
+    return count
+
+
+class MovingTables:
+    """Device tables of the translating circles for ``MovingCircleCells``/``Links``.
+
+    Band table (x0, y0, width, height, first_slot, body_id), immutable per run;
+    FP64 motion records (cx_old, cy_old, cx_new, cy_new, r, ux, uy), immutable
+    during a step; separate per-body diagnostics (see DIAG), never part of obs.
+    """
+    BAND_INTS, GEOM_DOUBLES = 6, 7
+    DIAG = ("cme_x", "cme_y", "wen_x", "wen_y", "link_dm", "n_links",
+            "cover_jx", "cover_jy", "cover_mass", "uncover_jx", "uncover_jy", "uncover_mass",
+            "eq_only", "no_donor", "crossings")
+    DIAG_FLOATS = len(DIAG)
+
+    def __init__(self, band, geom):
+        import numpy as np
+        import pycuda.driver as cuda
+        self.n_bodies = len(geom) // self.GEOM_DOUBLES
+        self.n_slots = int(band[-2] + band[-4] * band[-3] * 8)
+        self.n_cells = self.n_slots // 8
+        self.geom = list(geom)  # host copy of the uploaded records
+        self.band_gpu = cuda.to_device(np.asarray(band, np.int32))
+        self.geom_gpu = cuda.to_device(np.asarray(geom, np.float64))
+        # Owned page-locked staging for async record uploads (allocated once).
+        self._geom_pinned = cuda.pagelocked_empty(len(geom), np.float64)
+        self.diag_gpu = cuda.to_device(np.zeros(self.n_bodies * self.DIAG_FLOATS, np.float32))
+        self.slot_debug_gpu = None  # optional (q, fallback, rx, ry) per slot for parity tests
+
+    @classmethod
+    def upload(cls, circles):
+        """``circles``: [(MovingCircle, body_id)] in body order."""
+        return cls(*cls.layout(circles))
+
+    @staticmethod
+    def layout(circles):
+        """Host band and geometry tables (no device work)."""
+        band, geom, first = [], [], 0
+        for motion, body_id in circles:
+            x0, x1, y0, y1 = motion.band()
+            width, height = x1 - x0 + 1, y1 - y0 + 1
+            band += [x0, y0, width, height, first, int(body_id)]
+            first += width * height * 8
+            geom += MovingTables.record(motion)
+        return band, geom
+
+    @staticmethod
+    def record(motion):
+        """FP64 motion record of one circle for its next step (hold if nothing pending)."""
+        c_old, c_new, u = motion.record()
+        return [*c_old, *c_new, motion.radius, *u]
+
+    def upload_records(self, circles, stream) -> bool:
+        """Enqueue the records when they differ from the device copy; True if advancing.
+
+        The pinned staging buffer is reused: callers synchronize ``stream`` after every
+        upload before the next one (``Simulation._run_moving``), so no copy is in flight
+        when it is rewritten. ``geom`` is the enqueued state; after a failed copy the
+        owner is failed and the cache is never trusted again.
+        """
+        import pycuda.driver as cuda
+        geom = [v for motion, _ in circles for v in self.record(motion)]
+        if geom != self.geom:
+            self._geom_pinned[:] = geom
+            cuda.memcpy_htod_async(self.geom_gpu, self._geom_pinned, stream)
+            self.geom = geom
+        g = self.GEOM_DOUBLES
+        return any(geom[i:i + 2] != geom[i + 2:i + 4] for i in range(0, len(geom), g))
+
+    def enable_slot_debug(self):
+        import numpy as np
+        import pycuda.driver as cuda
+        self.slot_debug_gpu = cuda.to_device(np.zeros(self.n_slots * 4, np.float32))
+
+    def zero_diagnostics(self, stream):
+        import pycuda.driver as cuda
+        cuda.memset_d32_async(self.diag_gpu, 0, self.n_bodies * self.DIAG_FLOATS, stream)
+
+    def download_diagnostics(self):
+        import numpy as np
+        import pycuda.driver as cuda
+        out = np.empty((self.n_bodies, self.DIAG_FLOATS), np.float32)
+        cuda.memcpy_dtoh(out, self.diag_gpu)
+        return out
@@ -140,6 +140,7 @@ class SimObject:
            x=float(c[0]), y=float(c[1]), z=z0,
        )
        self.control = ObjectControl()
+        self.motion = None  # MovingCircle for a translating circle (body/motion.py)

    @property
    def surface_actuation(self):
@@ -18,7 +18,9 @@ all NQ populations from an owner core cell into another block's halo cell.
 Periodicity is explicit per axis. Across a non-periodic exterior, halo cells
 have no pair and stay invalid (NaN after packing): no wrap, no implied BC.

-Preparation only: no CUDA, no boundary conditions, no refinement levels.
+Host preparation plus thin launch wrappers: no CUDA import or allocation
+(launch_* only call kernels and device tables the caller owns), no refinement
+levels; level_config and static_circle_tables prepare one blocked level.
 The one-cell same-level halo here is unrelated to the four-cell fine halo of
 the fixed refinement patch (refinement.py), whose semantics are unchanged.

@@ -122,7 +122,7 @@ def leaf_partition(h):
    for level, kind in enumerate(h.kinds):
        if level:
            pair = h.pairs[level - 1]
-            py, px = (pair.fine_level_offset if isinstance(pair, FixedPatch2D) else (0, 0))
+            py, px = pair.fine_level_offset
            oy, ox = 2 * oy + py, 2 * ox + px
        y, x = np.nonzero(kind == ACTIVE)
        np.add.at(base, ((y + oy) >> level, (x + ox) >> level), 0.25 ** level)
@@ -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,
@@ -0,0 +1,183 @@
+// Translating circles (docs/moving_bodies.md): M1 MovingCircleCells reclassifies
+// band cells and refills uncovered ones; M2 MovingCircleLinks applies the cut links.
+// Geometry comes from an immutable FP64 record per body
+// {cx_old, cy_old, cx_new, cy_new, r, ux, uy}; the M2 operation order mirrors
+// CircleGeometry.build_cut_links / find_circle_intersection with explicit _rn
+// operations (bit parity with the host builder is a test target, not a claim).
+// Moving circles are never in the packed static links, so this kernel is the
+// only writer of their solid slots.
+#ifndef CELERIS_MOVING_CIRCLE_CU
+#define CELERIS_MOVING_CIRCLE_CU
+
+#if DIM == 2
+constexpr int MOVING_BAND_INTS = 6;   // x0, y0, width, height, first_slot, body_id
+constexpr int MOVING_GEOM_DOUBLES = 7;
+// cme_x, cme_y, wen_x, wen_y, link_dm, n_links, cover_jx, cover_jy, cover_mass,
+// uncover_jx, uncover_jy, uncover_mass, eq_only, no_donor, crossings
+constexpr int MOVING_DIAG_FLOATS = 15;
+
+__device__ __forceinline__ double mc_sq(double v) { return __dmul_rn(v, v); }
+
+__device__ __forceinline__ bool mc_inside(int x, int y, double cx, double cy, double r)
+{
+    return __dadd_rn(mc_sq(__dsub_rn((double)x, cx)), mc_sq(__dsub_rn((double)y, cy))) < mc_sq(r);
+}
+
+// Smallest root t in [0, 1] of |p + t c - centre| = r (find_circle_intersection order).
+__device__ __forceinline__ bool mc_hit(int x, int y, int ex, int ey,
+                                       double cx, double cy, double r, double& t)
+{
+    const double dx = (double)ex, dy = (double)ey;
+    const double px = __dsub_rn((double)x, cx), py = __dsub_rn((double)y, cy);
+    const double a = __dadd_rn(mc_sq(dx), mc_sq(dy));
+    const double b = __dmul_rn(2.0, __dadd_rn(__dmul_rn(dx, px), __dmul_rn(dy, py)));
+    const double c = __dsub_rn(__dadd_rn(mc_sq(px), mc_sq(py)), mc_sq(r));
+    const double det = __dsub_rn(mc_sq(b), __dmul_rn(__dmul_rn(4.0, a), c));
+    if (det < 0.0) return false;
+    const double s = __dsqrt_rn(det), a2 = __dmul_rn(2.0, a);
+    const double t1 = __ddiv_rn(__dsub_rn(-b, s), a2);
+    const double t2 = __ddiv_rn(__dadd_rn(-b, s), a2);
+    if (0.0 <= t1 && t1 <= 1.0) { t = t1; return true; }
+    if (0.0 <= t2 && t2 <= 1.0) { t = t2; return true; }
+    return false;
+}
+
+__global__ void MovingCircleLinks(
+    fpxx* fi, const uint16_t* flag,
+    const int* band, const double* geom,
+    unsigned int n_bodies, unsigned int n_slots,
+    const float* action, float* obs, float* diag, float* slot_debug)
+{
+    const unsigned int tid = threadIdx.x + blockIdx.x * blockDim.x;
+    if (tid >= n_slots) return;
+    unsigned int m = 0;
+    while (m + 1 < n_bodies && tid >= (unsigned int)band[(m + 1) * MOVING_BAND_INTS + 4]) ++m;
+    const int* b = band + m * MOVING_BAND_INTS;
+    const unsigned int local = tid - (unsigned int)b[4];
+    const unsigned int dir = local % 8u + 1u;
+    const int x = b[0] + (int)((local / 8u) % (unsigned int)b[2]);
+    const int y = b[1] + (int)((local / 8u) / (unsigned int)b[2]);
+    const double* g = geom + MOVING_GEOM_DOUBLES * m;
+    const double cx = g[2], cy = g[3], r = g[4];
+    const int ex = d_cx[dir], ey = d_cy[dir];
+    float* debug = slot_debug == nullptr ? nullptr : slot_debug + 4 * tid;
+    if (debug != nullptr) debug[0] = -1.0f;
+
+    double t;
+    if (mc_inside(x, y, cx, cy, r) || !mc_inside(x + ex, y + ey, cx, cy, r)
+            || !mc_hit(x, y, ex, ey, cx, cy, r, t))
+        return;
+    const double hx = __dadd_rn((double)x, __dmul_rn(t, (double)ex));
+    const double hy = __dadd_rn((double)y, __dmul_rn(t, (double)ey));
+    const double hdx = __dsub_rn(hx, (double)x), hdy = __dsub_rn(hy, (double)y);
+    const double norm = __dsqrt_rn((double)(ex * ex + ey * ey));
+    double qd = __ddiv_rn(__dsqrt_rn(__dadd_rn(mc_sq(hdx), mc_sq(hdy))), norm);
+    qd = fmax(0.01, fmin(qd, 1.0));
+    const float q = (float)qd;
+    const float rx = (float)__dsub_rn(hx, cx), ry = (float)__dsub_rn(hy, cy);
+
+    // Same two-stage order as build_cut_links + pack_cut_links_to_soa: the fallback
+    // class is chosen from the clamped FP64 q (q < 1/2 needs a fluid donor at x - c);
+    // only then are q, rx, ry rounded to FP32 for the BFL branch and spin velocity.
+    unsigned char fallback = CURVED_FALLBACK_BOUZIDI;
+    if (qd < 0.5) {
+        const int xd = x - ex, yd = y - ey;
+        if (xd < 0 || xd >= (int)NX || yd < 0 || yd >= (int)NY || mc_inside(xd, yd, cx, cy, r)
+                || (flag[linear_index((unsigned int)xd, (unsigned int)yd)] & FLAG_FLUID) == 0)
+            fallback = CURVED_FALLBACK_HALFWAY;
+    }
+    if (debug != nullptr) { debug[0] = q; debug[1] = (float)fallback; debug[2] = rx; debug[3] = ry; }
+
+    const int bid = b[5];
+    const float omega = action_omega(action, bid);
+    const float Uw = (float)g[5] - omega * ry, Vw = (float)g[6] + omega * rx;  // u + omega x r
+    float exchange[3], dm;
+    apply_bouzidi_link(dir, q, linear_index((unsigned int)x, (unsigned int)y), fi,
+                       Uw, Vw, rx, ry, fallback, obs, bid, exchange, &dm);
+    if (diag != nullptr) {
+        float* d = diag + m * MOVING_DIAG_FLOATS;
+        atomicAdd(&d[0], exchange[0]);
+        atomicAdd(&d[1], exchange[1]);
+        atomicAdd(&d[2], Uw * dm);   // Wen 2014 correction -u_w (f_t - f_r)
+        atomicAdd(&d[3], Vw * dm);
+        atomicAdd(&d[4], dm);
+        atomicAdd(&d[5], 1.0f);
+    }
+}
+// M1: one thread per band cell, launched only on steps with c_old != c_new.
+// Membership comes from the immutable old/new predicates, never from flags, so
+// classification, refill and the flag write of one cell cannot race with another.
+// ddf holds post-collision G^{n-1}; refilled cells get all nine slots at that level.
+__global__ void MovingCircleCells(
+    fpxx* fi, uint16_t* flag, const int* band, const double* geom,
+    unsigned int n_bodies, unsigned int n_cells, const float* action, float* diag)
+{
+    const unsigned int tid = threadIdx.x + blockIdx.x * blockDim.x;
+    if (tid >= n_cells) return;
+    unsigned int m = 0;
+    while (m + 1 < n_bodies && tid >= (unsigned int)band[(m + 1) * MOVING_BAND_INTS + 4] / 8u) ++m;
+    const int* b = band + m * MOVING_BAND_INTS;
+    const unsigned int local = tid - (unsigned int)b[4] / 8u;
+    const int x = b[0] + (int)(local % (unsigned int)b[2]);
+    const int y = b[1] + (int)(local / (unsigned int)b[2]);
+    const double* g = geom + MOVING_GEOM_DOUBLES * m;
+    const bool in_old = mc_inside(x, y, g[0], g[1], g[4]);
+    const bool in_new = mc_inside(x, y, g[2], g[3], g[4]);
+    if (in_old == in_new) return;
+    const unsigned long k = linear_index((unsigned int)x, (unsigned int)y);
+    float* d = diag + m * MOVING_DIAG_FLOATS;
+    float rho = 0.0f, jx = 0.0f, jy = 0.0f;
+    if (in_new) {  // covered: its fluid content leaves the active set
+        for (unsigned int i = 0; i < NQ; ++i) {
+            const float f = load_ddf(fi, index_f(k, i));
+            rho += f; jx += (float)d_cx[i] * f; jy += (float)d_cy[i] * f;
+        }
+        flag[k] = FLAG_SOLID | FLAG_OBSTACLE | FLAG_BC_CURVED;
+        atomicAdd(&d[6], jx); atomicAdd(&d[7], jy); atomicAdd(&d[8], rho);
+        atomicAdd(&d[14], 1.0f);
+        return;
+    }
+    // Uncovered: donors are fluid in BOTH geometries; NEQ from the neighbour with the
+    // largest unnormalized c.n > 0 (first direction wins ties), n from c_new.
+    const float nx = (float)__dsub_rn((double)x, g[2]), ny = (float)__dsub_rn((double)y, g[3]);
+    float rho_sum = 0.0f, best_dot = 0.0f;
+    int n_valid = 0, best = -1;
+    for (unsigned int i = 1; i < NQ; ++i) {
+        const int xn = x + d_cx[i], yn = y + d_cy[i];
+        if (mc_inside(xn, yn, g[0], g[1], g[4]) || mc_inside(xn, yn, g[2], g[3], g[4])) continue;
+        const unsigned long kn = linear_index((unsigned int)xn, (unsigned int)yn);
+        for (unsigned int q = 0; q < NQ; ++q) rho_sum += load_ddf(fi, index_f(kn, q));
+        ++n_valid;
+        const float dot = (float)d_cx[i] * nx + (float)d_cy[i] * ny;
+        if (dot > best_dot) { best_dot = dot; best = (int)i; }
+    }
+    const float omega = action_omega(action, b[5]);
+    const float uwx = (float)g[5] - omega * ny, uwy = (float)g[6] + omega * nx;  // at the cell
+    float feq[NQ], neq[NQ];
+    compute_feq(n_valid > 0 ? rho_sum / (float)n_valid : (float)RHO, uwx, uwy, feq);
+    for (unsigned int i = 0; i < NQ; ++i) neq[i] = 0.0f;
+    if (best >= 0) {
+        const unsigned long kd = linear_index((unsigned int)(x + d_cx[best]), (unsigned int)(y + d_cy[best]));
+        float gd[NQ], rd = 0.0f, mx = 0.0f, my = 0.0f, fd[NQ];
+        for (unsigned int i = 0; i < NQ; ++i) {
+            gd[i] = load_ddf(fi, index_f(kd, i));
+            rd += gd[i]; mx += (float)d_cx[i] * gd[i]; my += (float)d_cy[i] * gd[i];
+        }
+        compute_feq(rd, mx / rd, my / rd, fd);
+        for (unsigned int i = 0; i < NQ; ++i) neq[i] = gd[i] - fd[i];
+    } else {
+        atomicAdd(&d[12], 1.0f);
+    }
+    if (n_valid == 0) atomicAdd(&d[13], 1.0f);  // hard failure, checked before commit
+    for (unsigned int i = 0; i < NQ; ++i) {
+        const float f = feq[i] + neq[i];
+        store_ddf(fi, index_f(k, i), f);
+        rho += f; jx += (float)d_cx[i] * f; jy += (float)d_cy[i] * f;
+    }
+    flag[k] = FLAG_FLUID;
+    atomicAdd(&d[9], jx); atomicAdd(&d[10], jy); atomicAdd(&d[11], rho);
+    atomicAdd(&d[14], 1.0f);
+}
+#endif  // DIM == 2
+
+#endif  // CELERIS_MOVING_CIRCLE_CU
@@ -71,6 +71,7 @@ extern "C"
 #include "step/one_step_double.cu"
 #include "step/one_step_esopull.cu"
 #include "step/aux_kernels.cu"
+#include "step/moving_circle.cu"
 #include "step/esopull_macro.cu"
 #include "step/esopull_sync.cu"

@@ -120,6 +120,11 @@ class FixedPatch2D:
        (an index offset, not an (x, y) origin)."""
        return 2 * self.bounds[1] - HALO, 2 * self.bounds[0] - HALO

+    @property
+    def rectangles(self) -> tuple:
+        """The one parent rectangle, as a sequence like FixedPatchSet2D.rectangles."""
+        return (self.bounds,)
+
    @property
    def covered(self) -> np.ndarray:
        """Dense coarse mask of cells replaced by the patch (not leaves)."""
@@ -225,11 +230,18 @@ class FixedPatchSet2D:
    def bounds(self) -> tuple:
        return tuple(p.bounds for p in self.patches)

+    rectangles = bounds
+
    @property
    def fine_origin(self) -> tuple[float, float]:
        """Lower edge of the level-wide fine grid in coarse physical units."""
        return 0.0, 0.0

+    @property
+    def fine_level_offset(self) -> tuple[int, int]:
+        """The level-wide grid itself: zero (y, x) index offset."""
+        return 0, 0
+
    @property
    def covered(self) -> np.ndarray:
        mask = np.zeros(self.coarse_shape, bool)
@@ -111,6 +111,57 @@ def _check_levels(refinements, blocks, objects, jets):
            raise ValueError("bodies/Jets are supported only on the finest level of a nested hierarchy")


+# Host preparation (numpy only, before any device work); each returns one entry per level.
+
+def _level_layouts(cfg, pairs, level_blocks, wall_aware):
+    """Blocked storage layout of every level, checked against the launch limits."""
+    layouts = []
+    for level, block in enumerate(level_blocks):
+        if level == 0:
+            layout = BlockLayout2D.build(cfg.ny, cfg.nx, *block, periodic=(True, not wall_aware))
+        elif isinstance(pairs[level - 1], FixedPatchSet2D):
+            layout = pairs[level - 1].fine_layout(block)
+        else:
+            layout = BlockLayout2D.build(*pairs[level - 1].fine_shape, *block,
+                                         periodic=(False, False))
+        if layout.nblocks > MAX_GRID_DIM or layout.local_shape[0] > MAX_GRID_DIM:
+            raise ValueError(f"level {level} blocks exceed launch limits")
+        layouts.append(layout)
+    return layouts
+
+
+def _host_boundaries(hierarchy, level_cfgs, level_objects, level_jets, wall_aware, shifted):
+    """Host boundary geometry of every refined level in the root edge frame (None at L0)."""
+    hosts = [None]
+    for level in range(1, hierarchy.depth):
+        pair = hierarchy.pairs[level - 1]
+        if level == 1:  # the root config carries the wall-aware boundary fields
+            host = prepare_fine_boundary(level_cfgs[0], pair, level_objects[1],
+                                         wall_aware=wall_aware, jets=level_jets[1],
+                                         _shifted_blocked=shifted)
+        else:
+            origin, spacing = hierarchy.edge_frame(level - 1)
+            host = prepare_fine_boundary(level_cfgs[level - 1], pair, level_objects[level],
+                                         jets=level_jets[level], parent_origin=origin,
+                                         parent_spacing=spacing, _shifted_blocked=shifted)
+        hosts.append(host)
+    return hosts
+
+
+def _fluid_leaves(hierarchy, boundaries):
+    """ACTIVE cells of every level that no static circle or body makes solid."""
+    leaves = [kind == ACTIVE for kind in hierarchy.kinds]
+    for level, b in enumerate(boundaries):
+        if b is None:
+            continue
+        if b.host is None:  # single level: static root circles
+            leaves[level] &= ~b.tables["solid"]
+        else:
+            fine_shape = hierarchy.pairs[level - 1].fine_shape
+            leaves[level] &= (b.host.flags.reshape(fine_shape) & SOLID) == 0
+    return leaves
+
+
 class FixedPatchStepper:
    """D2Q9 SRT/MRT FP32 fixed 2:1 refinement hierarchy, built from level lists.

@@ -246,40 +297,15 @@ class FixedPatchStepper:
        self.hierarchy = build_hierarchy(cfg, refinements, wall_aware=self._wall_aware)
        pairs = self.hierarchy.pairs
        if self._wall_aware:
-            rects = (pairs[0].bounds,) if isinstance(pairs[0], FixedPatch2D) else pairs[0].bounds
-            if any(r[1] < 4 or r[3] > cfg.ny - 4 for r in rects):
+            if any(r[1] < 4 or r[3] > cfg.ny - 4 for r in pairs[0].rectangles):
                raise ValueError("slot_reservoir requires y0 >= 4 and y1 <= ny-4")
-        self._layouts = []
-        for level, block in enumerate(level_blocks):
-            if level == 0:
-                layout = BlockLayout2D.build(cfg.ny, cfg.nx, *block,
-                                             periodic=(True, not self._wall_aware))
-            elif isinstance(pairs[level - 1], FixedPatchSet2D):
-                layout = pairs[level - 1].fine_layout(block)
-            else:
-                layout = BlockLayout2D.build(*pairs[level - 1].fine_shape, *block,
-                                             periodic=(False, False))
-            if layout.nblocks > MAX_GRID_DIM or layout.local_shape[0] > MAX_GRID_DIM:
-                raise ValueError(f"level {level} blocks exceed launch limits")
-            self._layouts.append(layout)
+        self._layouts = _level_layouts(cfg, pairs, level_blocks, self._wall_aware)
        if self._channel:
            base = blocked_coarse_tables(pairs[0], self._layouts[0], np.argsort(pairs[0].coalesce_coarse))
            validate_channel_transfers(pairs[0], self._layouts[0], base[1])
        self._level_cfgs = _kernel_configs(cfg, pairs, self._wall_aware)
-        # Host boundary geometry per level in the root edge frame (Hierarchy.edge_frame).
-        self._boundaries = [None]
-        for level in range(1, depth):
-            pair = pairs[level - 1]
-            if level == 1:  # the root config carries the wall-aware boundary fields
-                host = prepare_fine_boundary(cfg, pair, level_objects[1],
-                                             wall_aware=self._wall_aware, jets=level_jets[1],
-                                             _shifted_blocked=self._ddf_shifting)
-            else:
-                origin, spacing = self.hierarchy.edge_frame(level - 1)
-                host = prepare_fine_boundary(self._level_cfgs[level - 1], pair, level_objects[level],
-                                             jets=level_jets[level], parent_origin=origin,
-                                             parent_spacing=spacing, _shifted_blocked=self._ddf_shifting)
-            self._boundaries.append(host)
+        hosts = _host_boundaries(self.hierarchy, self._level_cfgs, level_objects, level_jets,
+                                 self._wall_aware, self._ddf_shifting)
        if not isinstance(observe_raw_exchange, bool):
            raise ValueError("observe_raw_exchange must be bool")
        self._observe_raw_exchange = observe_raw_exchange
@@ -291,32 +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)
-        self._level_jets = level_jets
        root = None
        if depth == 1 and level_objects[0]:
            tables = static_circle_tables(self._layouts[0], level_objects[0])
            root = _BlockedBoundary(None, tables, (), False)
-        self._device_boundaries = [root] + [
-            _prepare_blocked_boundary(self._boundaries[level], self._layouts[level],
+        # The one boundary/Jet record per level; host geometry and Jets are views of it.
+        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 = [kind == ACTIVE for kind in self.hierarchy.kinds]
-        if root is not None:
-            self._leaf[0] &= ~root.tables["solid"]
-        for level in range(1, depth):
-            host = self._boundaries[level]
-            if host is not None:
-                self._leaf[level] &= (host.flags.reshape(pairs[level - 1].fine_shape) & SOLID) == 0
+        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"
@@ -371,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)
@@ -394,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)
@@ -556,19 +575,17 @@ 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 self._boundaries[1] is not None
-                    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,
                                           _shifted_blocked=self._ddf_shifting)
-            if (not same_physical_boundary(old.boundary, self._boundaries[2])
-                    or (self._boundaries[2] is not None and actual is None)):
+            if not same_physical_boundary(old.boundary, None if actual is None else actual.host):
                raise ValueError("objects must reproduce the original root-edge-frame level-2 boundary")
            if actual is not None:
                velocity = actual.tables.get("surface_velocity")
@@ -585,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():
@@ -603,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)
@@ -622,7 +639,7 @@ class FixedPatchStepper:
            return candidate

    def set_jet_amplitude(self, index, amplitude, *, level):
-        """Replace Jet ``index`` of ``level_jets[level]`` between steps; no geometry/table rebuild.
+        """Replace Jet ``index`` of level ``level`` (``jets[level]``) between steps; no table rebuild.

        Invalid input raises before any change. A device/sync error marks the
        instance failed (as step/initialize do); host state commits only after
@@ -633,27 +650,24 @@ 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]
-            jets = list(self._level_jets[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
-            velocity = None
-            if boundary is not None and boundary.jet_cache is not None:
-                shape, jet_id = boundary.jet_cache
-                velocity = boundary.tables["surface_velocity"].copy()
-                driven = jet_id >= 0
-                velocity[driven] = jet_link_velocity(jets, shape, jet_id)[driven]
+            shape, jet_id = boundary.jet_cache
+            velocity = boundary.tables["surface_velocity"].copy()
+            driven = jet_id >= 0
+            velocity[driven] = jet_link_velocity(jets, shape, jet_id)[driven]
+            # Device action: an uncertain write fails the instance closed.
            try:
-                if velocity is not None:
-                    self._stream.synchronize()  # no in-flight step reads the column
-                    cuda.memcpy_htod(boundary.device["surface_velocity"], velocity)
+                self._stream.synchronize()  # no in-flight step reads the column
+                cuda.memcpy_htod(boundary.device["surface_velocity"], velocity)
            except BaseException:
                self._status = "failed"
                raise
-            if velocity is not None:
-                boundary.tables["surface_velocity"] = velocity
-            if boundary is not None:
-                boundary.jets = tuple(jets)
-            self._level_jets[level] = tuple(jets)
+            # Commit: the one boundary record holds both the column and the Jets.
+            boundary.tables["surface_velocity"] = velocity
+            boundary.jets = tuple(jets)

    def _blocked_boundary(self, level, buffer, **kw):
        b = level.boundary
@@ -696,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)
@@ -740,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)
@@ -763,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:
@@ -776,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")
@@ -784,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."""
@@ -873,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
@@ -901,9 +915,20 @@ class FixedPatchStepper:
                raise

    # Read-only per-level state, indexed by level (0 = root). Pairs: hierarchy.pairs[l-1].
-    layouts = property(lambda self: tuple(self._layouts))
-    boundaries = property(lambda self: tuple(self._boundaries),
-                          doc="Host boundary geometry per level (None: no bodies/Jets).")
-    jets = property(lambda self: tuple(self._level_jets))
-    leaf = property(lambda self: tuple(m.copy() for m in self._leaf),
-                    doc="Per-level (ny, nx) masks of public (non-NaN) download cells.")
+    @property
+    def layouts(self):
+        return tuple(self._layouts)
+
+    @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._boundaries)
+
+    @property
+    def jets(self):
+        return tuple(() if b is None else b.jets for b in self._boundaries)
+
+    @property
+    def leaf(self):
+        """Per-level (ny, nx) masks of public (non-NaN) download cells."""
+        return tuple(m.copy() for m in self._leaf)
@@ -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):
@@ -20,6 +20,10 @@ class LBMStepper:
        self.cfg = cfg

        self._esopull = (cfg.streaming == "esopull")
+        # Set by the owner when a step would be wrong (translating circles
+        # excluded from the packed links before the moving kernels exist).
+        self.blocked: str | None = None
+        self.moving = None  # MovingTables attached by the owner (translating circles)

        if self._esopull:
            self.step_fn = module.get_function("EsoPullStep")
@@ -35,6 +39,8 @@ class LBMStepper:
        self.apply_esopull_curved_fn = module.get_function("ApplyEsoPullCurvedKernel")
        self.sensor_fn = module.get_function("SensorKernel")
        self.force_region_fn = module.get_function("ForceRegionKernel")
+        self.moving_fn = None if self._esopull else module.get_function("MovingCircleLinks")
+        self.moving_cells_fn = None if self._esopull else module.get_function("MovingCircleCells")

        tpb = cfg.threads_per_block
        self.block = (tpb, 1, 1)
@@ -70,6 +76,13 @@ class LBMStepper:
        stream: cuda.Stream | None = None,
    ):
        """Advance *n* time steps."""
+        if n > 0 and self.blocked is not None:
+            raise RuntimeError(self.blocked)
+        self._step(n, action_gpu=action_gpu, obs_gpu=obs_gpu, stream=stream)
+
+    def _step(self, n, *, action_gpu, obs_gpu, stream=None, advance=False):
+        """Unguarded steps; ``advance`` runs one step that moves translating circles."""
+        assert not advance or n == 1
        f = self.field
        launch_kw = {}
        if stream is not None:
@@ -81,6 +94,8 @@ class LBMStepper:
                    f.ddf_gpu, self.boundary_values_gpu, np.uint64(self._step_count),
                    block=self.block, grid=self.boundary_grid, **launch_kw,
                )
+            if self.moving is not None:
+                self._launch_moving(action_gpu, obs_gpu, advance, **launch_kw)
            self._launch_curved(action_gpu, obs_gpu, **launch_kw)

            if self._esopull:
@@ -106,6 +121,21 @@ class LBMStepper:
            # Any completed GPU step invalidates host DDF mirror.
            f.invalidate_host_ddf_cache()

+    def _launch_moving(self, action_gpu, obs_gpu, advance, **launch_kw):
+        m, tpb = self.moving, self.cfg.threads_per_block
+        if advance:  # M1 completes before M2 on the same stream
+            self.moving_cells_fn(
+                self.field.ddf_gpu, self.field.flag_gpu, m.band_gpu, m.geom_gpu,
+                np.uint32(m.n_bodies), np.uint32(m.n_cells), action_gpu, m.diag_gpu,
+                block=(tpb, 1, 1), grid=((m.n_cells + tpb - 1) // tpb, 1, 1), **launch_kw,
+            )
+        self.moving_fn(
+            self.field.ddf_gpu, self.field.flag_gpu, m.band_gpu, m.geom_gpu,
+            np.uint32(m.n_bodies), np.uint32(m.n_slots), action_gpu, obs_gpu,
+            m.diag_gpu, m.slot_debug_gpu if m.slot_debug_gpu is not None else np.uintp(0),
+            block=(tpb, 1, 1), grid=((m.n_slots + tpb - 1) // tpb, 1, 1), **launch_kw,
+        )
+
    def _launch_curved(self, action_gpu, obs_gpu, **launch_kw):
        f = self.field
        if f.n_curved == 0:
@@ -26,6 +26,13 @@ import pycuda.driver as cuda
 from .config import LBMConfig, BodyConfig, load_lbm_config, load_body_config
 from .cuda.context import CudaContext, context_bound
 from .cuda import compiler
+from .lbm.descriptors import FLUID
+
+# Translating circles step only through their moving tables (fixed centre in S2a).
+PENDING_BLOCKED = "a pending translation target is advanced only by Simulation.run(1)"
+MOVING_BLOCKED = ("translating circles need the moving kernels and their tables (attached "
+                  "at initialize); "
+                  "positive-step runs are rejected")
 from .lbm.field import LBMField
 from .lbm.stepper import LBMStepper
 from .body.objects import SimObject
@@ -97,7 +104,8 @@ class Simulation:

    # -- Object management ---------------------------------------------------
    def add_body(self, type: str = "circle", center=None,
-                 radius: float = 0.0, *, surface_actuation=None, **kwargs) -> int:
+                 radius: float = 0.0, *, surface_actuation=None, translation=None,
+                 **kwargs) -> int:
        """Add a simulation body. Returns body_id.

        Before ``initialize()`` the object is added to the formal registry
@@ -113,6 +121,9 @@ class Simulation:
                angle, half_width (radians), normal_velocity (positive outward).
                D2Q9/FP32/unshifted/DB only; set before initialization or stage
                a new body. No runtime motion or amplitude setter.
+            translation: ``body.motion.Translation`` declaring a translating
+                circle (before initialize only; docs/moving_bodies.md). Until
+                the moving kernels exist, positive-step runs are rejected.
            **kwargs: passed to the geometry constructor.

        Returns:
@@ -123,6 +134,10 @@ class Simulation:
        type_lower = str(type).strip().lower()
        if surface_actuation is not None and type_lower != "circle":
            raise ValueError("Surface actuation is only supported on circles")
+        if translation is not None and (type_lower != "circle" or surface_actuation is not None
+                                        or self._initialized):
+            raise ValueError("translation needs a circle without surface actuation, "
+                             "added before initialize()")

        if type_lower == "circle":
            from .body.geometry.circle import CircleGeometry
@@ -131,6 +146,11 @@ class Simulation:
            obj = SimObject(obj_id=-1, geometry=geom,
                            center=center, radius=radius,
                            is_sensor=False, surface_actuation=surface_actuation)
+            if translation is not None:
+                from .body.motion import MovingCircle
+                obj.motion = MovingCircle((cx, cy), radius, translation,
+                                          on_failure=self._fail_moving)
+                self.stepper.blocked = MOVING_BLOCKED
        elif type_lower == "sensor":
            from .body.geometry.circle import CircleGeometry
            cx, cy = float(center[0]), float(center[1])
@@ -190,6 +210,78 @@ class Simulation:
                raise
        self.bodies.jets[index] = jet  # later geometry syncs keep the new amplitude

+    def _fail_moving(self, exc: BaseException) -> None:
+        # Device geometry is unknown after a failed moving step: same lifecycle as Jets.
+        self._failed = f"moving circle step failed: {exc!r}"
+        self.stepper.blocked = self._failed  # direct stepper calls are refused as well
+
+    def _attach_moving(self) -> None:
+        # Direct stepper calls stay blocked until the moving tables exist.
+        tables = getattr(self, "_moving_tables", None)
+        self.stepper.moving = tables
+        moving = self._moving_bodies()
+        if moving and tables is None:
+            self.stepper.blocked = MOVING_BLOCKED
+        elif any(obj.motion.pending is not None for obj in moving):
+            self.stepper.blocked = PENDING_BLOCKED  # direct steps never consume a target
+        else:
+            self.stepper.blocked = None
+
+    def _build_moving_tables(self):
+        moving = self._moving_bodies()
+        if not moving:
+            return None
+        from .body.motion import MovingTables
+        return MovingTables.upload([(obj.motion, obj.obj_id) for obj in moving])
+
+    @context_bound
+    def moving_diagnostics(self) -> np.ndarray:
+        """Per translating circle (in body order), columns ``MovingTables.DIAG``,
+        accumulated since obs were last zeroed (diagnostics only, not a force RHS)."""
+        self._require_usable()
+        tables = self.stepper.moving
+        if tables is None:
+            return np.zeros((0, 6), np.float32)
+        self.stream.synchronize()
+        return tables.download_diagnostics()
+
+    def _moving_bodies(self):
+        return [obj for obj in self.bodies.objects if obj.motion is not None]
+
+    def _reject_moving_bodies(self, operation: str) -> None:
+        if self._moving_bodies():
+            raise ValueError(f"{operation} is not supported with translating circles")
+
+    def _validate_moving_bodies(self) -> None:
+        """Reject unsupported modes and swept-band placement before any mutation."""
+        moving = self._moving_bodies()
+        if not moving:
+            return
+        cfg = self.lbm_cfg
+        if (not cfg.is_d2q9 or cfg.streaming != "double_buffer" or cfg.ddf_shifting
+                or cfg.store_precision != "FP32"):
+            raise ValueError("translating circles need D2Q9/FP32/unshifted/double_buffer")
+        if self.bodies.jets:
+            raise ValueError("Jets are not supported together with translating circles")
+        nx, ny = self.field.nx, self.field.ny
+        base = self.field.build_channel_flags().reshape(ny, nx)
+        for obj in moving:
+            x0, x1, y0, y1 = obj.motion.band()
+            if x0 < 1 or y0 < 1 or x1 > nx - 2 or y1 > ny - 2:
+                raise ValueError(f"swept band {(x0, x1, y0, y1)} touches the domain edge")
+            if not np.all((base[y0:y1 + 1, x0:x1 + 1] & FLUID) != 0):
+                raise ValueError(f"swept band {(x0, x1, y0, y1)} contains non-fluid base cells")
+            for other in self.bodies.objects:
+                if other is obj:
+                    continue
+                (cx, cy), r = other.center[:2], float(other.radius)
+                if other.motion is not None:
+                    ox0, ox1, oy0, oy1 = other.motion.band()
+                else:
+                    ox0, ox1, oy0, oy1 = cx - r - 2, cx + r + 2, cy - r - 2, cy + r + 2
+                if ox0 <= x1 and x0 <= ox1 and oy0 <= y1 and y0 <= oy1:
+                    raise ValueError(f"swept band {(x0, x1, y0, y1)} overlaps body {other.obj_id}")
+
    def _require_usable(self):
        failed = getattr(self, "_failed", None)
        if failed is not None:
@@ -260,7 +352,7 @@ class Simulation:

    # -- Runtime control -----------------------------------------------------
    def set_body(self, id: int, omega: float = None,
-                 vx: float = None, vy: float = None) -> None:
+                 vx: float = None, vy: float = None, offset=None) -> None:
        """Set runtime body state (host action array only, no H2D).

        The GPU action buffer is updated at the next ``run()`` call.
@@ -270,8 +362,21 @@ class Simulation:
            omega: rotation speed (lattice units).  ``None`` = no change.
            vx: translation velocity x (future).  ``None`` = no change.
            vy: translation velocity y (future).  ``None`` = no change.
+            offset: target centre offset of a translating circle; with vx/vy
+                (default 0) it replaces the pending target after validation
+                against the applied offset (docs/moving_bodies.md, 4a). Only
+                ``run(1)`` advances it; direct stepper calls are blocked meanwhile.
        """
        obj = self.bodies.get(id)
+        motion = obj.motion
+        if motion is not None:
+            if offset is None and (vx is not None or vy is not None):
+                raise ValueError("a translating circle takes vx/vy only with an offset")
+            target = None if offset is None else motion.check_target(
+                offset, (0.0 if vx is None else vx, 0.0 if vy is None else vy))
+            vx = vy = None  # the motion target, not ObjectState, owns translation
+        elif offset is not None:
+            raise ValueError(f"body {id} was not declared with a translation")
        proposed = dict(vars(obj.state))
        proposed.update({k: float(v) for k, v in (("omega", omega), ("vx", vx), ("vy", vy))
                         if v is not None})
@@ -287,6 +392,9 @@ class Simulation:
        if vy is not None:
            obj.state.vy = float(vy)
            changed = True
+        if motion is not None and target is not None:
+            motion.pending = target
+            self._attach_moving()
        if changed:
            self.bodies._refresh_action_from_objects()

@@ -443,6 +551,7 @@ class Simulation:
            self.field, self._module, self.lbm_cfg,
        )
        self.stepper._step_count = _prev_step_count
+        self._attach_moving()

        self._assert_object_count_contract(expected_count=expected_n_objects)

@@ -468,6 +577,8 @@ class Simulation:
        # Validate the complete geometry before flags, populations or modules
        # change. Compact-list construction publishes no outputs on failure.
        next_flags = self.bodies.build_flags(self.field.build_channel_flags())
+        self._validate_moving_bodies()
+        self._attach_moving()
        if self.bodies.jets and self.bodies.count == 0:
            raise ValueError("Jets need bodies with cut links")
        self.bodies.build_compact_lists_for(self.bodies.objects, next_flags)
@@ -484,6 +595,8 @@ class Simulation:
            self.bodies.sync_to_gpu(self.field, rebuild_flags=False)
        else:
            self.bodies.ensure_packed_buffers()
+        self._moving_tables = self._build_moving_tables()
+        self._attach_moving()
        self._assert_runtime_contracts()
        self._initialized = True

@@ -499,6 +612,7 @@ class Simulation:
        self._require_usable()
        if not self._initialized:
            raise RuntimeError("Call initialize() before sync_bodies()")
+        self._reject_moving_bodies("sync_bodies()")

        self._assert_allocation_contract()
        # 1. Synchronize GPU
@@ -599,6 +713,13 @@ class Simulation:
        checkpoint_interval = operator.index(checkpoint_interval)
        if steps < 0 or checkpoint_interval < 0:
            raise ValueError("steps and checkpoint_interval must be nonnegative")
+        moving = self._moving_bodies()
+        if steps > 0 and moving and self.stepper.moving is None:
+            raise RuntimeError(MOVING_BLOCKED)
+        if moving and checkpoint_interval > 0:
+            raise ValueError("checkpoints are not supported with translating circles")
+        if steps > 1 and any(obj.motion.pending is not None for obj in moving):
+            raise ValueError("a pending translation target needs run(1)")
        self._assert_runtime_contracts()
        # Discard any uncommitted body edits before stepping.
        if self.bodies.has_pending_edit():
@@ -615,6 +736,8 @@ class Simulation:
        # Zero obs segments before step group
        if zero_obs:
            self.bodies.zero_obs_async(stream)
+            if self.stepper.moving is not None:
+                self.stepper.moving.zero_diagnostics(stream)

        if checkpoint_interval > 0:
            done = 0
@@ -631,6 +754,9 @@ class Simulation:
                if done < steps or done % checkpoint_interval == 0:
                    stream.synchronize()
                    self.save_checkpoint()
+        elif self.stepper.moving is not None and steps > 0:
+            self._run_moving(steps, stream, sync_obs)
+            return
        else:
            self.stepper.step(
                steps,
@@ -649,6 +775,34 @@ class Simulation:

        stream.synchronize()

+    def _run_moving(self, steps, stream, sync_obs):
+        """Steps with translating circles: hold, or one advancing step that commits
+        each pending target only after synchronization and a zero ``no_donor`` count
+        (docs/moving_bodies.md 4a, 5). Any failure leaves the Simulation unusable."""
+        tables = self.stepper.moving
+        circles = [(obj.motion, obj.obj_id) for obj in self._moving_bodies()]
+        try:
+            advancing = tables.upload_records(circles, stream)
+            self.stepper._step(steps, action_gpu=self.bodies.action_gpu,
+                               obs_gpu=self.bodies.obs_gpu, stream=stream, advance=advancing)
+            self.bodies.increment_obs_steps(steps)
+            if sync_obs:
+                self.bodies.download_obs_full_async(stream)
+            stream.synchronize()
+            if advancing:
+                no_donor = tables.download_diagnostics()[:, tables.DIAG.index("no_donor")]
+                if np.any(no_donor > 0):
+                    raise RuntimeError("translating circle refill found no donor cell")
+                for motion, _ in circles:
+                    motion.commit()
+                tables.upload_records(circles, stream)  # first hold: old = new = applied, u = 0
+                stream.synchronize()
+                self._attach_moving()
+        except BaseException as exc:
+            for motion, _ in circles:  # device geometry/state unknown after any failure
+                motion.fail(exc)
+            raise
+
    def step(self, n: int = 1):
        """Advance *n* steps (convenience for interactive use)."""
        self.run(n)
@@ -666,8 +820,16 @@ class Simulation:
        self.field.download_ddf(step_id=self.stepper.step_count)
        return self.field.ddf.copy()

+    @context_bound
    def get_flags(self) -> np.ndarray:
-        return self.field.flag.copy()
+        """Flag array; read back from the device when translating circles may have moved."""
+        if getattr(self.stepper, "moving", None) is None:
+            return self.field.flag.copy()
+        self._require_usable()
+        self.stream.synchronize()
+        flags = np.empty_like(self.field.flag)
+        cuda.memcpy_dtoh(flags, self.field.flag_gpu)
+        return flags

    # -- Runtime parameter updates -------------------------------------------
    @context_bound
@@ -698,6 +860,7 @@ class Simulation:
        """Save full simulation state to an HDF5 file."""
        self._require_usable()
        from .common.checkpoint import save_checkpoint
+        self._reject_moving_bodies("save_checkpoint()")
        if self.bodies.jets:
            raise NotImplementedError("Jets are not checkpointed yet")
        return save_checkpoint(self.field, self.stepper, self.lbm_cfg,
@@ -708,6 +871,7 @@ class Simulation:
        """Restore simulation state from an HDF5 checkpoint."""
        self._require_usable()
        from .common.checkpoint import load_checkpoint
+        self._reject_moving_bodies("load_checkpoint()")
        if not self._initialized:
            self.initialize()
        load_checkpoint(path, self.field, self.stepper, self.lbm_cfg,
@@ -726,6 +890,7 @@ class Simulation:
                cuda.Context.synchronize()
                # Release owned GPU resources while their context is current.
                self.stepper = None
+                self._moving_tables = None
                self.bodies = None
                self.field = None
                self._module = None
@@ -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
@@ -0,0 +1,536 @@
+"""Translating-circle host contract, S1 (docs/moving_bodies.md 4a, 5, 7).
+
+Host doubles only: no kernels exist yet, so every positive-step run with a
+moving body must be rejected before any device write.
+"""
+from types import SimpleNamespace
+
+import pytest
+
+from CelerisLab.body.motion import (MAX_STEP, R_MIN, MovingCircle, Translation,
+                                    enumerate_donor_coverage)
+
+BAD = (float("nan"), float("inf"))
+
+
+def circle(box=(-2., 2., -3., 3.), speed=.1):
+    return MovingCircle((50., 40.), 10., Translation(box, speed))
+
+
+def test_translation_and_radius_validation():
+    for box in ((0., 1., -1., 1.)[:3], (.5, 1., -1., 1.), (-1., 1., BAD[0], 1.)):
+        with pytest.raises(ValueError):
+            Translation(box, .1)
+    for speed in (0., MAX_STEP + 1e-9) + BAD:
+        with pytest.raises(ValueError):
+            Translation((0., 0., -1., 1.), speed)
+    with pytest.raises(ValueError):
+        MovingCircle((0., 0.), R_MIN - .1, Translation((0., 0., -1., 1.), .1))
+
+
+def test_targets_validate_against_applied_and_last_wins():
+    c = circle()
+    c.set_target((0., .3), (0., .05))
+    c.set_target((0., .5), (0., .05))  # second setter replaces, checked vs applied (0, 0)
+    assert c.pending == ((0., .5), (0., .05))
+    for offset, velocity in (((0., .51), (0., 0.)),           # Euclidean bound
+                             ((.36, .36), (0., 0.)),          # |dc|_2 = .509
+                             ((0., 3.1), (0., 0.)),           # outside box
+                             ((0., .1), (0., .11)),           # speed bound
+                             ((0., 0.), (0., .01)),           # hold with sliding velocity
+                             ((BAD[0], 0.), (0., 0.)), ((0., 0.), (0., BAD[1]))):
+        with pytest.raises(ValueError):
+            c.set_target(offset, velocity)
+        assert c.pending == ((0., .5), (0., .05))  # invalid replacement keeps the old target
+    c.set_target((.35, .35), (0., 0.))  # |dc|_2 = .495 inside the bound
+    assert c.applied == (0., 0.)
+
+
+def test_record_commit_and_hold():
+    c = circle()
+    assert c.record() == ((50., 40.), (50., 40.), (0., 0.))  # hold: no translation velocity
+    c.set_target((0., .25), (0., .05))
+    assert c.record() == ((50., 40.), (50., 40.25), (0., .05))
+    uploads = []
+    c.advance(uploads.append, lambda: None, lambda: None)
+    assert uploads == [((50., 40.), (50., 40.25), (0., .05))]
+    assert c.applied == (0., .25) and c.pending is None  # offset-only commit
+    c.advance(uploads.append, lambda: None, lambda: None)  # no new target: hold step
+    assert uploads[-1] == ((50., 40.25), (50., 40.25), (0., 0.)) and c.applied == (0., .25)
+    c.set_target((0., 0.), (0., -.05))  # return to zero offset is a real target
+    assert c.pending == ((0., 0.), (0., -.05))
+    c.advance(uploads.append, lambda: None, lambda: None)
+    assert c.applied == (0., 0.) and c.pending is None
+
+
+@pytest.mark.parametrize("stage", ("upload", "step", "synchronize"))
+def test_failure_after_upload_attempt_is_fail_closed(stage):
+    c = circle()
+    c.set_target((0., .25), (0., .05))
+
+    def call(name):
+        def run(*_):
+            if name == stage:
+                raise RuntimeError(f"{name} failed")
+        return run
+    with pytest.raises(RuntimeError, match=f"{stage} failed"):
+        c.advance(call("upload"), call("step"), call("synchronize"))
+    assert c.failed and c.applied == (0., 0.)
+    for later in (c.record, lambda: c.set_target((0., .1), (0., 0.)),
+                  lambda: c.advance(*(lambda *_: None,) * 3)):
+        with pytest.raises(RuntimeError, match="unknown"):
+            later()
+
+
+def test_band_covers_offset_box_with_margin():
+    assert circle().band() == (50 - 2 - 10 - 2, 50 + 2 + 10 + 2, 40 - 3 - 10 - 2, 40 + 3 + 10 + 2)
+
+
+@pytest.mark.parametrize("radius", (R_MIN, 4.3, 7.5, 12.0))
+def test_refill_donor_sets_nonempty_within_bounds(radius):
+    # Evidence for R >= R_MIN and |dc|_2 <= MAX_STEP on a fixed sample, not a proof.
+    assert enumerate_donor_coverage(radius) == 0
+
+
+pytest.importorskip("pycuda.driver")
+from CelerisLab import simulation as module  # noqa: E402
+from tests.unit.test_simulation_lifecycle import Owner, make_simulation  # noqa: E402
+
+
+class ClosableOwner(Owner):
+    """Owner double that refuses activation after close, like CudaContext.activate."""
+    def activate(self):
+        if self.is_closed:
+            raise RuntimeError("CUDA context is closed")
+        return super().activate()
+
+
+def moving_simulation(monkeypatch, *, moving=True):
+    sim = make_simulation(monkeypatch)
+    obj = SimpleNamespace(motion=circle() if moving else None, obj_id=0,
+                          center=(50., 40.), radius=10.)
+    writes = []
+    sim._failed, sim._initialized = None, True
+    sim.bodies = SimpleNamespace(objects=[obj], count=1, jets=[], action_gpu=None, obs_gpu=None,
+                                 has_pending_edit=lambda: False,
+                                 _upload_action_async=lambda s: writes.append("action"),
+                                 zero_obs_async=lambda s: writes.append("obs"),
+                                 increment_obs_steps=lambda n: writes.append("count"),
+                                 download_obs_full_async=lambda s: writes.append("download"))
+    sim.stepper = SimpleNamespace(step=lambda n, **k: writes.append(("step", n)), blocked=None,
+                                  moving=None)
+    sim._assert_runtime_contracts = lambda: None
+    return sim, writes
+
+
+STREAM = SimpleNamespace(synchronize=lambda: None)
+
+
+def test_positive_steps_rejected_before_device_writes(monkeypatch):
+    sim, writes = moving_simulation(monkeypatch)
+    sim.bodies.objects[0].motion.set_target((0., .1), (0., .01))
+    sim.bodies.objects[0].motion.applied = (0., 0.)
+    for call in (lambda: sim.run(1, stream=STREAM), lambda: sim.run(5, stream=STREAM),
+                 lambda: sim.step(), lambda: sim.run(3, stream=STREAM, checkpoint_interval=1)):
+        with pytest.raises(RuntimeError, match="moving kernels"):
+            call()
+    with pytest.raises(RuntimeError, match="moving kernels"):  # nominal at rest as well
+        sim.bodies.objects[0].motion.pending = None
+        sim.run(1)
+    sim.bodies.objects[0].state = SimpleNamespace(omega=.02)  # spin does not lift the guard
+    with pytest.raises(RuntimeError, match="moving kernels"):
+        sim.run(2, stream=STREAM)
+    assert writes == [] and sim.bodies.objects[0].motion.applied == (0., 0.)
+
+
+def test_run_zero_keeps_pending_and_static_runs_proceed(monkeypatch):
+    sim, writes = moving_simulation(monkeypatch)
+    motion = sim.bodies.objects[0].motion
+    motion.set_target((0., .1), (0., .01))
+    sim.run(0, stream=STREAM)  # existing action/obs work only; no motion upload or commit
+    assert motion.pending == ((0., .1), (0., .01)) and motion.applied == (0., 0.)
+    assert ("step", 0) in writes
+    static, writes = moving_simulation(monkeypatch, moving=False)
+    static.save_checkpoint = lambda *a: writes.append("checkpoint")
+    static.run(3, stream=STREAM)
+    static.run(4, stream=STREAM, checkpoint_interval=2)
+    assert writes.count(("step", 3)) == 1 and writes.count(("step", 2)) == 2
+
+
+def test_direct_stepper_entry_is_guarded():
+    from CelerisLab.lbm.stepper import LBMStepper
+    launches = []
+    module = SimpleNamespace(get_function=lambda name: lambda *a, **k: launches.append(name))
+    cfg = SimpleNamespace(streaming="double_buffer", threads_per_block=128, is_d2q9=True,
+                          nx=64, ny=32, nz=1)
+    field = SimpleNamespace(ddf_gpu="ddf", temp_gpu="temp", n_curved=0)
+    stepper = LBMStepper(field, module, cfg)
+    assert stepper.blocked is None
+    stepper.blocked = "translating circles need the moving kernels"
+    for n in (1, 3):
+        with pytest.raises(RuntimeError, match="moving kernels"):
+            stepper.step(n, action_gpu=None, obs_gpu=None)
+    assert launches == [] and (field.ddf_gpu, field.temp_gpu) == ("ddf", "temp")
+    assert stepper.step_count == 0
+    stepper.step(0, action_gpu=None, obs_gpu=None)  # zero steps stay a no-op
+    assert launches == [] and stepper.step_count == 0
+
+
+def test_moving_body_rejects_topology_and_checkpoint(monkeypatch):
+    sim, _ = moving_simulation(monkeypatch)
+    for call in (sim.sync_bodies, sim.save_checkpoint, lambda: sim.load_checkpoint("x.h5")):
+        with pytest.raises(ValueError, match="translating circles"):
+            call()
+
+
+# -- Real-object routing, initialize rejections, packing, lifecycle -------------
+import numpy as np  # noqa: E402
+
+from CelerisLab.body.geometry.circle import CircleGeometry  # noqa: E402
+from CelerisLab.body.manager import ObjectManager  # noqa: E402
+from CelerisLab.body.objects import SimObject  # noqa: E402
+from CelerisLab.lbm.descriptors import FLUID, SOLID  # noqa: E402
+
+
+class Cfg:
+    dim, is_d3q19, is_d2q9 = 2, False, True
+    streaming, ddf_shifting, store_precision, y_wall_bc = "double_buffer", False, "FP32", "bounce_back"
+
+
+def body(x, y, r, translation=None, **kw):
+    obj = SimObject(obj_id=-1, geometry=CircleGeometry(x, y, r), center=(x, y), radius=r, **kw)
+    if translation is not None:
+        obj.motion = MovingCircle((x, y), r, translation)
+    return obj
+
+
+def routed_simulation(monkeypatch):
+    sim = make_simulation(monkeypatch)
+    sim._failed, sim._initialized, sim.lbm_cfg = None, False, Cfg()
+    sim.stepper = SimpleNamespace(blocked=None)
+    sim.bodies = ObjectManager(nx=120, ny=80, nz=1, nq=9, cfg=Cfg())
+    return sim
+
+
+def test_add_body_and_set_body_route_targets(monkeypatch):
+    sim = routed_simulation(monkeypatch)
+    bid = sim.add_body("circle", (50., 40.), 10., translation=Translation((-2., 2., -3., 3.), .1))
+    still = sim.add_body("circle", (100., 40.), 5.)
+    motion = sim.bodies.get(bid).motion
+    sim._moving_tables = "tables"
+    sim._attach_moving()
+    assert sim.stepper.blocked is None
+    sim.set_body(bid, offset=(0., .3), vy=.05)
+    assert sim.stepper.blocked == module.PENDING_BLOCKED  # direct steps never consume a target
+    sim.set_body(bid, offset=(0., .4), vy=.05)  # last wins
+    for kwargs in (dict(offset=(0., .6)), dict(vy=.05), dict(vx=.05),
+                   dict(offset=(0., .1), vy=.2), dict(offset=(0., 0.), vy=.01)):
+        with pytest.raises(ValueError):
+            sim.set_body(bid, **kwargs)
+        assert motion.pending == ((0., .4), (0., .05))  # invalid input keeps the target
+    sim.set_body(bid, omega=.01)  # omega-only keeps the pending target
+    obj = sim.bodies.get(bid)
+    assert obj.state.omega == .01 and (obj.state.vx, obj.state.vy) == (0., 0.)
+    assert motion.pending == ((0., .4), (0., .05)) and motion.applied == (0., 0.)
+    with pytest.raises(ValueError, match="not declared"):
+        sim.set_body(still, offset=(0., .1))
+    sim.set_body(still, vx=.01)  # stationary API unchanged
+    assert sim.bodies.get(still).state.vx == .01 and sim.bodies.get(still).motion is None
+    for kwargs in (dict(type="sensor"), dict(surface_actuation=dict(angle=0., half_width=.2,
+                                                                    normal_velocity=.01))):
+        kind = kwargs.pop("type", "circle")
+        with pytest.raises(ValueError, match="translation"):
+            sim.add_body(kind, (20., 20.), 5., translation=Translation((0., 0., 0., 0.), .1),
+                         **kwargs)
+
+
+class Tables:
+    """MovingTables double: records uploads, steps and diagnostics reads."""
+    DIAG = ("no_donor",)
+
+    def __init__(self, log, no_donor=0., fail=None):
+        self.log, self.no_donor, self.fail = log, no_donor, fail
+
+    def upload_records(self, circles, stream):
+        records = [m.record() for m, _ in circles]
+        self.log.append(("upload", records))
+        uploads = sum(1 for w in self.log if isinstance(w, tuple) and w[0] == "upload")
+        if self.fail == "upload" or (self.fail == "hold_upload" and uploads == 2):
+            raise RuntimeError("upload failed")
+        return any(r[0] != r[1] for r in records)
+
+    def zero_diagnostics(self, stream):
+        self.log.append("zero_diag")
+
+    def download_diagnostics(self):
+        return np.array([[self.no_donor]], np.float32)
+
+
+def advancing_simulation(monkeypatch, **tables):
+    sim, writes = moving_simulation(monkeypatch)
+    sim.bodies.objects[0].motion = circle()
+    sim.stepper = SimpleNamespace(
+        moving=Tables(writes, **tables), blocked=None,
+        _step=lambda n, advance=False, **k: writes.append(("step", n, advance)))
+    from CelerisLab.lbm.stepper import LBMStepper
+    sim.stepper.step = LBMStepper.step.__get__(sim.stepper)  # real public entry
+    sim._moving_tables = sim.stepper.moving
+    return sim, writes, sim.bodies.objects[0].motion
+
+
+def test_run_one_advances_and_commits_after_sync(monkeypatch):
+    sim, writes, motion = advancing_simulation(monkeypatch)
+    motion.set_target((0., .25), (0., .05))
+    with pytest.raises(ValueError, match="run\\(1\\)"):  # multi-step with a target
+        sim.run(2, stream=STREAM)
+    assert writes == [] and motion.pending is not None
+    sim.run(0, stream=STREAM)  # zero steps never consume the target
+    assert motion.pending is not None
+    sim.run(1, stream=STREAM)
+    uploads = [w for w in writes if isinstance(w, tuple) and w[0] == "upload"]
+    assert uploads[0][1] == [((50., 40.), (50., 40.25), (0., .05))]
+    assert uploads[1][1] == [((50., 40.25), (50., 40.25), (0., 0.))]  # first hold after move
+    assert ("step", 1, True) in writes
+    assert motion.applied == (0., .25) and motion.pending is None
+    writes.clear()
+    sim.run(3, stream=STREAM)  # hold steps: no advance
+    assert ("step", 3, False) in writes and motion.applied == (0., .25)
+
+
+@pytest.mark.parametrize("failure", ("upload", "first_sync", "no_donor", "hold_upload",
+                                     "hold_sync"))
+def test_failed_advance_fails_simulation(monkeypatch, failure):
+    kwargs = dict(no_donor=1.) if failure == "no_donor" else dict(fail=failure)
+    sim, writes, motion = advancing_simulation(monkeypatch, **kwargs)
+    motion._on_failure = sim._fail_moving
+    motion.set_target((0., .25), (0., .05))
+    syncs = []
+
+    def synchronize():  # stands in for an asynchronous failure surfacing at sync
+        syncs.append(1)
+        if (failure, len(syncs)) in (("first_sync", 1), ("hold_sync", 2)):
+            raise RuntimeError("sync failed")
+    with pytest.raises(RuntimeError):
+        sim.run(1, stream=SimpleNamespace(synchronize=synchronize))
+    # Before the step-group sync and no_donor check nothing is committed; the hold
+    # record is staged after the commit (applied already matches the moved device).
+    committed = failure in ("hold_upload", "hold_sync")
+    assert motion.applied == ((0., .25) if committed else (0., 0.))
+    assert (motion.pending is None) == committed
+    assert motion.failed and sim._failed is not None
+    assert sim.stepper.blocked == sim._failed  # direct stepper reuse refused too
+    with pytest.raises(RuntimeError, match="Simulation failed"):
+        sim.run(1, stream=STREAM)
+
+
+def test_failed_hold_step_blocks_the_direct_stepper(monkeypatch):
+    # No target is pending, so only the failure latch can refuse direct steps.
+    from CelerisLab.lbm.stepper import LBMStepper
+    sim, writes, motion = advancing_simulation(monkeypatch, fail="upload")
+    motion._on_failure = sim._fail_moving
+    launches = []
+    kernels = SimpleNamespace(get_function=lambda name: lambda *a, **k: launches.append(name))
+    cfg = SimpleNamespace(streaming="double_buffer", threads_per_block=128, is_d2q9=True,
+                          nx=64, ny=32, nz=1)
+    field = SimpleNamespace(ddf_gpu="ddf", temp_gpu="temp", flag_gpu="flag", n_curved=0)
+    stepper = LBMStepper(field, kernels, cfg)
+    stepper.moving = sim.stepper.moving
+    sim.stepper = stepper
+    assert motion.pending is None and stepper.blocked is None
+    with pytest.raises(RuntimeError, match="upload failed"):
+        sim.run(2, stream=STREAM)  # hold step
+    assert sim._failed is not None and stepper.blocked == sim._failed
+    with pytest.raises(RuntimeError, match="moving circle step failed"):
+        stepper.step(1, action_gpu=None, obs_gpu=None)
+    assert launches == [] and stepper.step_count == 0
+    assert (field.ddf_gpu, field.temp_gpu) == ("ddf", "temp")
+
+
+def test_moving_tables_layout_and_attachment(monkeypatch):
+    from CelerisLab.body.motion import MovingTables
+    a = MovingCircle((50., 40.), 10., Translation((-1., 1., -3., 3.), .1))
+    b = MovingCircle((90.25, 30.5), 5., Translation((0., 0., 0., 0.), .1))
+    band, geom = MovingTables.layout([(a, 2), (b, 5)])
+    wa, ha = 63 - 37 + 1, 55 - 25 + 1
+    assert band == [37, 25, wa, ha, 0, 2, 83, 23, 16, 16, wa * ha * 8, 5]
+    assert geom == [50., 40., 50., 40., 10., 0., 0., 90.25, 30.5, 90.25, 30.5, 5., 0., 0.]  # hold
+    sim = routed_simulation(monkeypatch)
+    sim.add_body("circle", (50., 40.), 10., translation=Translation((0., 0., 0., 0.), .1))
+    sim._moving_tables = "tables"
+    sim._attach_moving()
+    assert sim.stepper.moving == "tables" and sim.stepper.blocked is None
+    sim._moving_tables = None
+    sim._attach_moving()
+    assert sim.stepper.moving is None and sim.stepper.blocked is not None
+
+
+class CountingOwner(Owner):
+    """Owner double that counts nested activations like CudaContext."""
+    def __init__(self):
+        super().__init__()
+        self.active = 0
+
+    def activate(self):
+        from contextlib import contextmanager
+
+        @contextmanager
+        def scope():
+            self.active += 1
+            try:
+                yield self
+            finally:
+                self.active -= 1
+        return scope()
+
+
+def test_public_diagnostics_keep_owner_active_through_sync_and_download(monkeypatch):
+    sim = make_simulation(monkeypatch)
+    sim.ctx, sim._failed = CountingOwner(), None
+    seen = []
+    sim._stream = SimpleNamespace(synchronize=lambda: seen.append(("sync", sim.ctx.active)))
+    expected = np.arange(6, dtype=np.float32).reshape(1, 6)
+
+    def download():
+        seen.append(("download", sim.ctx.active))
+        return expected
+    sim.stepper = SimpleNamespace(moving=SimpleNamespace(download_diagnostics=download))
+    assert sim.moving_diagnostics() is expected
+    assert [name for name, _ in seen] == ["sync", "download"]
+    assert all(active > 0 for _, active in seen) and sim.ctx.active == 0
+
+
+def test_stepper_launches_moving_links_before_static_curved():
+    from CelerisLab.lbm.stepper import LBMStepper
+    launches = []
+    module = SimpleNamespace(get_function=lambda name: lambda *a, **k: launches.append((name, a, k)))
+    cfg = SimpleNamespace(streaming="double_buffer", threads_per_block=128, is_d2q9=True,
+                          nx=64, ny=32, nz=1)
+    curved = SimpleNamespace(**{f"{name}_gpu": name for name in (
+        "fluid_idx", "dir", "q", "rx", "ry", "rz", "fallback_class", "body_id")},
+        surface_velocity=None, surface_velocity_gpu=None, link_exchange_gpu=None)
+    field = SimpleNamespace(ddf_gpu="ddf", temp_gpu="temp", flag_gpu="flag", n_curved=5,
+                            curved=curved, invalidate_host_ddf_cache=lambda: None)
+    stepper = LBMStepper(field, module, cfg)
+    stepper._launch_force_region = stepper._launch_sensor = lambda *a, **k: None
+    stepper.moving = SimpleNamespace(band_gpu="band", geom_gpu="geom", diag_gpu="diag",
+                                     slot_debug_gpu=None, n_bodies=1, n_slots=300)
+    stepper.step(2, action_gpu="act", obs_gpu="obs")
+    names = [name for name, _, _ in launches]
+    assert names == ["MovingCircleLinks", "CurvedBoundaryKernel", "OneStep"] * 2
+    assert stepper.step_count == 2
+    assert launches[1][1][0] == "ddf" and launches[3][1][0] == "temp"  # swap between steps
+    args, kw = launches[0][1], launches[0][2]
+    assert args[:4] == ("ddf", "flag", "band", "geom") and int(args[5]) == 300
+    assert kw["grid"] == (3, 1, 1) and kw["block"] == (128, 1, 1)
+
+
+def test_recompile_installs_a_guarded_stepper(monkeypatch):
+    from CelerisLab.lbm.stepper import LBMStepper
+    sim = routed_simulation(monkeypatch)
+    sim.add_body("circle", (50., 40.), 10., translation=Translation((0., 0., -1., 1.), .1))
+    launches = []
+    kernel_module = SimpleNamespace(get_function=lambda name: lambda *a, **k: launches.append(name))
+    monkeypatch.setattr(module.compiler, "build_module", lambda *a, **k: kernel_module)
+    sim.lbm_cfg = SimpleNamespace(validate_runtime=lambda: None, viscosity=.01, omega=1.,
+                                  streaming="double_buffer", threads_per_block=128,
+                                  is_d2q9=True, nx=120, ny=80, nz=1)
+    sim._compiled_viscosity = .01
+    sim._assert_allocation_contract = lambda: None
+    sim._resolve_compile_arch = lambda: "sm_70"
+    sim._assert_object_count_contract = lambda **k: None
+    sim.field = SimpleNamespace(module=None, invalidate_module_constant_cache=lambda: None,
+                                _upload_params=lambda: None, _runtime_params={},
+                                ddf_gpu="ddf", temp_gpu="temp", n_curved=0)
+    sim.stepper = SimpleNamespace(_step_count=7, blocked=None)
+    sim._recompile()
+    assert isinstance(sim.stepper, LBMStepper) and sim.stepper.step_count == 7
+    with pytest.raises(RuntimeError, match="moving kernels"):
+        sim.stepper.step(1, action_gpu=None, obs_gpu=None)
+    assert launches == [] and (sim.field.ddf_gpu, sim.field.temp_gpu) == ("ddf", "temp")
+    assert sim.stepper.step_count == 7
+
+
+def validate(monkeypatch, objects, flags=None, **cfg):
+    sim = routed_simulation(monkeypatch)
+    for obj in objects:
+        sim.bodies.add(obj)
+    for name, value in cfg.items():
+        setattr(sim.lbm_cfg, name, value)
+    base = np.full(120 * 80, FLUID, np.uint16) if flags is None else flags
+    sim.field = SimpleNamespace(nx=120, ny=80, build_channel_flags=lambda: base.copy())
+    sim._validate_moving_bodies()
+
+
+def test_initialize_rejects_band_and_modes(monkeypatch):
+    t = Translation((-1., 1., -3., 3.), .1)
+    validate(monkeypatch, [body(50., 40., 10., t), body(100., 40., 5.)])  # accepted
+    walls = np.full(120 * 80, FLUID, np.uint16)
+    walls[60 + 30 * 120] = SOLID
+    cases = [([body(10., 40., 10., t)], None, {}, "domain edge"),
+             ([body(50., 40., 10., t)], walls, {}, "non-fluid"),
+             ([body(50., 40., 10., t), body(68., 40., 4.)], None, {}, "overlaps"),
+             ([body(50., 40., 10., t), body(50., 60., 3., is_sensor=True)], None, {}, "overlaps"),
+             ([body(50., 40., 10., t), body(64., 40., 3., is_force_region=True)], None, {}, "overlaps"),
+             ([body(50., 40., 10., t)], None, {"streaming": "esopull"}, "double_buffer"),
+             ([body(50., 40., 10., t)], None, {"ddf_shifting": True}, "unshifted"),
+             ([body(50., 40., 10., t)], None, {"is_d2q9": False}, "D2Q9")]
+    for objects, flags, cfg, match in cases:
+        with pytest.raises(ValueError, match=match):
+            validate(monkeypatch, objects, flags, **cfg)
+
+
+def test_initialize_validates_before_mutation(monkeypatch):
+    sim = routed_simulation(monkeypatch)
+    sim.bodies.add(body(10., 40., 10., Translation((0., 0., 0., 0.), .1)))
+    uploads = []
+    sim.field = SimpleNamespace(nx=120, ny=80, build_channel_flags=lambda: np.full(9600, FLUID, np.uint16),
+                                upload_flags=lambda: uploads.append("flags"))
+    sim._assert_allocation_contract = lambda: None
+    with pytest.raises(ValueError, match="domain edge"):
+        sim.initialize()
+    assert uploads == [] and not sim._initialized
+
+
+def test_moving_circle_excluded_from_static_links():
+    def manager(*objects):
+        mgr = ObjectManager(nx=120, ny=80, nz=1, nq=9, cfg=Cfg())
+        ids = [mgr.add(obj) for obj in objects]
+        return mgr, ids
+    reference_mgr, _ = manager(body(100., 40., 5.), body(20., 40., 6.))
+    mixed_mgr, ids = manager(body(100., 40., 5.), body(20., 40., 6.),
+                             body(55., 40., 10., Translation((0., 0., -1., 1.), .1)))
+    reference = reference_mgr.build_compact_lists_for(reference_mgr.objects)
+    mixed = mixed_mgr.build_compact_lists_for(mixed_mgr.objects)
+    for a, b in zip(reference[:8], mixed[:8]):  # static link columns byte-identical
+        assert np.asarray(a).tobytes() == np.asarray(b).tobytes()
+    assert len(reference[0]) > 0 and ids[2] not in set(np.asarray(mixed[3]).tolist())
+
+
+def test_failed_moving_step_marks_simulation_failed(monkeypatch):
+    sim = routed_simulation(monkeypatch)
+    sim.ctx = ClosableOwner()
+    bid = sim.add_body("circle", (50., 40.), 10., translation=Translation((0., 0., -1., 1.), .1))
+    motion = sim.bodies.get(bid).motion
+    motion.set_target((0., .2), (0., .05))
+
+    def fail(*_):
+        raise RuntimeError("sync failed")
+    with pytest.raises(RuntimeError, match="sync failed"):
+        motion.advance(lambda *_: None, lambda: None, fail)
+    assert sim._failed is not None and motion.applied == (0., 0.)
+    for call in (lambda: sim.run(0), lambda: sim.run(1), sim.step, sim.sync_bodies,
+                 sim.save_checkpoint, lambda: sim.load_checkpoint("x.h5")):
+        with pytest.raises(RuntimeError, match="Simulation failed"):
+            call()
+    with pytest.raises(RuntimeError, match="unknown"):  # no retry on the failed circle
+        motion.advance(lambda *_: None, lambda: None, lambda: None)
+    with sim.ctx.activate():  # a refused close keeps the failed state and resources
+        with pytest.raises(RuntimeError, match="active CUDA"):
+            sim.close()
+    assert sim._failed is not None and sim.bodies is not None and not sim._closed
+    sim.close()  # retry outside the activation succeeds
+    assert sim._closed and sim._failed is not None
+    with pytest.raises(RuntimeError, match="CUDA context is closed"):  # no re-entry after close
+        sim.run(1)
+    with pytest.raises(RuntimeError, match="unknown"):  # the failed circle stays refused
+        motion.record()
@@ -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):
@@ -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()

@@ -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()
@@ -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()
@@ -1,5 +1,6 @@
 """CPU mocks execute production remap lifecycle; no GPU numerical claim."""
 from dataclasses import replace
+from types import SimpleNamespace
 import unittest
 from unittest.mock import patch

@@ -69,10 +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)
-        hosts, pairs = run._boundaries, run.hierarchy.pairs
+        records, pairs = run._boundaries, run.hierarchy.pairs
        for attribute, value in (('_wall_aware', True),
-                                 ('_boundaries', [*hosts[:2], object()]),        # child boundary
-                                 ('_boundaries', [hosts[0], object(), hosts[2]]),  # fine boundary
+                                 # level-2 record whose host geometry differs from objects=()
+                                 ('_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)
@@ -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()


@@ -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)

@@ -294,8 +294,11 @@ def test_passive_overlap_failure_keeps_live_state(device, entry, second_x):
    sim = SimpleNamespace(bodies=bodies, field=field, lbm_cfg=field.cfg,
                          _initialized=entry == "runtime_edit",
                          _assert_allocation_contract=lambda: None,
-                          stream=SimpleNamespace(synchronize=device.synchronize))
-    sim._require_usable = Simulation._require_usable.__get__(sim)
+                          stream=SimpleNamespace(synchronize=device.synchronize),
+                          stepper=SimpleNamespace(blocked=None))
+    for name in ("_require_usable", "_moving_bodies", "_reject_moving_bodies",
+                 "_validate_moving_bodies", "_attach_moving"):
+        setattr(sim, name, getattr(Simulation, name).__get__(sim))
    with pytest.raises(ValueError, match="active fluid receivers|unique cut-link writers"):
        if entry == "sync":
            bodies.sync_to_gpu(field)
@@ -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')
Author	SHA1	Message	Date
ClaudeandClaude Opus 5.5	8c1f3c4d38	docs: Dütsch protocol — adopt Grok's cycle-index labels (K = first scored cycle) Same cycles and pairs as `8505a723`, relabelled: proof pairs 3→4..5→6 give K_case = 7, window K..K+2 with remain-settled pairs, last scored cycle 11. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:50:22 +00:00
ClaudeandClaude Opus 5.5	8505a723a5	docs: Dütsch protocol — worked cycle-index example, sparse sampling, driver unit tests Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:48:30 +00:00
ClaudeandClaude Opus 5.5	10cea29e8e	docs: Dütsch protocol — fold the closed domain order, settling rule, window and criteria Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:44:38 +00:00
ClaudeandClaude Opus 5.5	cd73fcd689	docs: Dütsch protocol — LDA profile comparison qualitative until data exist Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:43:48 +00:00
ClaudeandClaude Opus 5.5	ee6e68bbfd	docs: Dütsch protocol — consolidate verified source facts and close protocol choices Source facts with page provenance (motion, quarter-cycle start, stations, Morison c_d/added-mass c_i definition, Table 1 as computed reference, no recovered uncertainty); lab-frame kinematics with base at the right turning point, interval-average wall velocity and force attribution, explicit cycle indexing; T multiple of 360 steps; native-channel domain sensitivity order; fluid-only interpolation; periodic gate with fixed-scale waveform/profile settling, N=3, cap 30, common window; engineering criteria labelled as such. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:43:10 +00:00
ClaudeandClaude Opus 5.5	f68f5264b2	docs: Dütsch protocol — mandatory periodic-regime gate ledger rows Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:30:38 +00:00
ClaudeandClaude Opus 5.5	03b23cb511	docs: Dütsch protocol — duration/periodic gate, phase convention, reference types Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:30:05 +00:00
ClaudeandClaude Opus 5.5	8ce810cc99	docs: Dütsch Re100/KC5 prescribed-oscillation protocol draft Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 10:28:24 +00:00
ClaudeandClaude Opus 5.5	653fd997a6	body: stage motion-record uploads in an owned page-locked buffer memcpy_htod_async needs page-locked host memory. MovingTables allocates one pinned FP64 buffer at construction (manager.py action pattern) and rewrites it only when the records change; every upload is followed by a stream sync in Simulation._run_moving before the next one, so no copy is in flight when it is reused. The advance test fixture now binds the real LBMStepper.step. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 09:16:28 +00:00
ClaudeandClaude Opus 5.5	25cc5cfd68	tests: a failed hold step refuses a real direct stepper call Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 09:05:34 +00:00
ClaudeandClaude Opus 5.5	6f9e40caa3	simulation: a failed moving step also blocks direct stepper calls Tests cover failures at the first upload, the step-group sync, the no_donor check, the hold-record upload and its asynchronous sync: before the commit point the target stays pending, after it applied matches the moved device; in every case run and direct stepper reuse are refused. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 09:01:28 +00:00
ClaudeandClaude Opus 5.5	c21a365eca	docs/tests: state the S3a commit point; cover a failed hold-record upload Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 08:59:45 +00:00
ClaudeandClaude Opus 5.5	d50b8ba1be	lbm: translating circles move with refill (S3a) MovingCircleCells reclassifies swept-band cells from immutable old/new FP64 predicates: covered cells add their raw momentum and mass to diagnostics and become solid; uncovered cells are refilled in all nine slots from donors fluid in both geometries (mean density, NEQ from the best-aligned donor, wall velocity at the cell), counting eq_only and no_donor. MovingCircleLinks adds the translation velocity to the wall velocity. Simulation.run(1) uploads the motion record, runs M1 then M2, synchronizes, rejects a nonzero no_donor count and only then commits the target and uploads the hold record; any failure marks the Simulation failed. Targets block direct stepper calls and multi-step runs. get_flags reads device flags back when translating circles exist. Contract section 10b. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 08:51:20 +00:00
ClaudeandClaude Opus 5.5	b920a303a4	simulation: keep the owner context active through moving diagnostics readback moving_diagnostics is context_bound, so stream synchronization and the DTOH copy run inside one owner activation. Tests: owner stays active through sync and download; launch order moving -> static curved -> OneStep with a real static link count and the buffer swap between steps. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 07:58:40 +00:00
ClaudeandClaude Opus 5.5	e25828968b	lbm: choose the moving-link fallback from FP64 q, as the host builder does Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 07:35:45 +00:00
ClaudeandClaude Opus 5.5	0f8a428f07	lbm: fixed-centre moving-circle cut links (S2a) MovingCircleLinks computes the cut links of each translating circle from an FP64 {cx, cy, r} table over its swept band (same q/hit/clamp/fallback order as build_cut_links, explicit _rn operations) and applies the unchanged apply_bouzidi_link with spin wall velocity. obs keeps its contract; a separate per-body diagnostics buffer records CME, the Wen term u_w*dm, link_dm and the link count. Moving tables are built at initialize and attached to the stepper, which then launches the kernel before the static curved kernel; until then direct and run steps stay blocked. Displacement is disabled: set_body offset/vx/vy on a translating circle raises before mutation. Contract section 10a documents S2a and the S3 pending rule. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 07:35:02 +00:00
ClaudeandClaude Opus 5.5	983020b625	tests: align closed-owner semantics; route moving-body checks in state regressions test_failed_moving_step_marks_simulation_failed uses an owner double that refuses activation after close like CudaContext and expects 'CUDA context is closed'; test_passive_overlap_failure_keeps_live_state binds the real moving-body validation methods and a stepper double. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 07:11:06 +00:00
ClaudeandClaude Opus 5.5	e1ed27b438	simulation: revert latch-only refresh (`e682d16`, `36eac713`) No reachable path leaves a stale latch (Grok review): initialize and _recompile recompute it, and post-init removal of a translating circle cannot commit. simulation.py is back to b448bc3f; the latch-only test is dropped. The recompile, close-retry and tightened after-close tests stay. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:57:02 +00:00
ClaudeandClaude Opus 5.5	71c12aedce	tests: recompile installs a guarded stepper; refused close keeps failed state Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:55:26 +00:00
ClaudeandClaude Opus 5.5	36eac713f4	simulation: refresh the stepper latch only when a moving body is added or removed Keeps the stationary add/remove path unchanged (test_surface_actuation builds a Simulation without a stepper). Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:53:35 +00:00
ClaudeandClaude Opus 5.5	e682d16c3a	simulation: derive the stepper latch from moving bodies on add/remove One helper sets stepper.blocked from _moving_bodies() after add_body, pre-init remove_body, initialize and stepper recreation, so the latch never disagrees with motion ownership. Tests: latch follows add/remove; the after-close assertion checks the actual refusal. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:53:12 +00:00
ClaudeandClaude Opus 5.5	b448bc3f71	lbm: guard direct stepper entry while translating circles exist LBMStepper.step raises before any launch or buffer swap when the owner has set 'blocked'; Simulation sets it when a translating circle is declared, at initialize and after stepper recreation. Tests: direct-stepper guard (zero launches, no swap, unchanged step count), spin/multi-step/checkpoint run variants, run(0) keeps the pending target, static runs proceed, mixed moving/static packing parity with the moving id absent, failure refusal of checkpoint calls, no retry and no re-entry after close. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:45:13 +00:00
ClaudeandClaude Opus 5.5	8a8e8a03ff	body: moving-step failure fails the Simulation; S1 acceptance tests A failed MovingCircle.advance now calls the owner hook, so the Simulation is marked failed (Jet lifecycle) and refuses run/step/sync until close(). Adds host tests for set_body routing on real objects, initialize band and mode rejections before mutation, exclusion from static link packing, and the Simulation failure lifecycle. Contract status line updated for S1. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:43:16 +00:00
ClaudeandClaude Opus 5.5	305ab2409b	body: host-only translating-circle state (S1, no kernels) MovingCircle keeps the applied offset and an optional pending target apart, validates targets against the applied offset (box, Euclidean 0.5-cell step, speed, no sliding hold) and commits only after upload, step and sync succeed; any failure after the upload attempt is fail-closed. Simulation.add_body (translation=...) declares one; moving circles are excluded from the packed static links, unsupported modes and swept-band placement are rejected at initialize, and every positive-step run raises before any device write until the moving kernels exist. Host tests include a refill-donor enumeration. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:34:06 +00:00
ClaudeandClaude Opus 5.5	c986343a15	docs: moving-bodies contract (from S0 `26f29b09`) with applied/pending fixes Carries the reviewed S0 draft into the S1 branch and fixes: applied stores only the offset, explicit 'pending is not None', reject a hold target with nonzero translational velocity, single Euclidean 0.5-cell bound. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-30 06:31:26 +00:00
Frank14f	b2034f87e0	docs: distinguish fresh and carried PTX evidence	2026-10-28 20:20:01 +08:00
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
ClaudeandClaude Opus 5.5	65339ab54e	test_runtime_child_remap: inject remap body rejects through the boundary record The removed _boundaries list made the two body-reject patches no-ops. Inject them via _device_boundaries (level-2 record with foreign host geometry; non-None level-1 record). set_jet_amplitude docstring names the public jets view. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-28 10:22:29 +00:00
ClaudeandClaude Opus 5.5	ae766f7112	FixedPatchStepper: one boundary/Jet record per level, host preparation helpers - Jets and host boundary geometry are views of the per-level _BlockedBoundary record (drop duplicate _level_jets/_boundaries registration). - Host preparation split into _level_layouts, _host_boundaries, _fluid_leaves (numpy only, before any device work); constructor order unchanged. - set_jet_amplitude: explicit prepare -> device write -> commit; behaviour unchanged (invalid input raises before any change, uncertain write fails closed). - Pair accessor: FixedPatch2D.rectangles, FixedPatchSet2D.rectangles and fine_level_offset=(0, 0); removes two incidental isinstance branches. Geometry/representation unchanged; scalar/set types stay distinct. - @property read-only views; blocks.py responsibility docstring corrected. - Migrate test_runtime_child_remap to the public boundaries view. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>	2026-10-28 10:09:25 +00:00