Compare commits

..
Author SHA1 Message Date
Frank14f 4f17c3c586 motion: exercise real step ownership and move host test helpers out of src 2026-10-31 11:24:09 +08:00
3f7f303e17 research: VIV driver reads diagnostics by MovingTables.DIAG name; coupling-contract README
Index values are unchanged (1, 3, 7, 10, 13); the run recorded at fb507aa1
used the same columns. Research-only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-31 11:20:57 +08:00
ClaudeandClaude Opus 5.5 d5b13d7f31 motion: mark MovingCircle.advance as a test helper, not the production path
Docstring only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-31 01:24:35 +00:00
ClaudeandClaude Opus 5.5 353a7b03c4 stepper: bind moving-circle kernels only for D2Q9 double buffer
Both kernels are compiled under DIM == 2, but LBMStepper requested them for
every double-buffer module, so stationary D3Q19 double-buffer construction
would fail at get_function. Moving circles already reject D3Q19 and EsoPull
at initialize. Stale latch comment updated. Parametrized CPU test with a
module double that lacks the symbols.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-31 01:23:54 +00:00
ClaudeandClaude Opus 5.5 0d6df6e6b8 simulation: consume a zero-motion translation target; empty diagnostics shape
A pending target equal to the applied offset uploaded the same record, so
_run_moving saw no advance and never committed it: pending stayed latched,
blocking run(n > 1) and direct steps. Commit whenever a target was pending
(no_donor check and hold-record re-upload still only on advancing steps).
moving_diagnostics() without moving circles now returns (0, DIAG_FLOATS)
instead of the stale (0, 6). Tests for both.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-31 01:23:04 +00:00
ClaudeandClaude Opus 5.5 7886906517 docs: moving-circle docstrings — tables-attached guard, vx/vy meaning, center keyword
Docstring and documentation wording only; no code change.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-30 16:33:12 +00:00
ClaudeandClaude Opus 5.5 48ee0982c9 tests/docs: unique-name driver loading; moving-path lifecycle and settings notes; VIV predeclare renamed
Tests load research drivers via importlib under unique module names
instead of a bare 'import run' (shared-process collision between the
Dütsch and VIV drivers). lifecycle.md qualifies immutable geometry and the
host-only get_flags for translating circles; supported_settings.md points to
the dense-only moving-circle feature; leontini_viv_predeclare.md is renamed
viv_re150_predeclare.md (Re 150 Kolahdouz/Qin case) with references updated.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-30 16:12:05 +00:00
ClaudeandClaude Opus 5.5 307468b671 docs: moving circles — current status and user guide, BC capability, stage-1 evidence
docs/moving_bodies.md gains a current status/user-guide section; the staged
contract below it is kept unchanged as the design record.
configs/CONFIG.md states the single dense outer-boundary package.
docs/validation/moving_bodies_stage1.md summarizes the Dütsch and VIV
evidence with SHAs and binding limitations. One stale source comment fixed
(comment only).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-30 16:10:32 +00:00
15 changed files with 351 additions and 145 deletions
+2 -2
View File
@@ -6,8 +6,8 @@ it is; it adds no framework. Line references are to candidate 0f8fdbc3.
| Class | Contents | How it changes | Owner rule |
|---|---|---|---|
| 1. Module and its globals | Compile-time (visible to the compiler, fixed per module): `config.h` macros incl. `NX`/`NY`, `COLLISION_MODEL`, `BLOCK_ROOT_BC`, `WALL_*`. Runtime per-module globals (not visible to the compiler; written by host code): `__constant__` `d_params`, `d_w`, `d_cx`, `d_cy`, `d_block_channel`, `d_block_wall`; plus bound kernel handles | Globals are written once by the owner right after `load_module`. Uniform grid: `Simulation._recompile` (construction, body-topology sync) loads a new module and re-uploads; `update_runtime_params` rewrites `d_params` (omega, force, n_objects) between runs. Steppers: never rewritten after construction, so they are immutable for the instance's lifetime although they are runtime values. The blocked slot/reservoir root reads its slot from the `d_block_wall` constant (the former dense route used `WALL_*` macros; both were bit-identical in the parity gates). | Every owner loads its own module (`compiler.load_module` = `module_from_file`, no cache). A loaded module is never shared between live owners: globals are per module, so sharing would share mutable `d_params`. |
| 2. Geometry, layout, tables | flags, curved-link tables, block layouts and halo pairs, transfer (explode/coalesce) tables, block kinds | Immutable for the life of an instance. Uniform grid: rebuilt by body sync. `FixedPatchStepper`: a geometry change is a new instance (`remap_finest`), committed only after a successful transfer. | No in-place table edits. |
| 3. Controls and dependent boundary columns | Jet amplitude -> per-link `surface_velocity`; body actions (uniform grid) | Between steps, without a table or module rebuild (`set_jet_amplitude`). | Input (index, level, amplitude) is validated before any device write; invalid input changes nothing and the instance stays usable. The host copy commits only after the upload returned. Any error after validation leaves the device column undefined and fails the instance: `FixedPatchStepper` status `failed`, `Simulation._failed`; every later device-touching public method (`Simulation`: `stream`, `set_jet_amplitude`, `initialize`, `sync_bodies`, `run`/`step`, `read_surface_exchange`, `read_body(s)`, `get_macroscopic`, `get_ddf`, `update_runtime_params`, `snapshot`, `restore`, save/load checkpoint) raises; cached host readers and host-side staging (e.g. `get_flags`, `add_body`/`remove_body` before a sync) do not touch the device and stay callable. `close()` still releases resources and stays retryable. |
| 2. Geometry, layout, tables | flags, curved-link tables, block layouts and halo pairs, transfer (explode/coalesce) tables, block kinds | Immutable for the life of an instance, except translating circles on the uniform grid: their swept-band flags and cut links are rewritten on the device every advancing step (`MovingCircleCells`/`MovingCircleLinks`, docs/moving_bodies.md); the host flag mirror is then stale for those cells. Uniform grid otherwise: rebuilt by body sync. `FixedPatchStepper`: a geometry change is a new instance (`remap_finest`), committed only after a successful transfer. | No in-place table edits. |
| 3. Controls and dependent boundary columns | Jet amplitude -> per-link `surface_velocity`; body actions (uniform grid) | Between steps, without a table or module rebuild (`set_jet_amplitude`). | Input (index, level, amplitude) is validated before any device write; invalid input changes nothing and the instance stays usable. The host copy commits only after the upload returned. Any error after validation leaves the device column undefined and fails the instance: `FixedPatchStepper` status `failed`, `Simulation._failed`; every later device-touching public method (`Simulation`: `stream`, `set_jet_amplitude`, `initialize`, `sync_bodies`, `run`/`step`, `read_surface_exchange`, `read_body(s)`, `get_macroscopic`, `get_ddf`, `update_runtime_params`, `snapshot`, `restore`, save/load checkpoint) raises; cached host readers and host-side staging (e.g. `get_flags` without translating circles, `add_body`/`remove_body` before a sync) do not touch the device and stay callable; with translating circles `get_flags` reads the device flags back and refuses after a failure. Translating-circle targets follow the same rule: `set_body` validates before mutation, the applied offset commits only after a successful synchronized step, and a failed moving step fails the instance (docs/moving_bodies.md). `close()` still releases resources and stays retryable. |
| 4. Per-step state | population double buffers and their parity, raw-exchange ledgers, observations | Device-only during `step(n)` / `run(n)`; buffer swaps are host-handle swaps. | Valid on the host only after the single wait at the end of the call. |
Storage routes are separate concerns:
+56 -12
View File
@@ -1,17 +1,61 @@
# Moving circles on the uniform dense grid — S0 DRAFT contract
# Moving circles on the uniform dense grid
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.
## Status and user guide (current)
Items marked **OPEN** are undecided; items marked **CANDIDATE** are proposals.
Development-stage feature, integrated from the reviewed slices S1–S3a
(production code identical to `653fd997`). Physical evidence and its limits:
`docs/validation/moving_bodies_stage1.md`. The force used as a structural RHS
is a declared stage choice, not an established result.
*What it does.* Translating circles (with optional spin) on the dense
double-buffer D2Q9 path: per-step cut links from an FP64 geometry record
(M2 `MovingCircleLinks`), reclassification and refill of covered/uncovered
cells (M1 `MovingCircleCells`), moving-wall Bouzidi with wall velocity
`u + omega x r`, and separate per-body diagnostics.
*API.*
- `Simulation.add_body("circle", center=(x, y), radius=r,
translation=Translation(offset_box, max_speed))` before `initialize()`;
`offset_box = (x_min, x_max, y_min, y_max)` in lattice cells around the
base centre (must contain 0), `0 < max_speed ≤ .5`.
- `set_body(id, offset=(dx, dy), vx=..., vy=..., omega=...)` sets a pending
target (validated against the applied offset: box, Euclidean step ≤ .5
cell, speed, no sliding hold); `vx`/`vy` are the interval wall velocity.
- `run(1)` advances a pending target and commits it only after the step
synchronizes with `no_donor == 0`; `run(n > 1)` with a pending target
raises; without a target `run(n)` holds the geometry (spin allowed).
Direct `LBMStepper.step` refuses while a target is pending.
- `moving_diagnostics()` returns per body the columns of
`MovingTables.DIAG` (CME, Wen term, link mass, cover/uncover momentum and
mass, eq_only, no_donor, crossings) accumulated since obs were zeroed;
`obs`/`read_force` keep their contract (CME on the body).
- `get_flags()` reads device flags back when moving circles exist.
*Supported:* dense double buffer, D2Q9, FP32 unshifted, SRT/MRT, circles,
spin. *Rejected before mutation:* EsoPull, ddf shifting, D3Q19, refinement /
blocked steppers, Jets or surface actuation on a moving circle,
`sync_bodies` and checkpoints while a moving circle exists, swept bands that
touch the domain edge, other bodies, sensors or force regions.
*Failure lifecycle.* Device-side failures (upload, step, sync, no_donor)
mark the Simulation failed; later run/step/sync/checkpoint calls and direct
stepper calls refuse until `close()`.
*Known limitations.* No mass correction at curved walls (static or moving);
cover/uncover and Wen terms are diagnostics, not a selected force law; the
host loop costs a fixed per-step synchronization (≈ 1 ms/step measured on a
7.2 M-cell grid, run(1) class); outer boundaries are the single dense
channel package (see `configs/CONFIG.md`).
## Design record (staged contract, kept as written)
The sections below are the contract the slices were reviewed against
(S0 → S3a). They are kept unchanged for provenance; where they say "draft",
"proposed" or "until S2/S3", the status above supersedes them.
Original scope (Frank via TRON, thread C0C78UFKGRM/1793335557.460509):
uniform-grid prescribed translation and single-cylinder free transverse VIV;
base main `b2034f87e03656c5f67ab88562a70d6c338b9918`.
## 1. Working standards
+7
View File
@@ -82,6 +82,13 @@ Tests: `tests/unit/test_jets.py`, `test_jets_review.py`, `test_shifted_jets.py`,
- DFG 2D-1 Re 20 (V4, SRT, t* 60, pressure-NEQ outlet and local Zou-He inlet): refined and uniform
runs were NOT_STEADY under the frozen steady rule, so no interval verdict was given.
## Moving circles (dense `Simulation`, uniform grid only)
Translating circles with optional spin are supported on the dense double
buffer D2Q9 FP32 path only (not on `FixedPatchStepper` or refined levels).
API, supported/rejected modes and limits: `docs/moving_bodies.md`; stage-1
evidence: `docs/validation/moving_bodies_stage1.md`.
## Not supported in this release
3D, dynamic/automatic AMR (tagging, time-varying hierarchies beyond `remap_finest`),
+46
View File
@@ -0,0 +1,46 @@
# Moving circles — stage-1 evidence (development level)
Owner criterion (Yanqi via TRON): prescribed motion and one free VIV case both
"reasonably close" to a predeclared reference = development-stage success.
This is practical closeness, not a formal validation claim. Every number
below is from Astra's receipts as posted in Slack C0C78UFKGRM thread
1793335557.460509; raw outputs stay on Horsea under `/tmp/celeris-results-*`
and the Mac evidence directory (not in git).
## Code
- Production: `src/` identical to `653fd997` (slices S1 `983020b6`,
S2a `b920a303`, S3a `653fd997`). CPU suites and GPU hold/crossing gates
passed on those slices (G0a/G0b population-bitwise; X1 crossing ledger).
- Research drivers: `research/duetsch_re100_kc5/` (`dd37883d`, capture
`a0aea3b9`) and `research/viv_re150/` (`fb507aa1`), run as committed.
## Prescribed in-line oscillation (Dütsch et al. 1998, Re 100, KC 5)
- Protocol: `docs/validation/duetsch_re100_kc5_protocol.md`.
- Frozen settling gate (.002 F_ref waveform): **FAILED** and kept on record:
D40 domain stage (20/40/60 D, 11.25 T) and the 60 D duration diagnostic
(31.25 T, 0/28 pairs; last waveform drift .00613 F_ref).
- Unscored late-cycle CME fits c_d ≈ 2.02–2.05, c_i ≈ 1.46–1.48 vs the
paper's computed (2.09, 1.45); fields clean and symmetric (owner view).
The owner accepts ~0.61 % cycle drift at D40 for this stage; the gate is
not re-scored.
- Fig. 8 LDA profiles: qualitative only (no digitized ordinates).
## Free transverse VIV (Kolahdouz §4.1.1 setup, Qin §4.5 target)
- Predeclaration: `docs/validation/viv_re150_predeclare.md` (agreed at
`1abfe954` under its former name `leontini_viv_predeclare.md`); run at `fb507aa1`, single case, D = 40, 75 D × 60 D.
- Result: status `scored`, 336,384 steps (t U/D = 336.4); A_max = .5464 D vs
.520 D (+5.1 %, inside the declared ±10 % band); mean half peak-to-peak
.5448 D; f/f_n = .955; block drifts .78/.29/.18 %, last-five spreads
.15/.11/.10 %; finite, no_donor = 0, no seed; 337.7 s, ~1.0 ms/step.
## Limitations (binding)
- One resolution (D40) per case; no grid or time convergence.
- Structural RHS = CME + Wen is a theory-motivated stage choice; CME-only
and node-impulse variants were not run as RHS; no force law is validated.
- Boundary mismatches: NEQ outlet with backflow clamp vs traction outlet;
finite box; the strict Dütsch settling failure is unexplained.
- No mass correction at curved walls; ~1 ms/step host-coupled overhead.
@@ -1,4 +1,6 @@
# Single-cylinder free transverse VIV — PREDECLARATION (agreed at 1abfe954)
# Single-cylinder free transverse VIV, Re 150 (Kolahdouz/Qin) — PREDECLARATION (agreed at 1abfe954)
(File renamed from `leontini_viv_predeclare.md`; the Leontini Re 200 case is only the queued second case below.)
Stage goal (Yanqi via TRON): prescribed Dütsch case + one VIV case both
"reasonably close" to a predeclared reference = development-stage success, then
+23
View File
@@ -0,0 +1,23 @@
# Research VIV driver — coupling contract (not production)
`run.py` couples a 1-DOF undamped spring–mass (`viv.py`) to the production
moving-circle API. Everything specific to fluid–structure coupling lives here,
not in `src/`.
- Per step: known (y_n, v_n) and the last completed force F_last → exact
harmonic candidate (y_{n+1}, v_{n+1}) → `set_body(offset=(0, y_{n+1}),
vy=y_{n+1} − y_n)` (interval wall velocity, lattice dt = 1) → `run(1)` →
commit y, v and F_last only after the step succeeded. One step of explicit
force lag (weak coupling); no sub-iteration.
- Structural RHS = CME + Wen (y), read from `moving_diagnostics()` by the
names in `MovingTables.DIAG`; raw CME, the `obs` force and node terms are
logged as diagnostics. This is a stage choice, not a validated force law.
- Production guarantees used: fail-closed target commit (synchronized step,
`no_donor == 0`), input rejection before mutation, a failed device step
marks the Simulation failed.
- Not provided, by design: device-resident structure state or batched
(multi-step) coupling, strong/sub-iterated coupling (would need state
rollback), more than one DOF. These are future work and would be a separate
execution path, not an abstraction over this driver.
- Recorded result and limits: `docs/validation/moving_bodies_stage1.md`;
predeclaration: `docs/validation/viv_re150_predeclare.md`.
+6 -5
View File
@@ -1,4 +1,4 @@
"""One-shot predeclared VIV run (docs/validation/leontini_viv_predeclare.md).
"""One-shot predeclared VIV run (docs/validation/viv_re150_predeclare.md).
Re 150, U* = 5, m/(rho D^2) = 2, zeta = 0, D = 40 on a 75 D x 60 D native
channel (uniform inlet, NEQ outlet, free-slip y walls), cylinder at
@@ -65,7 +65,8 @@ def main():
step_cap = args.steps
centre = (30.0 * Dn, 0.5 * (ny - 1))
from CelerisLab.body.motion import Translation
from CelerisLab.body.motion import MovingTables, Translation
col = {k: MovingTables.DIAG.index(k) for k in ("cme_y", "wen_y", "cover_jy", "uncover_jy", "no_donor")}
from CelerisLab.simulation import Simulation
with tempfile.NamedTemporaryFile("w", suffix=".json", delete=False) as fh:
json.dump(cfg, fh)
@@ -95,9 +96,9 @@ def main():
sim.run(1) # device failures raise and mark the Simulation failed
d = sim.moving_diagnostics()[0]
f_obs = sim.read_force(body)[1]
f_rhs = float(d[1] + d[3]) # CME + Wen (y)
row = (n, y1, v1, f_rhs, float(d[1]), float(d[3]), float(f_obs),
float(d[7] - d[10]), float(d[13]))
f_rhs = float(d[col["cme_y"]] + d[col["wen_y"]]) # CME + Wen (y)
row = (n, y1, v1, f_rhs, float(d[col["cme_y"]]), float(d[col["wen_y"]]), float(f_obs),
float(d[col["cover_jy"]] - d[col["uncover_jy"]]), float(d[col["no_donor"]]))
if not viv.finite(row): # every recorded output, not only the RHS
status = "fail_nonfinite_output"
break
+1 -1
View File
@@ -1,6 +1,6 @@
"""Pure logic of the predeclared Re 150 free transverse VIV case (stdlib only).
Contract: docs/validation/leontini_viv_predeclare.md. Oscillator update with a
Contract: docs/validation/viv_re150_predeclare.md. Oscillator update with a
constant force over one lattice step, response-cycle tracking with fixed y = 0
crossings, eligibility, 10-cycle blocks, settling/scoring with deterministic
restart, and the conditional seed.
+3 -61
View File
@@ -1,10 +1,10 @@
# CelerisLab/body/motion.py
"""Host state of a translating circle (dense double buffer, D2Q9 only).
"""Translation targets and device tables (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.
separate facts; centres are always ``base + offset``. MovingCircle owns host
targets, MovingTables owns device records, and Simulation commits completed steps.
"""
from __future__ import annotations
@@ -18,10 +18,6 @@ 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)
@@ -98,23 +94,6 @@ class MovingCircle:
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:
@@ -138,43 +117,6 @@ class MovingCircle:
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``.
+14
View File
@@ -86,6 +86,20 @@ step/one_step_*.cu → kernel 编排
- 转速 `omega` 属于运行时状态,可在步进期间更新;该更新不改变 `N_OBJS` 或 `OBS_*` 内存布局,因此**不需要**因 `omega` 变化触发重新编译。
- 力矩观测从 `obs` 的 torque segment 读取;force / torque / sensor 段统一由 `config_obs.h` 宏描述。
### 外边界能力(Outer boundaries, dense `Simulation`)
The dense path supports one channel package, not four independently
selectable sides:
- left x = 0: velocity inlet (`method.inlet.scheme`, `profile`, `physics.velocity`);
- right x = nx−1: outlet (`method.outlet.mode`; `neq_extrap` holds rho = 1 with
a `max(u, 0)` backflow clamp when `backflow_clamp` is set);
- bottom and top: one shared `y_wall_bc` (`bounce_back` halfway or
`free_slip`); corners are walls.
`y_wall_bc = "on_site_velocity"` switches to x-periodic with on-site y planes
and requires an object-free domain. There is no x/y-periodic or symmetry
option for domains with bodies; adding outer BC types is a separate design
question, not covered by these keys.
### 稳定性默认窗口(Stability defaults)
- 默认 `method.omega_guard = {min: 0.01, max: 1.99}`,用于约束碰撞频率与 LES 有效粘度路径。
+5 -3
View File
@@ -21,7 +21,7 @@ class LBMStepper:
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).
# without attached moving tables, a pending target, or a failed step).
self.blocked: str | None = None
self.moving = None # MovingTables attached by the owner (translating circles)
@@ -39,8 +39,10 @@ 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")
# Moving-circle kernels exist only in D2Q9 double-buffer modules.
moving = not self._esopull and cfg.is_d2q9
self.moving_fn = module.get_function("MovingCircleLinks") if moving else None
self.moving_cells_fn = module.get_function("MovingCircleCells") if moving else None
tpb = cfg.threads_per_block
self.block = (tpb, 1, 1)
+16 -8
View File
@@ -28,7 +28,7 @@ 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).
# Translating circles step only through their moving tables (attached at initialize).
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); "
@@ -122,8 +122,9 @@ class Simulation:
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.
circle (before initialize only; docs/moving_bodies.md).
Positive-step runs are rejected until ``initialize()`` has
attached its moving tables.
**kwargs: passed to the geometry constructor.
Returns:
@@ -241,7 +242,8 @@ class Simulation:
self._require_usable()
tables = self.stepper.moving
if tables is None:
return np.zeros((0, 6), np.float32)
from .body.motion import MovingTables
return np.zeros((0, MovingTables.DIAG_FLOATS), np.float32)
self.stream.synchronize()
return tables.download_diagnostics()
@@ -360,8 +362,11 @@ class Simulation:
Args:
id: body_id from add_body().
omega: rotation speed (lattice units). ``None`` = no change.
vx: translation velocity x (future). ``None`` = no change.
vy: translation velocity y (future). ``None`` = no change.
vx: translation velocity x. For a translating circle: the
interval wall velocity of the target, only together with
``offset``; otherwise stored in the body state (not applied
to stationary walls). ``None`` = no change.
vy: translation velocity y; same rules as ``vx``.
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
@@ -781,6 +786,7 @@ class Simulation:
(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()]
pending = any(motion.pending is not None for motion, _ in circles)
try:
advancing = tables.upload_records(circles, stream)
self.stepper._step(steps, action_gpu=self.bodies.action_gpu,
@@ -793,10 +799,12 @@ class Simulation:
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")
if pending: # a target equal to the applied offset is consumed as a hold step
for motion, _ in circles:
motion.commit()
tables.upload_records(circles, stream) # first hold: old = new = applied, u = 0
stream.synchronize()
if advancing:
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
+16 -13
View File
@@ -11,6 +11,20 @@ _spec = importlib.util.spec_from_file_location(
p = importlib.util.module_from_spec(_spec)
_spec.loader.exec_module(p)
def _load(rel, name):
"""Load a research driver module under a unique name (no bare ``import run``)."""
import sys
path = Path(__file__).resolve().parents[2] / rel
sys.path.insert(0, str(path.parent)) # for the driver's own sibling imports
try:
spec = importlib.util.spec_from_file_location(name, path)
mod = importlib.util.module_from_spec(spec)
spec.loader.exec_module(mod)
finally:
sys.path.remove(str(path.parent))
return mod
T = 2520
@@ -93,12 +107,7 @@ def test_exact_point_inside_body_is_missing_even_with_fluid_corners():
def test_generated_config_passes_the_real_loader(tmp_path):
import json
import sys
sys.path.insert(0, str(Path(__file__).resolve().parents[2] / "research/duetsch_re100_kc5"))
try:
import run
finally:
sys.path.pop(0)
run = _load("research/duetsch_re100_kc5/run.py", "duetsch_run")
from CelerisLab.config import load_lbm_config
cfg_dict, u_max, n = run.config(20, 10.0, 2520)
path = tmp_path / "config_lbm.json"
@@ -111,13 +120,7 @@ def test_generated_config_passes_the_real_loader(tmp_path):
def _analyze():
import sys
sys.path.insert(0, str(Path(__file__).resolve().parents[2] / "research/duetsch_re100_kc5"))
try:
import analyze
finally:
sys.path.pop(0)
return analyze
return _load("research/duetsch_re100_kc5/analyze.py", "duetsch_analyze")
def test_hf_energy_is_parseval_weighted():
+138 -33
View File
@@ -1,14 +1,13 @@
"""Translating-circle host contract, S1 (docs/moving_bodies.md 4a, 5, 7).
"""Translating-circle host contract (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.
Real Simulation entrypoints with device doubles; no CUDA context or GPU work.
"""
import math
from types import SimpleNamespace
import pytest
from CelerisLab.body.motion import (MAX_STEP, R_MIN, MovingCircle, Translation,
enumerate_donor_coverage)
from CelerisLab.body.motion import MAX_STEP, R_MIN, MovingCircle, Translation
BAD = (float("nan"), float("inf"))
@@ -46,46 +45,96 @@ def test_targets_validate_against_applied_and_last_wins():
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))
def test_record_commit_and_hold(monkeypatch):
sim, writes, c = routed_advancing_simulation(monkeypatch)
assert c.record() == ((50., 40.), (50., 40.), (0., 0.))
sim.set_body(0, offset=(0., .25), vy=.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))]
seen_at_sync = []
sim.run(1, stream=SimpleNamespace(synchronize=lambda: seen_at_sync.append(c.applied)))
assert seen_at_sync == [(0., 0.), (0., .25)] # commit between step sync and hold sync
uploads = [w[1][0] for w in writes if isinstance(w, tuple) and w[0] == "upload"]
assert uploads == [((50., 40.), (50., 40.25), (0., .05)),
((50., 40.25), (50., 40.25), (0., 0.))]
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
sim.run(1, stream=STREAM) # no new target: hold step
uploads = [w[1][0] for w in writes if isinstance(w, tuple) and w[0] == "upload"]
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
sim.set_body(0, offset=(0., 0.), vy=-.05) # zero offset is a real target
assert c.pending == ((0., 0.), (0., -.05))
c.advance(uploads.append, lambda: None, lambda: None)
sim.run(1, stream=STREAM)
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 test_failure_after_upload_attempt_is_fail_closed(monkeypatch, stage):
sim, writes, c = routed_advancing_simulation(monkeypatch, fail="upload" if stage == "upload" else None)
sim.set_body(0, offset=(0., .25), vy=.05)
def call(name):
def run(*_):
if name == stage:
raise RuntimeError(f"{name} failed")
return run
def fail(*a, **k):
raise RuntimeError(f"{stage} failed")
if stage == "step":
sim.stepper._step = fail
stream = SimpleNamespace(synchronize=fail) if stage == "synchronize" else STREAM
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)):
sim.run(1, stream=stream)
assert c.failed and c.applied == (0., 0.) and c.pending == ((0., .25), (0., .05))
assert sim._failed is not None and sim.stepper.blocked == sim._failed
for later in (c.record, lambda: sim.set_body(0, offset=(0., .1), vy=0.)):
with pytest.raises(RuntimeError, match="unknown"):
later()
with pytest.raises(RuntimeError, match="Simulation failed"):
sim.run(1, stream=STREAM)
with pytest.raises(RuntimeError, match="moving circle step failed"):
sim.stepper.step(1, action_gpu=None, obs_gpu=None)
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)
# 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 _inside(x: int, y: int, c, r: float) -> bool:
return (x - c[0]) ** 2 + (y - c[1]) ** 2 < r * r
def empty_donor_cells(c_old, c_new, 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
@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.
@@ -278,6 +327,26 @@ def advancing_simulation(monkeypatch, **tables):
return sim, writes, sim.bodies.objects[0].motion
def routed_advancing_simulation(monkeypatch, **tables):
"""Real body registration, target routing, run and owner hook; double device I/O."""
sim = routed_simulation(monkeypatch)
bid = sim.add_body("circle", (50., 40.), 10., translation=Translation((-2., 2., -3., 3.), .1))
sim._initialized = True
sim._assert_runtime_contracts = lambda: None
writes = []
monkeypatch.setattr(sim.bodies, "_upload_action_async", lambda s: writes.append("action"))
monkeypatch.setattr(sim.bodies, "zero_obs_async", lambda s: writes.append("obs"))
monkeypatch.setattr(sim.bodies, "download_obs_full_async", lambda s: writes.append("download"))
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)
sim._moving_tables = sim.stepper.moving
sim._attach_moving()
return sim, writes, sim.bodies.get(bid).motion
def test_run_one_advances_and_commits_after_sync(monkeypatch):
sim, writes, motion = advancing_simulation(monkeypatch)
motion.set_target((0., .25), (0., .05))
@@ -297,6 +366,44 @@ def test_run_one_advances_and_commits_after_sync(monkeypatch):
assert ("step", 3, False) in writes and motion.applied == (0., .25)
def test_zero_motion_target_is_consumed_by_run_one(monkeypatch):
sim, writes, motion = advancing_simulation(monkeypatch)
motion.set_target((0., 0.), (0., 0.)) # equal to applied: a hold target
sim._attach_moving()
assert sim.stepper.blocked == module.PENDING_BLOCKED
sim.run(1, stream=STREAM)
assert motion.pending is None and motion.applied == (0., 0.)
assert sim.stepper.blocked is None and ("step", 1, False) in writes
writes.clear()
sim.run(3, stream=STREAM) # no stale latch
assert ("step", 3, False) in writes
@pytest.mark.parametrize("streaming,d2q9,bound", [("double_buffer", True, True),
("double_buffer", False, False),
("esopull", True, False)])
def test_moving_kernels_bound_only_for_d2q9_double_buffer(monkeypatch, streaming, d2q9, bound):
from CelerisLab.lbm import stepper as stepper_module
monkeypatch.setattr(stepper_module.cuda, "mem_alloc", lambda n: object())
moving = {"MovingCircleLinks", "MovingCircleCells"}
def get_function(name):
if name in moving and not (d2q9 and streaming == "double_buffer"):
raise KeyError(name) # absent from D3Q19 / EsoPull modules
return lambda *a, **k: None
cfg = SimpleNamespace(streaming=streaming, threads_per_block=128, is_d2q9=d2q9,
nx=16, ny=8, nz=4, nq=9 if d2q9 else 19)
stepper = stepper_module.LBMStepper(SimpleNamespace(), SimpleNamespace(get_function=get_function), cfg)
assert (stepper.moving_fn is not None) == bound and (stepper.moving_cells_fn is not None) == bound
def test_moving_diagnostics_empty_shape_matches_layout(monkeypatch):
from CelerisLab.body.motion import MovingTables
sim = make_simulation(monkeypatch)
sim._failed, sim.stepper = None, SimpleNamespace(moving=None)
assert sim.moving_diagnostics().shape == (0, MovingTables.DIAG_FLOATS)
@pytest.mark.parametrize("failure", ("upload", "first_sync", "no_donor", "hold_upload",
"hold_sync"))
def test_failed_advance_fails_simulation(monkeypatch, failure):
@@ -507,23 +614,21 @@ def test_moving_circle_excluded_from_static_links():
def test_failed_moving_step_marks_simulation_failed(monkeypatch):
sim = routed_simulation(monkeypatch)
sim, writes, motion = routed_advancing_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))
sim.set_body(0, offset=(0., .2), vy=.05)
def fail(*_):
raise RuntimeError("sync failed")
with pytest.raises(RuntimeError, match="sync failed"):
motion.advance(lambda *_: None, lambda: None, fail)
sim.run(1, stream=SimpleNamespace(synchronize=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)
sim.set_body(0, offset=(0., .1), vy=.05)
with sim.ctx.activate(): # a refused close keeps the failed state and resources
with pytest.raises(RuntimeError, match="active CUDA"):
sim.close()
+15 -6
View File
@@ -8,6 +8,20 @@ _spec = importlib.util.spec_from_file_location(
viv = importlib.util.module_from_spec(_spec)
_spec.loader.exec_module(viv)
def _load(rel, name):
"""Load a research driver module under a unique name (no bare ``import run``)."""
import sys
path = Path(__file__).resolve().parents[2] / rel
sys.path.insert(0, str(path.parent)) # for the driver's own sibling imports
try:
spec = importlib.util.spec_from_file_location(name, path)
mod = importlib.util.module_from_spec(spec)
spec.loader.exec_module(mod)
finally:
sys.path.remove(str(path.parent))
return mod
D, U = 40, 0.04
@@ -72,12 +86,7 @@ def test_seed_reset_discards_incomplete_cycle_and_proof():
def test_viv_config_passes_the_real_loader(tmp_path):
import json
import sys
sys.path.insert(0, str(Path(__file__).resolve().parents[2] / "research/viv_re150"))
try:
import run
finally:
sys.path.pop(0)
run = _load("research/viv_re150/run.py", "viv_run")
from CelerisLab.config import load_lbm_config
cfg_dict, nx, ny = run.config()
path = tmp_path / "config_lbm.json"