CelerisLab/docs/moving_bodies.md

Moving circles on the uniform dense grid — S0 DRAFT contract

Status: DRAFT for review, not frozen. No implementation exists. 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.

4. State, masks and authority

• Host ObjectState is the single authority for offset = (dx, dy), velocity = (vx, vy) and omega. Velocity is supplied, never differenced from offsets.
• 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.)
• 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.
• 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.
• No valid donor for the NEQ part → equilibrium only (counted). No valid neighbour at all → cannot happen when the per-step move is < 1 cell (validated on host); counted anyway.
• 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. To be checked analytically in review: q < 1/2 branch (reads k_f - c), halfway fallback, and that every fluid pull from a solid slot is BFL-written.

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.

Accumulators (separate, per body, per step window; all reported raw): link_cme, link_wen_correction (= sum u_w dm), node_cover, node_uncover, plus counters of crossings and refill fallbacks, 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 ≥ 1 cell, speed above the declared bound.
• Host state is committed only after a successful upload (same contract as Jet amplitude updates); an uncertain device-write failure makes the simulation unusable.
• Device-side soft diagnostics (refill fallbacks) are counters, not errors.

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 et al. JFS 2006 under verification; vacuum vs in-fluid natural frequency must not be mapped silently). Stability/convergence is gated by a halved-time-step run, not assumed from the mass ratio.

9. Gates (map; tolerances OPEN, declared before runs)

G0 zero-motion regression (stationary path unchanged; moving path at rest reproduces stationary link set and DDF) · 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.

10. 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.