Compare commits

...
Author SHA1 Message Date
ClaudeandClaude Opus 5.5 0f8fdbc30e refinement_stepper: compatibility views are kept on purpose (comment only)
The views and the PeriodicPatchStepper alias stay for the existing tests
and docs (about 130 references); the comment no longer says REMOVE AFTER.
No code change.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 02:26:24 +00:00
ClaudeandClaude Opus 5.5 e4538219b2 docs: supported settings cite guards by file and symbol at 888c8038
The file:line citations pointed at c95a7ba and went stale after P1, the
compile-time root modes and the kernel entry/step split. Each guard is now
named by file plus function (or the kernel #error guard) at source
revision 888c8038. No claim, scope, observation or evidence revision
changed.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 02:24:37 +00:00
ClaudeandClaude Opus 5.5 888c80381c kernels (K10): refinement.cu becomes an include-only entry; patch kernels under step/
refinement.cu keeps its header, includes, the curved_double.cu entry
block, the #error guard and the optional CELERIS_REFINEMENT_AUDIT
counters, then includes step/patch_step.cu, which holds patch_collide,
PatchExplode, PatchFineStep and PatchCoarseStep moved verbatim (checked
against the previous refinement.cu) in the same order.
test_refinement_mrt reads refinement.cu + step/patch_step.cu; its
assertions are unchanged. The audit-build tests include refinement.cu
from the build root, so step/patch_step.cu still resolves.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 02:03:28 +00:00
ClaudeandClaude Opus 5.5 bf5088f84c kernels (K9): blocks.cu becomes an include-only entry; kernels under step/
blocks.cu keeps its header, includes and #error guards and then includes
the kernels, as kernel_v2.cu does:
- step/block_transfer.cu: BlockHaloCopy, BlockExplode,
  BlockExplodeProjected, BlockCoalesce
- step/block_boundary.cu: BlockCurvedBoundary, BlockAccumulateExchange
- step/block_step.cu: BlockStep (rejects BLOCK_ROOT_BC == 2)
- step/block_channel.cu: unchanged (mode 2 only)
Kernel text moved verbatim (checked against the previous blocks.cu) and in
the same definition order, so the translation unit is the same sequence of
kernels. The two source-reading tests follow the move:
test_projected_explosion reads step/block_transfer.cu, test_refinement_mrt
reads blocks.cu + step/block_step.cu; their assertions are unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 01:46:49 +00:00
ClaudeandClaude Opus 5.5 af406b44bf Merge K4 helper withdrawal into K5
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 01:32:42 +00:00
ClaudeandClaude Opus 5.5 78c0f1811b kernels (K4): withdraw the address helpers; keep only fill_invalid
Horsea's gates showed block_ddf/block_slot change real PTX in
BlockHaloCopy, BlockExplode and BlockExplodeProjected (6/10 modules;
BlockExplodeProjected 38 -> 40 registers) although the results stayed
exact. The stack's bar is byte-identical PTX for naming-only steps, so
those three kernels return to their K3 expressions, the helpers leave
core/layout.cuh, and test_projected_explosion goes back to its K3 form.
fill_invalid stays (its kernels were byte-identical). K4 vs K3 is now
only the NaN-fill helper.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 01:32:40 +00:00
ClaudeandClaude Opus 5.5 0ddc6c3d62 Merge K4 test harness fix into K5
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:35:09 +00:00
ClaudeandClaude Opus 5.5 43e226c59a test_projected_explosion: give the host harness the real block address helpers
K4 made BlockExplodeProjected address storage through block_ddf and
block_slot (core/layout.cuh). This host test compiles the kernel text
extracted from blocks.cu, so it now also extracts those two helper
definitions from core/layout.cuh (device qualifiers -> inline) instead of
failing to compile. Assertions and tolerances unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:35:07 +00:00
ClaudeandClaude Opus 5.5 93e692a40f kernels (K5): reuse the existing bounce_back_swap instead of a new copy
operators/helpers.cuh already had bounce_back_swap (the uniform steps use
it). It moves verbatim, #pragma unroll included, into its own
operators/bounce_back_swap.cuh, which helpers.cuh includes; the channel
root includes that header and calls it. bounce_back_swap_d2q9 from the
previous commit is removed. bounce_back.cuh keeps only the half-way rows
helper.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:32:07 +00:00
ClaudeandClaude Opus 5.5 33837697da kernels (K5): channel root reuses boundary and streaming helpers
- boundary/bounce_back.cuh: apply_wall_bb_d2q9_rows(y, top_row, ...) is
  the existing half-way y-wall correction with an explicit top row;
  apply_wall_bb_d2q9_y_pull keeps its signature and calls it with NY - 2.
  bounce_back_swap_d2q9(f) names the full opposite-pair swap.
- streaming/pull_double_buffer.cuh: load_node(fi, n, f) loads one node's
  stored populations.
- step/block_channel.cu uses all three instead of its inline copies: the
  wall swap, the outlet donor gather, and rows 1 / ny-2 (now one call with
  the physical top row; same y == 1 / y == top_row order and loads).
Same statements, branch order, indices and arithmetic.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:31:26 +00:00
ClaudeandClaude Opus 5.5 b47af12c9d kernels (K4): blocked storage address helpers and one invalid-cell fill
- core/layout.cuh: block_ddf(f, block_cell[, cells]) returns the start of a
  block's populations and block_slot(block_cell, cells, q) the offset of one
  population; BlockHaloCopy, BlockExplode and BlockExplodeProjected use them
  instead of repeating (c / cells) * NQ * cells [+ q * cells + c % cells].
  BlockExplode keeps its loop-hoisted parent base; types stay uint64_t /
  unsigned long; the expressions are the same ones, now named.
- core/cell_kind.cuh: fill_invalid(f) writes the NaN marker; BlockStep,
  BlockChannelStep and PatchFineStep call it instead of three copies of the
  loop.
No change to operation order, types, launches or arithmetic.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:29:03 +00:00
ClaudeandClaude Opus 5.5 a3df5855a2 kernels (K3): cell kinds in core/, capability guard first, self-contained block collision
- core/cell_kind.cuh holds KIND_ACTIVE/INNER/OUTER (was a constexpr in
  blocks.cu that step/block_channel.cu and the channel_root.cpp shim
  relied on or redefined); blocks.cu, block_channel.cu and the shim now
  include it.
- blocks.cu checks the lattice/model #error guard before the
  BLOCK_ROOT_BC guard and mode includes, so an unsupported build reports
  the real reason first.
- operators/block_collision.cuh includes what it uses (params, macro,
  SRT/MRT collision) instead of relying on blocks.cu's include order.
Preprocessor-only: same kernels, same statements, same arithmetic.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-22 00:26:41 +00:00
Frank14f 11e96a188b Test historical build fingerprints with their pinned renderer 2026-10-21 18:50:15 +08:00
Frank14f 3a374faace Specialize blocked root boundary mode at compile time 2026-10-21 18:35:27 +08:00
Grok fce24a64e8 test: run history-independent identity checks before the history lookup 2026-10-21 16:49:41 +08:00
Grok d437faa501 test: check historical source manifests against their recorded commits
Manifests and receipts unchanged. Each of the 8 group-A identity tests now
extracts its manifest's recorded commit and verifies hashes there instead of
against the current checkout. Tamper assertions and runner source refusals are
kept; missing recorded history is an explicit failure, not a skip or repin.
2026-10-21 16:47:22 +08:00
Frank14f 26e49451ec Merge reviewed checkpoint additions for release candidate 2026-10-21 15:20:04 +08:00
Frank14f e214480452 Merge reviewed docs additions for release candidate 2026-10-21 15:20:04 +08:00
Frank14f c414f6c7f6 Merge reviewed fixtures additions for release candidate 2026-10-21 15:20:04 +08:00
ClaudeandClaude Opus 5.5 decd419e94 docs: cite only the recorded MRT parity pass (K1 7813334)
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 07:08:31 +00:00
ClaudeandClaude Opus 5.5 1c47cb8570 docs: supported settings review fixes (exact parity test IDs, divisibility guard)
Names the exact SRT/MRT parity tests and their sub-case scope (SRT: one
of 18 sub-cases, unidentified; MRT: 12 sub-cases, default outlet settings,
with the runs that recorded the pass), cites blocks.py:134 for tile
divisibility, labels refinement.py:457 as the fine-boundary re-check, and
says 'measured at 36ffc8b1' for the Taylor-Green observation.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 07:05:34 +00:00
ClaudeandClaude Opus 5.5 ee7074a53d tests: restore the negative checks masked by stale fixtures (WP-B groups B-E)
Test-side only; no production, manifest or pin change. Group A identity
checks (test_source_pin and the jets.py pins) are unchanged and stay
failing as recorded.

B test_refined_dfg: receipt/plan fixtures load the frozen manifest JSON
  that verify_source returns for a matching tree (frozen_manifest), and
  the two manifest-selection assertions call verify_source with only its
  working-tree hash comparison neutralised (selected_manifest), so the
  malformed/collision/observation-mode/checkout receipt checks run again.
  test_source_pin remains the single identity check.
C test_channel_root: int(owner)-1 (uint64-int became float on NumPy 1);
  test_shifted_storage: FP32 NaN fill so the arrays stay FP32. Both reach
  the intended immutable/representation/FP32/density checks.
D test_block_stepper fake gets _ddf_shifting; test_fine_surface_gpu sets
  _boundary_hosts/_level_jets instead of the read-only fine_boundary view;
  test_refinement_followup fake_backend sets _layouts/_ddf_shifting/
  _level_cfgs instead of the read-only fine_layout view (as reviewed in
  32da2d73).
E test_kernel_contracts_reject_before_upload: an inconsistent pair tau is
  now first rejected by the hierarchy chain check; the test asserts that
  rejection, then bypasses _check_chain to prove the per-level omega/tau
  guard still rejects, and asserts no mem_alloc/memcpy_htod in either case.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 07:04:47 +00:00
Codex b97fc79af9 test: stage existing same-scene body checkpoint roundtrip
Copy tests/integration/test_checkpoint_bodies_gpu.py byte-for-byte from 34b30c434f. Base c95a7baae4. No production or dependency changes. GPU execution is deferred to the assigned executor.
2026-10-21 15:04:08 +08:00
ClaudeandClaude Opus 5.5 8e9fc72e19 docs: curved-wall observation states only the recorded matrix outcome
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:59:14 +00:00
ClaudeandClaude Opus 5.5 d2cd1114ed docs: supported settings for the 2D Jet + multilevel candidate
One page listing supported settings, each with the guard enforcing it
(file:line at c95a7ba) and the existing tests exercising it. Measured
observations are labelled by their exact tested case (SRT strict-parity
1/4320 values ~2 ULP with unresolved cause; V3 TG constant explosion
aggregate; curved-wall r1 pilot; V4 NOT_STEADY), not stated as general
bounds. Documentation only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:59:02 +00:00
ClaudeandClaude Opus 5.5 c95a7baae4 test_projected_pairs: regenerate pinned build sources from NUMERICAL_REVISION kernels
The pinned BUILD_SOURCES are the a135252 qualification build inputs.
generate_config hashes the whole kernel snapshot, so after K1/K2 the
working-tree kernels (legitimately) no longer reproduce them. The test
now archives src/CelerisLab/lbm/kernels at NUMERICAL_REVISION and points
compiler_v2._KERNEL_DIR there; config headers are still rendered by
today's compiler. BUILD_SOURCES and the equality assertion are unchanged;
each size is a subTest; a missing NUMERICAL_REVISION fails (no skip).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:26:45 +00:00
ClaudeandClaude Opus 5.5 f4372a84c0 Merge a454f023 (claude/jet-shapes-once) for integration review
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:03:36 +00:00
ClaudeandClaude Opus 5.5 e48498d089 Merge 580fcf65 (grok/p1-docs-wording) for integration review
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:03:36 +00:00
ClaudeandClaude Opus 5.5 531739454e Merge 1cbdc479 (claude/step-single-sync) for integration review
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:03:36 +00:00
ClaudeandClaude Opus 5.5 b93de2ce24 Merge 1fa498c7 (claude/k2-block-coords) for integration review
Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 06:03:36 +00:00
ClaudeandClaude Opus 5.5 a454f02355 Jets: compute jet link shapes once per boundary preparation
surface_link_velocity already runs jet_link_shapes; CurvedLinkSoA.
assign_surface and _prepare_blocked_boundary then ran it again with the
same jets and hits to fill the host jet cache. A private _surface_links
returns (velocity, aperture, shape, jet_id) from the one call; both
callers use it. surface_link_velocity keeps its two-value public
contract (now a thin wrapper). Same values, same validation order, same
cache contents; host geometry preparation only.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 05:57:11 +00:00
ClaudeandClaude Opus 5.5 1cbdc479cf test_block_step_sync: locate block_stepper without importing it
The module-level import of CelerisLab.lbm.block_stepper loaded
pycuda.driver into the shared test process, unlike the other CPU host
tests (which load mocked copies via AST). find_spec gives the source path
without the import. Test content and assertions unchanged.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:39:52 +00:00
ClaudeandClaude Opus 5.5 1fa498c7d7 kernels: shared block-window to global coordinate helper
BlockStep (slot-reservoir wall rows) and BlockChannelStep computed the
global cell of a block-window cell with the same four lines.
block_global_xy (core/layout.cuh, 2D, __forceinline__) now holds them;
both kernels call it with their own domain width. Same integer
operations: no numerical, ABI or scheduling change. The slot-reservoir
kernel layout is unchanged (no split).

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:24:27 +00:00
ClaudeandClaude Opus 5.5 781333407e kernels: one shared unforced MRT collision for blocked, channel and patch paths
BlockStep, BlockChannelStep and refinement.cu patch_collide each carried
the same compute_rho_u -> omega clamp -> collide_mrt_no_force sequence.
collide_mrt_unforced (collision_mrt.cuh, __forceinline__) now holds it
once; the three sites call it. Same operations in the same order: no
numerical, ABI or SRT change. The host source contract now requires one
helper definition and exactly one call per file, and forbids the old
inline sequence.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-21 01:15:47 +00:00
43 changed files with 881 additions and 389 deletions
+88
View File
@@ -0,0 +1,88 @@
# Supported settings: 2D Jet + fixed multilevel refinement (candidate 888c8038)
Every entry names the guard that enforces it (file and function or kernel guard,
at source revision 888c80381c006e976be2af8a452288729fc17c1b) and, where one
exists, an existing test that exercises it. A setting not listed here is not
claimed. Observations at the end are tied to the exact case in which they were
measured; they are not general bounds.
## Lattice, storage, collision
| Setting | Supported | Guard |
|---|---|---|
| Dimension / lattice | 2D D2Q9 only (3D refinement is future work) | `lbm/refinement.py` `FixedPatch2D.from_config`; `#error` guards in `lbm/kernels/blocks.cu` and `lbm/kernels/refinement.cu` |
| Precision / streaming | FP32 storage, double buffer | same guards |
| Collision on refinement levels | SRT or MRT; TRT and LES rejected | `lbm/refinement.py` `FixedPatch2D.from_config` (re-checked for fine boundaries in `prepare_fine_boundary`) |
| Shifted (ddf_shifting) storage | only when every level is blocked; dense or mixed rejected | `lbm/refinement_stepper.py` `FixedPatchStepper.__init__`, `FixedPatchStepper._check_storage`; `#error` guard in `lbm/kernels/refinement.cu` (dense patch kernel unshifted only) |
| Fine relaxation | 1/tau_fine must lie inside the configured omega guards | `lbm/refinement.py` `FixedPatch2D.from_config` |
| Uniform blocked backend `PeriodicBlockStepper` | doubly periodic, SRT only, stationary passive circles | `lbm/block_stepper.py` `PeriodicBlockStepper.__init__`, `_static_circles` |
Tests: `tests/unit/test_refinement_mrt.py`, `tests/validation/multigrid/test_refinement_mrt_gpu.py`,
`tests/unit/test_shifted_storage.py`, `tests/unit/test_shifted_block_operators.py`.
## Refinement hierarchy (`FixedPatchStepper(cfg, *, context, refinements, ...)`)
| Setting | Supported | Guard |
|---|---|---|
| Depth | at least one refinement level; one entry per level in `level_blocks`/`level_objects`/`level_jets` | `lbm/refinement_stepper.py` `_check_levels` |
| Single L1 rectangle | two-coarse-cell margin to the domain edge | `lbm/refinement.py` `FixedPatch2D.from_config` |
| Several L1 rectangles (union) | require `level_blocks[0]` and `level_blocks[1]`; transfers must close (dependency checks, no fixed gap rule) | `lbm/refinement_stepper.py` `_check_levels`; `lbm/refinement.py` `_check_parent_pull_closure` |
| Nested levels | depth >= 3 must be all blocked; a blocked level requires blocked levels below; a scalar child needs a scalar parent | `lbm/refinement_stepper.py` `_check_levels`; `lbm/hierarchy.py` `build_hierarchy` |
| Nested transfer support | child only over inherited ACTIVE cells, donors/receivers uncovered ACTIVE, no child over the parent restriction footprint | `lbm/refinement.py` `_nested_kind` |
| Tile shape | positive `(by, bx)`; the level allocation must be divisible by the tile (no partial tiles yet) | `lbm/refinement_stepper.py` `_check_levels`, `lbm/refinement.py` `FixedPatchSet2D.fine_layout` (shape); `lbm/blocks.py` `BlockLayout2D.build` (divisibility) |
| Coarse boundary | `periodic`, `slot_reservoir` or `channel` | `lbm/refinement_stepper.py` `FixedPatchStepper.__init__` |
| `slot_reservoir` | on-site velocity bottom slot, reservoir top, velocity 0, 2 <= slot width < nx, patch rows 4..ny-4 | `lbm/refinement_stepper.py` `FixedPatchStepper.__init__`; `lbm/refinement.py` `prepare_fine_boundary` |
| `channel` | all levels blocked; Zou-He inlet, pressure-NEQ outlet, stationary bounce-back walls; transfers away from the closure | `lbm/channel.py` `validate_channel_config`, `validate_channel_transfers` |
| Explosion | `constant` (default); `projected_linear` only periodic, all-blocked SRT, no bodies/Jets, no LES (experimental) | `lbm/refinement_stepper.py` `FixedPatchStepper.__init__` |
| L2 runtime remap (`remap_child`) | periodic three-level L2 sets only, no L0/L1 bodies, identical pinned boundary | `lbm/refinement_stepper.py` `FixedPatchStepper.remap_child` |
Tests: `tests/unit/test_refinement_three_level_multipatch.py`, `tests/validation/multigrid/test_refinement_multipatch.py`,
`test_refinement_three_level.py`, `test_refinement_three_level_multipatch_gpu.py`, `test_refinement_blocked_base.py`,
`test_refinement_blocked_fine.py`, `test_channel_root_gpu.py`, `tests/unit/test_channel_root.py`,
`tests/unit/test_projected_pair_stepper.py`, `test_runtime_child_remap_gpu.py`, `test_pinned_body_remap_gpu.py`.
## Bodies and Jets on refined levels
| Setting | Supported | Guard |
|---|---|---|
| Bodies | stationary circles without sensor/force-region roles, entirely inside the patch union, fully resolved, disjoint, stencil away from the two-hop transfer footprint | `lbm/refinement.py` `prepare_fine_boundary` |
| Level | not on L0; only on the finest level of a nested hierarchy | `lbm/refinement_stepper.py` `_check_levels` |
| Jet condition | velocity only; `pressure` is reserved and raises NotImplementedError | `body/jets.py` `Jet.__post_init__` |
| Jet geometry | each Jet selects resolved links of exactly one stationary body; no Jet/Jet or Jet/SurfaceActuation overlap; abs(amplitude) < 1 lattice | `body/jets.py` `Jet.__post_init__`, `jet_link_shapes`, `_surface_links` |
| Jets with shifted storage | only in the all-blocked stepper with SRT and no LES | `body/jets.py` `_surface_links` |
| SurfaceActuation | unshifted storage only | `body/objects.py` `SimObject.validate_surface_actuation` |
| Raw boundary exchange (force) | blocked body levels only; the window ledger requires `observe_raw_exchange` | `lbm/refinement_stepper.py` `FixedPatchStepper.__init__` |
| Checkpoint with Jets | refused | `simulation.py` `Simulation.save_checkpoint` |
Raw exchange is the cut-link momentum reaction (wall + aperture), not a validated
nozzle thrust or stress load. Jet amplitude is a lattice velocity, not a flow rate.
Tests: `tests/unit/test_jets.py`, `test_jets_review.py`, `test_shifted_jets.py`, `test_raw_exchange.py`,
`tests/validation/multigrid/test_refinement_jets.py`, `test_shifted_jets_gpu.py`, `test_raw_exchange_gpu.py`,
`test_exchange_window.py`, `test_fine_surface_gpu.py`.
## Observations from specific tested cases (not general bounds)
- SRT channel root vs uniform OneStep, strict-parity test
`tests/validation/multigrid/test_channel_root_gpu.py::ChannelRootGPU::test_uniform_matched_root`:
in one of its 18 sub-cases (2 storage modes x 3 profile/clamp/damp settings x 3 blocks; which
one is not recorded), 1 of the 4320 compared values differs, max abs 2.33e-10, max rel 1.77e-7
(about 2 ULP). The test remains failed; its cause and any general bound are unresolved.
- MRT root parity test
`tests/validation/multigrid/test_refinement_mrt_gpu.py::RefinementMRTGPU::test_uniform_matched_mrt_root`
(12 sub-cases: 2 storage modes x 2 profiles x 3 blocks, default outlet settings) passed bitwise
in the Horsea focused GPU run at K1 7813334 (2026-10-21).
- Periodic Taylor-Green Re 100, measured at 36ffc8b1 (V3 matched pair, N128 depth 2):
every per-sample screen passed; the constant explosion's aggregate worst velocity L2,
0.004418, exceeds its preregistered 0.002 budget, while projected_linear passed that budget.
- Curved-wall pressure (pilot harness, exact-q analytic annulus, tau 0.8): the coarse swirl r1
case exceeded its terminal density-gain cap (1.105e-4 vs 1e-4), so the frozen 9-case matrix
is NOT_ESTABLISHED. No validated fix exists. This is a pilot harness, not production code.
- 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.
## Not supported in this release
3D, dynamic/automatic AMR (tagging, time-varying hierarchies beyond `remap_child`),
moving bodies on refined levels, pressure Jets, Jet checkpoints, LES or TRT on refined
levels, partial (non-dividing) tiles.
+14 -3
View File
@@ -207,19 +207,30 @@ def surface_link_velocity(objects, cfg, flags, fluid_idx, dir, body_id, rx, ry,
jets=(), hits=None, *, _shifted_blocked=False):
"""Per-link prescribed surface velocity (n, 2) FP32 and aperture mask.
See _surface_links, which also returns the unit-amplitude jet shapes.
"""
return _surface_links(objects, cfg, flags, fluid_idx, dir, body_id, rx, ry, jets, hits,
_shifted_blocked=_shifted_blocked)[:2]
def _surface_links(objects, cfg, flags, fluid_idx, dir, body_id, rx, ry,
jets=(), hits=None, *, _shifted_blocked=False):
"""(velocity, aperture, jet shape, jet id); the last two are None without jets.
Host-only and storage-neutral: rows follow the given physical cut-link
order. Returns (None, None) when neither an object carries surface
order. Returns all None when neither an object carries surface
actuation nor a jet is given. Jets need true FP64 ``hits`` (n, 2) and may
not share a link with each other or with a legacy aperture.
Shifted storage is privately enabled only by the all-blocked stepper.
"""
enabled = [obj for obj in objects if obj.surface_actuation is not None]
if not enabled and not jets:
return None, None
return None, None, None, None
from .preprocess.flag_overlay import validate_cut_link_ownership
validate_cut_link_ownership(flags, fluid_idx, dir, cfg.nq)
surface_velocity = np.zeros((len(dir), 2), np.float32)
surface_aperture = np.zeros(len(dir), bool)
shape = jet_id = None
for obj in enabled:
obj.validate_surface_actuation(cfg)
selected = body_id == obj.obj_id
@@ -252,4 +263,4 @@ def surface_link_velocity(objects, cfg, flags, fluid_idx, dir, body_id, rx, ry,
raise ValueError("Jet overlaps a SurfaceActuation aperture")
surface_velocity[driven] = jet_link_velocity(jets, shape, jet_id)[driven]
surface_aperture |= driven
return surface_velocity, surface_aperture
return surface_velocity, surface_aperture, shape, jet_id
+12 -4
View File
@@ -54,14 +54,18 @@ def _write(path: str, content: str):
def generate_config(cfg: LBMConfig, n_objects: int = 0,
*, build_root: Optional[str] = None) -> str:
*, build_root: Optional[str] = None, block_root_bc: int = 0) -> str:
"""Snapshot kernels and render headers in a fresh private directory.
Returns its path. Pass it explicitly to compile_kernel for interleaved
builds. Legacy generate_config(); compile_kernel() remains context-local.
The caller owns cleanup_build_dir() once the loaded module is no longer
needed for recompilation/debugging. No package files are changed.
block_root_bc specializes a blocked root module; fine/dense builds use 0.
"""
if (isinstance(block_root_bc, bool) or not isinstance(block_root_bc, int)
or block_root_bc not in (0, 1, 2)):
raise ValueError("block_root_bc must be an integer mode 0, 1 or 2")
macros = cfg.to_macros()
if isinstance(n_objects, bool) or not isinstance(n_objects, int) or n_objects < 0:
raise ValueError("n_objects must be a nonnegative integer")
@@ -71,14 +75,15 @@ def generate_config(cfg: LBMConfig, n_objects: int = 0,
try:
shutil.copytree(_KERNEL_DIR, build_dir, dirs_exist_ok=True,
ignore=shutil.ignore_patterns("*.ptx", "*.cubin", "*.tmp"))
_render_config(macros, n_objects, os.path.join(build_dir, "config"))
_render_config(macros, n_objects, os.path.join(build_dir, "config"),
block_root_bc=block_root_bc)
digest = hashlib.sha256()
for path in sorted(Path(build_dir).rglob("*")):
if path.is_file():
digest.update(str(path.relative_to(build_dir)).encode())
digest.update(path.read_bytes())
metadata = {"source_sha256": digest.hexdigest(), "config": asdict(cfg),
"n_objects": n_objects}
"n_objects": n_objects, "block_root_bc": block_root_bc}
_write(os.path.join(build_dir, ".celeris-build"), json.dumps(metadata, sort_keys=True))
except BaseException:
shutil.rmtree(build_dir)
@@ -101,7 +106,7 @@ def cleanup_build_dir(build_dir: Optional[str]) -> None:
_last_ptx_path.set(None)
def _render_config(m: dict, n_objects: int, cdir: str):
def _render_config(m: dict, n_objects: int, cdir: str, *, block_root_bc: int = 0):
_write(os.path.join(cdir, "config_grid.h"), f"""{_HEADER}\
// Layer: Global / Grid
#ifndef CELERIS_CONFIG_GRID_H
@@ -165,6 +170,9 @@ def _render_config(m: dict, n_objects: int, cdir: str):
#define STORE_PRECISION {m['STORE_PRECISION']}
#define USE_DDF_SHIFTING {m['USE_DDF_SHIFTING']}
// Internal blocked root mode: 0=periodic/fine, 1=slot/reservoir, 2=channel.
#define BLOCK_ROOT_BC {block_root_bc}
#define USE_LES {m['USE_LES']}
#define LES_CS {m['LES_CS']}
#define LES_CLOSED_FORM {m['LES_CLOSED_FORM']}
+2 -4
View File
@@ -94,12 +94,10 @@ class CurvedLinkSoA:
"""
self.surface_velocity = self.aperture = None
self.jets, self.jet_shape, self.jet_id = tuple(jets), None, None
self.surface_velocity, self.aperture = surface_link_velocity(
from ..body.jets import _surface_links
self.surface_velocity, self.aperture, self.jet_shape, self.jet_id = _surface_links(
objects, cfg, flags, self.fluid_idx, self.dir, self.body_id, self.rx, self.ry,
jets=self.jets, hits=hits)
if self.jets:
from ..body.jets import jet_link_shapes
self.jet_shape, self.jet_id = jet_link_shapes(self.jets, hits)
def set_jet_amplitude(self, index, amplitude, stream=None):
"""Rewrite jet rows from cached shapes; re-upload only the velocity column.
+22 -179
View File
@@ -1,5 +1,6 @@
// CelerisLab/lbm/kernels/blocks.cu
// One periodic level of fixed-size blocks (lbm/blocks.py layout). Separate
// One fixed-size blocked level (lbm/blocks.py layout), specialized at compile
// time for periodic/fine, slot/reservoir or channel stepping. Separate
// compilation entry: NX, NY are the block window including halo, so index_f
// and linear_index address one block; blockIdx.z selects the block, whose
// populations start at b*NQ*TOTAL_CELLS (block-contiguous, Q-major inside).
@@ -12,6 +13,7 @@
#include "core/params.cuh"
#include "core/flags.cuh"
#include "core/obs.cuh"
#include "core/cell_kind.cuh"
#include "operators/macro.cuh"
#include "operators/collision_srt.cuh"
#if COLLISION_MODEL == 2
@@ -19,188 +21,29 @@
#endif
#include "streaming/pull_double_buffer.cuh"
#include "boundary/curved_boundary.cuh"
#include "boundary/planar_velocity.cuh"
#if DIM != 2 || NQ != 9 || (COLLISION_MODEL != 0 && COLLISION_MODEL != 2) || STORE_PRECISION != 0 || USE_LES || STREAMING_MODEL != 0
#error "Blocks require FP32 D2Q9 SRT or MRT double buffers without LES"
#endif
// Fixed-patch cell kinds, matching refinement.py ACTIVE/INNER_GHOST/OUTER_GHOST.
constexpr unsigned char KIND_ACTIVE = 0u, KIND_INNER = 1u, KIND_OUTER = 2u;
// Zero-initialized for ordinary blocks/fine levels. Set only for parent walls;
// NX/NY remain local storage strides, never physical wall coordinates.
struct BlockWall { unsigned int nx, ny, slot_width; float amplitude; };
extern "C" { __constant__ BlockWall d_block_wall; }
// Halo refresh: block cell ids b*TOTAL_CELLS+c, sources are always core cells.
extern "C" __global__ void BlockHaloCopy(
fpxx* f, const uint64_t* dst, const uint64_t* src, unsigned int count)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t d = dst[row], s = src[row];
fpxx* to = f + (d / TOTAL_CELLS) * NQ * TOTAL_CELLS;
const fpxx* from = f + (s / TOTAL_CELLS) * NQ * TOTAL_CELLS;
for (unsigned int q = 0; q < NQ; ++q)
to[index_f(d % TOTAL_CELLS, q)] = from[index_f(s % TOTAL_CELLS, q)];
}
// Fixed-patch explosion into blocked fine storage: target = owner block cell
// id of a patch-ghost cell; donors are parent block cell ids. Dense parent
// storage is the one-window case with coarse_cells = full domain size.
extern "C" __global__ void BlockExplode(
const fpxx* coarse, fpxx* fine, const uint64_t* target, const uint64_t* donor,
unsigned int count, unsigned long coarse_cells)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t t = target[row];
fpxx* to = fine + (t / TOTAL_CELLS) * NQ * TOTAL_CELLS;
const uint64_t d = donor[row];
const fpxx* from = coarse + (d / coarse_cells) * NQ * coarse_cells;
for (unsigned int q = 0; q < NQ; ++q)
to[index_f(t % TOTAL_CELLS, q)] = from[q * coarse_cells + d % coarse_cells];
}
// Experimental Schornbaum--Ruede (arXiv:1508.07982, Algorithm 2)
// and Chen et al. (doi:10.1016/j.physa.2005.09.036) population interpolation.
// Independent implementation of the mathematical stencil; no upstream source copied.
// Stencil cells are owner cores; a disabled axis repeats its center twice.
extern "C" __global__ void BlockExplodeProjected(
const fpxx* coarse, fpxx* fine, const uint64_t* target, const uint64_t* donor,
unsigned int count, unsigned long coarse_cells,
const uint64_t* neighbors, const float* offsets)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t t = target[row], d = donor[row];
fpxx* to = fine + (t / TOTAL_CELLS) * NQ * TOTAL_CELLS;
for (unsigned int q = 0; q < NQ; ++q) {
float v[4];
for (unsigned int a = 0; a < 4; ++a) {
const uint64_t n = neighbors[4ul * row + a];
v[a] = coarse[(n / coarse_cells) * NQ * coarse_cells
+ q * coarse_cells + n % coarse_cells];
}
float gx = .5f * (v[1] - v[0]), gy = .5f * (v[3] - v[2]);
const float ex = d_cx[q], ey = d_cy[q];
if (q != 0) {
const float along = (ex * gx + ey * gy) / (ex * ex + ey * ey);
gx -= ex * along;
gy -= ey * along;
}
to[index_f(t % TOTAL_CELLS, q)] =
coarse[(d / coarse_cells) * NQ * coarse_cells + q * coarse_cells + d % coarse_cells]
+ offsets[2ul * row] * gx + offsets[2ul * row + 1] * gy;
}
}
// Targets are receiver-local SOURCE slots for the immediately following pull.
// Keep child summation order identical to PatchCoarseStep; no collision here.
extern "C" __global__ void BlockCoalesce(
const fpxx* fine, fpxx* coarse, const uint64_t* target, const uint64_t* children,
unsigned int count)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t* s = children + 4ul * row;
coarse[target[row]] = 0.25f * (fine[s[0]] + fine[s[1]] + fine[s[2]] + fine[s[3]]);
}
// Static cut links, receiver-local staging: cell is the receiver's local core
// cell, so the existing BFL writes the solid-source window image that this
// block's next pull reads (no halo copy between this kernel and BlockStep).
// surface_velocity (nullable, 2 floats per link): prescribed fluid velocity on
// fixed geometry, same convention as CurvedBoundaryKernel. Optional raw rows
// (jx,jy,torque,dm) are unique physical links, never body-observation slots.
extern "C" __global__ void BlockCurvedBoundary(
fpxx* f, const uint32_t* block, const uint32_t* cell, const unsigned char* dir,
const float* q, const unsigned char* fallback, unsigned int count,
const float* surface_velocity, const float* rx, const float* ry, float* raw)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
float Uw = 0.0f, Vw = 0.0f;
if (surface_velocity != nullptr) {
const float u = surface_velocity[2 * row], v = surface_velocity[2 * row + 1];
if (u != 0.0f || v != 0.0f) { Uw = u; Vw = v; }
}
apply_bouzidi_link(dir[row], q[row], cell[row], f + (unsigned long)block[row] * NQ * TOTAL_CELLS,
Uw, Vw, rx ? rx[row] : 0.0f, ry ? ry[row] : 0.0f,
fallback[row], nullptr, -1, raw ? raw + 4ul * row : nullptr,
raw ? raw + 4ul * row + 3 : nullptr);
}
// Same-stream BFL -> this row-owned accumulation -> pull step. No atomics,
// no cross-row reduction and no history allocation proportional to substeps.
extern "C" __global__ void BlockAccumulateExchange(
const float* raw, double* integral, unsigned int count, double h, double dt)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const unsigned long k = 4ul * row;
const double impulse = h * dt;
integral[k] += impulse * (double)raw[k];
integral[k + 1] += impulse * (double)raw[k + 1];
integral[k + 2] += h * impulse * (double)raw[k + 2];
integral[k + 3] += h * h * (double)raw[k + 3];
}
// Pull/collide/store core cells only; halo cells of out are left untouched.
// flags (block cell ids, nullable): solid cores store NaN and never collide.
// kind (block cell ids, nullable = all ACTIVE): OUTER skips the pull and
// stores NaN, INNER pulls without collision, ACTIVE pulls and collides.
extern "C" __global__ void BlockStep(const fpxx* in, fpxx* out, unsigned int halo,
const uint16_t* flags, const unsigned char* kind) {
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x < halo || y < halo || x >= NX - halo || y >= NY - halo) return;
const unsigned long base = (unsigned long)blockIdx.z * NQ * TOTAL_CELLS;
const unsigned long cell = (unsigned long)blockIdx.z * TOTAL_CELLS + k;
const unsigned char cell_kind = kind != nullptr ? kind[cell] : KIND_ACTIVE;
float f[NQ];
if (cell_kind == KIND_OUTER || (flags != nullptr && is_solid(flags[cell]))) {
for (int q = 0; q < NQ; ++q) f[q] = NAN;
} else {
unsigned long j[NQ];
compute_neighbors(x, y, j); // Core neighbors stay inside the window: no wrap.
stream_pull_load(k, f, in + base, j);
if (cell_kind == KIND_ACTIVE) {
if (d_block_wall.nx != 0) {
const unsigned int bx = NX - 2 * halo, by = NY - 2 * halo;
const unsigned int nbx = d_block_wall.nx / bx;
const unsigned int gx = (blockIdx.z % nbx) * bx + x - halo;
const unsigned int gy = (blockIdx.z / nbx) * by + y - halo;
if (gy == 1)
reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(
gx, d_block_wall.nx, d_block_wall.slot_width, d_block_wall.amplitude));
else if (gy == d_block_wall.ny - 2) reconstruct_top_reservoir(f, RHO);
}
#if COLLISION_MODEL == 2
// Uniform MRT sequence (collide_dispatch, zero forcing); SRT path below unchanged.
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
#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];
const float jx = f[1]-f[2]+f[5]-f[6]+f[7]-f[8];
const float jy = f[3]-f[4]+f[5]-f[6]-f[7]+f[8];
compute_feq_from_momentum(drho, jx, jy, ux, uy, feq);
#else
compute_feq(rho, ux, uy, feq);
#ifndef BLOCK_ROOT_BC
#error "Blocks require a generated BLOCK_ROOT_BC mode"
#elif BLOCK_ROOT_BC < 0 || BLOCK_ROOT_BC > 2
#error "BLOCK_ROOT_BC must be 0, 1 or 2"
#endif
collide_srt_no_force(f, feq, d_params.omega);
#endif
}
}
stream_pull_store(k, f, out + base);
}
// Root channel closure; fine and periodic levels keep BlockStep.
#if BLOCK_ROOT_BC == 1
#include "boundary/block_wall.cuh"
#endif
#if BLOCK_ROOT_BC != 2
#include "operators/block_collision.cuh"
#endif
// Kernels, in definition order: transfers, boundary staging, then the one full
// step selected for this module (BlockStep, or BlockChannelStep for a channel root).
#include "step/block_transfer.cu"
#include "step/block_boundary.cu"
#if BLOCK_ROOT_BC != 2
#include "step/block_step.cu"
#else
#include "step/block_channel.cu"
#endif
@@ -0,0 +1,20 @@
// Physical parent walls; NX/NY remain local storage strides.
#ifndef CELERIS_BOUNDARY_BLOCK_WALL_CUH
#define CELERIS_BOUNDARY_BLOCK_WALL_CUH
#include "planar_velocity.cuh"
struct BlockWall { unsigned int nx, ny, slot_width; float amplitude; };
extern "C" { __constant__ BlockWall d_block_wall; }
__device__ __forceinline__ void block_wall(float* f, unsigned int x, unsigned int y,
unsigned int halo) {
unsigned int gx, gy;
block_global_xy(x, y, halo, d_block_wall.nx, gx, gy);
if (gy == 1)
reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(
gx, d_block_wall.nx, d_block_wall.slot_width, d_block_wall.amplitude));
else if (gy == d_block_wall.ny - 2) reconstruct_top_reservoir(f, RHO);
}
#endif
@@ -15,10 +15,13 @@
// D2Q9 y-wall half-way bounce-back for pull double-buffer streaming.
// ---------------------------------------------------------------------------
#if LATTICE_MODEL == LATTICE_D2Q9
__device__ inline void apply_wall_bb_d2q9_y_pull(unsigned int y,
float* __restrict__ f,
const fpxx* __restrict__ fi_in,
unsigned long k)
// Rows 1 and top_row (= physical ny - 2) next to walls at rows 0 and ny - 1.
// y is the physical row; k addresses the node in fi_in (a block window may
// have NX/NY strides different from the physical extent).
__device__ inline void apply_wall_bb_d2q9_rows(unsigned int y, unsigned int top_row,
float* __restrict__ f,
const fpxx* __restrict__ fi_in,
unsigned long k)
{
if (y == 1) {
// Directions sourced from y=0 wall in pull step: +y, +x+y, -x+y.
@@ -26,14 +29,22 @@ __device__ inline void apply_wall_bb_d2q9_y_pull(unsigned int y,
f[5] = load_ddf(fi_in, index_f(k, 6u));
f[8] = load_ddf(fi_in, index_f(k, 7u));
}
else if (y == (unsigned int)(NY - 2)) {
// Directions sourced from y=NY-1 wall in pull step: -y, -x-y, +x-y.
else if (y == top_row) {
// Directions sourced from the y = ny - 1 wall in pull step: -y, -x-y, +x-y.
f[4] = load_ddf(fi_in, index_f(k, 3u));
f[6] = load_ddf(fi_in, index_f(k, 5u));
f[7] = load_ddf(fi_in, index_f(k, 8u));
}
}
__device__ inline void apply_wall_bb_d2q9_y_pull(unsigned int y,
float* __restrict__ f,
const fpxx* __restrict__ fi_in,
unsigned long k)
{
apply_wall_bb_d2q9_rows(y, (unsigned int)(NY - 2), f, fi_in, k);
}
// ---------------------------------------------------------------------------
// D2Q9 y-wall free-slip (specular reflection) for pull double-buffer streaming.
//
@@ -8,6 +8,9 @@
#define STORE_PRECISION 0
#define USE_DDF_SHIFTING 0
// Internal blocked root mode: 0=periodic/fine, 1=slot/reservoir, 2=channel.
#define BLOCK_ROOT_BC 0
#define USE_LES 0
#define LES_CS 0.160000f
#define LES_CLOSED_FORM 1
@@ -0,0 +1,16 @@
// CelerisLab – core/cell_kind.cuh
// Fixed-patch cell kinds stored per block cell, matching lbm/refinement.py
// ACTIVE/INNER_GHOST/OUTER_GHOST: OUTER is invalid (stores NaN, no pull),
// INNER streams without collision, ACTIVE pulls and collides.
#ifndef CELERIS_CORE_CELL_KIND_CUH
#define CELERIS_CORE_CELL_KIND_CUH
constexpr unsigned char KIND_ACTIVE = 0u, KIND_INNER = 1u, KIND_OUTER = 2u;
// NaN populations mark a cell without valid state (OUTER ghost or solid core).
__device__ __forceinline__ void fill_invalid(float* f)
{
for (int q = 0; q < NQ; ++q) f[q] = NAN;
}
#endif
@@ -53,6 +53,18 @@ __device__ __forceinline__ void index_from_thread(unsigned int& x,
k = (unsigned long)y * NX + x;
}
// Global cell of window cell (x, y) in block blockIdx.z of a row-major,
// full-domain block layout domain_nx cells wide. NX/NY are the local strides.
__device__ __forceinline__ void block_global_xy(unsigned int x, unsigned int y,
unsigned int halo, unsigned int domain_nx,
unsigned int& gx, unsigned int& gy)
{
const unsigned int bx = NX - 2 * halo, by = NY - 2 * halo;
const unsigned int nbx = domain_nx / bx;
gx = (blockIdx.z % nbx) * bx + x - halo;
gy = (blockIdx.z / nbx) * by + y - halo;
}
// ---------------------------------------------------------------------------
// Neighbor list (periodic BCs via modular arithmetic)
// ---------------------------------------------------------------------------
@@ -0,0 +1,31 @@
// Collision used by periodic/fine and slot/reservoir BlockStep.
#ifndef CELERIS_OPERATORS_BLOCK_COLLISION_CUH
#define CELERIS_OPERATORS_BLOCK_COLLISION_CUH
#include "../core/params.cuh"
#include "macro.cuh"
#include "collision_srt.cuh"
#if COLLISION_MODEL == 2
#include "collision_mrt.cuh"
#endif
__device__ __forceinline__ void block_collide(float* f) {
#if COLLISION_MODEL == 2
collide_mrt_unforced(f); // SRT path below unchanged.
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
#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];
const float jx = f[1]-f[2]+f[5]-f[6]+f[7]-f[8];
const float jy = f[3]-f[4]+f[5]-f[6]-f[7]+f[8];
compute_feq_from_momentum(drho, jx, jy, ux, uy, feq);
#else
compute_feq(rho, ux, uy, feq);
#endif
collide_srt_no_force(f, feq, d_params.omega);
#endif
}
#endif
@@ -0,0 +1,14 @@
// CelerisLab – operators/bounce_back_swap.cuh
// Full-way bounce-back at a wall node: swap every opposite pair (paired
// order q, q+1). Needs only NQ; shared by the uniform and blocked steps.
#ifndef CELERIS_OPERATORS_BOUNCE_BACK_SWAP_CUH
#define CELERIS_OPERATORS_BOUNCE_BACK_SWAP_CUH
__device__ __forceinline__ void bounce_back_swap(float* f) {
#pragma unroll
for (int i = 1; i < NQ; i += 2) {
float t = f[i]; f[i] = f[i+1]; f[i+1] = t;
}
}
#endif
@@ -114,6 +114,16 @@ __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,
@@ -5,12 +5,7 @@
#ifndef CELERIS_OPERATORS_HELPERS_CUH
#define CELERIS_OPERATORS_HELPERS_CUH
__device__ __forceinline__ void bounce_back_swap(float* f) {
#pragma unroll
for (int i = 1; i < NQ; i += 2) {
float t = f[i]; f[i] = f[i+1]; f[i+1] = t;
}
}
#include "bounce_back_swap.cuh"
__device__ __forceinline__ void collide_dispatch(
float* f, float rho_n,
+3 -93
View File
@@ -9,6 +9,7 @@
#include "core/params.cuh"
#include "core/flags.cuh"
#include "core/obs.cuh"
#include "core/cell_kind.cuh"
#include "operators/macro.cuh"
#include "operators/collision_srt.cuh"
#if COLLISION_MODEL == 2
@@ -35,96 +36,5 @@ __device__ unsigned int incoming_counts[NQ * TOTAL_CELLS];
}
#endif
__device__ __forceinline__ void patch_collide(float* f, unsigned long k) {
#if COLLISION_MODEL == 2
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
compute_feq(rho, ux, uy, feq);
collide_srt_no_force(f, feq, d_params.omega);
#endif
#ifdef CELERIS_REFINEMENT_AUDIT
++collision_counts[k];
#endif
}
extern "C" __global__ void PatchExplode(
const float* coarse, float* fine, const uint32_t* targets,
const uint32_t* donors, unsigned int count, unsigned long coarse_cells)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
for (unsigned int q = 0; q < NQ; ++q)
fine[index_f(targets[row], q)] = coarse[q * coarse_cells + donors[row]];
#ifdef CELERIS_REFINEMENT_AUDIT
++explosion_counts[targets[row]];
#endif
}
extern "C" __global__ void PatchFineStep(const float* in, float* out, const uint16_t* flags) {
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x >= NX || y >= NY) return;
float f[NQ];
// Outer halo is invalid after EACH hop, including when reusing the buffer
// that held explosion data. Inner halo only streams; it never collides.
if (x < 2 || y < 2 || x >= NX - 2 || y >= NY - 2
|| (flags != nullptr && is_solid(flags[k]))) {
for (int q = 0; q < NQ; ++q) f[q] = NAN;
} else {
unsigned long j[NQ];
compute_neighbors(x, y, j); // These interior targets never wrap.
stream_pull_load(k, f, in, j);
if (x >= 4 && y >= 4 && x < NX - 4 && y < NY - 4)
patch_collide(f, k);
}
stream_pull_store(k, f, out);
}
extern "C" __global__ void PatchCoarseStep(
const float* in, float* out, const float* fine,
const uint64_t* targets, const uint64_t* children, unsigned int count,
unsigned int x0, unsigned int y0, unsigned int x1, unsigned int y1,
unsigned int wall_aware)
{
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x >= NX || y >= NY) return;
// Both buffers were poisoned at initialization; covered cells stay invalid.
if (x >= x0 && x < x1 && y >= y0 && y < y1) return;
if (wall_aware && (y == 0 || y == NY - 1)) return;
unsigned long j[NQ];
compute_neighbors(x, y, j);
float f[NQ];
stream_pull_load(k, f, in, j);
if (x >= x0 - 1 && x <= x1 && y >= y0 - 1 && y <= y1) {
for (unsigned int q = 1; q < NQ; ++q) {
uint64_t slot = index_f(k, q);
unsigned int lo = 0, hi = count;
while (lo < hi) {
unsigned int mid = lo + (hi - lo) / 2;
if (targets[mid] < slot) lo = mid + 1;
else hi = mid;
}
if (lo < count && targets[lo] == slot) {
const uint64_t* s = children + 4ul * lo;
f[q] = 0.25f * (fine[s[0]] + fine[s[1]] + fine[s[2]] + fine[s[3]]);
#ifdef CELERIS_REFINEMENT_AUDIT
++incoming_counts[slot];
#endif
}
}
}
// Includes the eight grazing-corner links with UNCOVERED coarse donors.
// All replacement arrivals must be in place before computing rho/u.
if (wall_aware) {
if (y == 1) reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(x));
else if (y == NY - 2) reconstruct_top_reservoir(f, RHO);
}
patch_collide(f, k);
stream_pull_store(k, f, out);
}
// Kernels: patch collision helper, PatchExplode, PatchFineStep, PatchCoarseStep.
#include "step/patch_step.cu"
@@ -0,0 +1,42 @@
// CelerisLab/lbm/kernels/step/block_boundary.cu
// Receiver-local curved-boundary staging and its raw-exchange accumulation.
// Part of the blocks.cu compilation entry; relies on its includes and guards.
// Static cut links, receiver-local staging: cell is the receiver's local core
// cell, so the existing BFL writes the solid-source window image that this
// block's next pull reads (no halo copy between this kernel and BlockStep).
// surface_velocity (nullable, 2 floats per link): prescribed fluid velocity on
// fixed geometry, same convention as CurvedBoundaryKernel. Optional raw rows
// (jx,jy,torque,dm) are unique physical links, never body-observation slots.
extern "C" __global__ void BlockCurvedBoundary(
fpxx* f, const uint32_t* block, const uint32_t* cell, const unsigned char* dir,
const float* q, const unsigned char* fallback, unsigned int count,
const float* surface_velocity, const float* rx, const float* ry, float* raw)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
float Uw = 0.0f, Vw = 0.0f;
if (surface_velocity != nullptr) {
const float u = surface_velocity[2 * row], v = surface_velocity[2 * row + 1];
if (u != 0.0f || v != 0.0f) { Uw = u; Vw = v; }
}
apply_bouzidi_link(dir[row], q[row], cell[row], f + (unsigned long)block[row] * NQ * TOTAL_CELLS,
Uw, Vw, rx ? rx[row] : 0.0f, ry ? ry[row] : 0.0f,
fallback[row], nullptr, -1, raw ? raw + 4ul * row : nullptr,
raw ? raw + 4ul * row + 3 : nullptr);
}
// Same-stream BFL -> this row-owned accumulation -> pull step. No atomics,
// no cross-row reduction and no history allocation proportional to substeps.
extern "C" __global__ void BlockAccumulateExchange(
const float* raw, double* integral, unsigned int count, double h, double dt)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const unsigned long k = 4ul * row;
const double impulse = h * dt;
integral[k] += impulse * (double)raw[k];
integral[k + 1] += impulse * (double)raw[k + 1];
integral[k + 2] += h * impulse * (double)raw[k + 2];
integral[k + 3] += h * h * (double)raw[k + 3];
}
@@ -1,5 +1,12 @@
#if BLOCK_ROOT_BC != 2
#error "BlockChannelStep requires BLOCK_ROOT_BC == 2"
#endif
// Root-only blocked counterpart of OneStep's D2Q9/SRT channel closure.
// NX/NY are LOCAL strides; all physical boundary decisions use global coordinates.
#include "../core/cell_kind.cuh"
#include "../operators/bounce_back_swap.cuh"
#include "../boundary/bounce_back.cuh"
#include "../boundary/inlet/common.cuh"
#include "../boundary/inlet/zou_he_local.cuh"
#include "../boundary/outlet/pressure_neq.cuh"
@@ -15,15 +22,13 @@ extern "C" __global__ void BlockChannelStep(
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x < halo || y < halo || x >= NX - halo || y >= NY - halo) return;
const unsigned int bx = NX - 2 * halo, by = NY - 2 * halo;
const unsigned int nbx = d_block_channel.nx / bx;
const unsigned int gx = (blockIdx.z % nbx) * bx + x - halo;
const unsigned int gy = (blockIdx.z / nbx) * by + y - halo;
unsigned int gx, gy;
block_global_xy(x, y, halo, d_block_channel.nx, gx, gy);
const unsigned long base = (unsigned long)blockIdx.z * NQ * TOTAL_CELLS;
const unsigned long cell = (unsigned long)blockIdx.z * TOTAL_CELLS + k;
float f[NQ];
if (kind != nullptr && kind[cell] == KIND_OUTER) {
for (int q = 0; q < NQ; ++q) f[q] = NAN;
fill_invalid(f);
} else {
unsigned long j[NQ];
compute_neighbors(x, y, j);
@@ -31,35 +36,21 @@ extern "C" __global__ void BlockChannelStep(
const bool wall = gy == 0 || gy == d_block_channel.ny - 1;
const bool outlet = gx == d_block_channel.nx - 1;
if (wall) {
// Wall owns all four corners. Ghosts pull and swap, never collide.
for (int q = 1; q < NQ; q += 2) {
const float t = f[q]; f[q] = f[q + 1]; f[q + 1] = t;
}
bounce_back_swap(f); // Wall owns all four corners; ghosts never collide.
} else if (gx == 0) {
apply_zou_he_left_velocity_inlet_d2q9(
f, inlet_target_u((float)gy, d_block_channel.ny, d_block_channel.velocity), 0.0f);
} else if (outlet) {
// Immutable input donor, receiver-local (also valid across a block seam).
float donor[NQ];
for (int q = 0; q < NQ; ++q)
donor[q] = load_ddf(in + base, index_f(k - 1, q));
load_node(in + base, k - 1, donor);
apply_pressure_outlet_d2q9(f, donor);
} else if (gy == 1) {
f[3] = load_ddf(in + base, index_f(k, 4));
f[5] = load_ddf(in + base, index_f(k, 6));
f[8] = load_ddf(in + base, index_f(k, 7));
} else if (gy == d_block_channel.ny - 2) {
f[4] = load_ddf(in + base, index_f(k, 3));
f[6] = load_ddf(in + base, index_f(k, 5));
f[7] = load_ddf(in + base, index_f(k, 8));
} else {
apply_wall_bb_d2q9_rows(gy, d_block_channel.ny - 2, f, in + base, k);
}
if (!wall && !outlet) {
#if COLLISION_MODEL == 2
// Match uniform MRT arithmetic (collide_dispatch, zero forcing).
float rho, ux, uy;
compute_rho_u(f, rho, ux, uy);
collide_mrt_no_force(f, rho, ux, uy,
fminf(OMEGA_COLLISION_MAX, fmaxf(OMEGA_COLLISION_MIN, d_params.omega)));
collide_mrt_unforced(f);
#else
// Match uniform SRT arithmetic, including its shifted equilibrium path.
float rho, ux, uy, feq[NQ], forcing[NQ] = {};
@@ -0,0 +1,35 @@
// CelerisLab/lbm/kernels/step/block_step.cu
// Generic blocked pull/collide/store for periodic/fine and slot/reservoir levels.
// Part of the blocks.cu compilation entry; relies on its includes and guards.
#if BLOCK_ROOT_BC == 2
#error "BlockStep is not built for the channel root (BLOCK_ROOT_BC == 2)"
#endif
// Pull/collide/store core cells only; halo cells of out are left untouched.
// flags (block cell ids, nullable): solid cores store NaN and never collide.
// kind (block cell ids, nullable = all ACTIVE): OUTER skips the pull and
// stores NaN, INNER pulls without collision, ACTIVE pulls and collides.
extern "C" __global__ void BlockStep(const fpxx* in, fpxx* out, unsigned int halo,
const uint16_t* flags, const unsigned char* kind) {
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x < halo || y < halo || x >= NX - halo || y >= NY - halo) return;
const unsigned long base = (unsigned long)blockIdx.z * NQ * TOTAL_CELLS;
const unsigned long cell = (unsigned long)blockIdx.z * TOTAL_CELLS + k;
const unsigned char cell_kind = kind != nullptr ? kind[cell] : KIND_ACTIVE;
float f[NQ];
if (cell_kind == KIND_OUTER || (flags != nullptr && is_solid(flags[cell]))) {
fill_invalid(f);
} else {
unsigned long j[NQ];
compute_neighbors(x, y, j); // Core neighbors stay inside the window: no wrap.
stream_pull_load(k, f, in + base, j);
if (cell_kind == KIND_ACTIVE) {
#if BLOCK_ROOT_BC == 1
block_wall(f, x, y, halo);
#endif
block_collide(f);
}
}
stream_pull_store(k, f, out + base);
}
@@ -0,0 +1,78 @@
// CelerisLab/lbm/kernels/step/block_transfer.cu
// Same-level halo refresh and 2:1 explosion/coalescence between blocked levels.
// Part of the blocks.cu compilation entry; relies on its includes and guards.
// Halo refresh: block cell ids b*TOTAL_CELLS+c, sources are always core cells.
extern "C" __global__ void BlockHaloCopy(
fpxx* f, const uint64_t* dst, const uint64_t* src, unsigned int count)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t d = dst[row], s = src[row];
fpxx* to = f + (d / TOTAL_CELLS) * NQ * TOTAL_CELLS;
const fpxx* from = f + (s / TOTAL_CELLS) * NQ * TOTAL_CELLS;
for (unsigned int q = 0; q < NQ; ++q)
to[index_f(d % TOTAL_CELLS, q)] = from[index_f(s % TOTAL_CELLS, q)];
}
// Fixed-patch explosion into blocked fine storage: target = owner block cell
// id of a patch-ghost cell; donors are parent block cell ids. Dense parent
// storage is the one-window case with coarse_cells = full domain size.
extern "C" __global__ void BlockExplode(
const fpxx* coarse, fpxx* fine, const uint64_t* target, const uint64_t* donor,
unsigned int count, unsigned long coarse_cells)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t t = target[row];
fpxx* to = fine + (t / TOTAL_CELLS) * NQ * TOTAL_CELLS;
const uint64_t d = donor[row];
const fpxx* from = coarse + (d / coarse_cells) * NQ * coarse_cells;
for (unsigned int q = 0; q < NQ; ++q)
to[index_f(t % TOTAL_CELLS, q)] = from[q * coarse_cells + d % coarse_cells];
}
// Experimental Schornbaum--Ruede (arXiv:1508.07982, Algorithm 2)
// and Chen et al. (doi:10.1016/j.physa.2005.09.036) population interpolation.
// Independent implementation of the mathematical stencil; no upstream source copied.
// Stencil cells are owner cores; a disabled axis repeats its center twice.
extern "C" __global__ void BlockExplodeProjected(
const fpxx* coarse, fpxx* fine, const uint64_t* target, const uint64_t* donor,
unsigned int count, unsigned long coarse_cells,
const uint64_t* neighbors, const float* offsets)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t t = target[row], d = donor[row];
fpxx* to = fine + (t / TOTAL_CELLS) * NQ * TOTAL_CELLS;
for (unsigned int q = 0; q < NQ; ++q) {
float v[4];
for (unsigned int a = 0; a < 4; ++a) {
const uint64_t n = neighbors[4ul * row + a];
v[a] = coarse[(n / coarse_cells) * NQ * coarse_cells
+ q * coarse_cells + n % coarse_cells];
}
float gx = .5f * (v[1] - v[0]), gy = .5f * (v[3] - v[2]);
const float ex = d_cx[q], ey = d_cy[q];
if (q != 0) {
const float along = (ex * gx + ey * gy) / (ex * ex + ey * ey);
gx -= ex * along;
gy -= ey * along;
}
to[index_f(t % TOTAL_CELLS, q)] =
coarse[(d / coarse_cells) * NQ * coarse_cells + q * coarse_cells + d % coarse_cells]
+ offsets[2ul * row] * gx + offsets[2ul * row + 1] * gy;
}
}
// Targets are receiver-local SOURCE slots for the immediately following pull.
// Keep child summation order identical to PatchCoarseStep; no collision here.
extern "C" __global__ void BlockCoalesce(
const fpxx* fine, fpxx* coarse, const uint64_t* target, const uint64_t* children,
unsigned int count)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
const uint64_t* s = children + 4ul * row;
coarse[target[row]] = 0.25f * (fine[s[0]] + fine[s[1]] + fine[s[2]] + fine[s[3]]);
}
@@ -0,0 +1,95 @@
// CelerisLab/lbm/kernels/step/patch_step.cu
// Dense fixed-patch collision, explosion and the fine/coarse steps.
// Part of the refinement.cu compilation entry; relies on its includes, guards
// and the optional CELERIS_REFINEMENT_AUDIT counters declared there.
__device__ __forceinline__ void patch_collide(float* f, unsigned long k) {
#if COLLISION_MODEL == 2
collide_mrt_unforced(f);
#else
float rho, ux, uy, feq[NQ];
compute_rho_u(f, rho, ux, uy);
compute_feq(rho, ux, uy, feq);
collide_srt_no_force(f, feq, d_params.omega);
#endif
#ifdef CELERIS_REFINEMENT_AUDIT
++collision_counts[k];
#endif
}
extern "C" __global__ void PatchExplode(
const float* coarse, float* fine, const uint32_t* targets,
const uint32_t* donors, unsigned int count, unsigned long coarse_cells)
{
unsigned int row = blockIdx.x * blockDim.x + threadIdx.x;
if (row >= count) return;
for (unsigned int q = 0; q < NQ; ++q)
fine[index_f(targets[row], q)] = coarse[q * coarse_cells + donors[row]];
#ifdef CELERIS_REFINEMENT_AUDIT
++explosion_counts[targets[row]];
#endif
}
extern "C" __global__ void PatchFineStep(const float* in, float* out, const uint16_t* flags) {
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x >= NX || y >= NY) return;
float f[NQ];
// Outer halo is invalid after EACH hop, including when reusing the buffer
// that held explosion data. Inner halo only streams; it never collides.
if (x < 2 || y < 2 || x >= NX - 2 || y >= NY - 2
|| (flags != nullptr && is_solid(flags[k]))) {
fill_invalid(f);
} else {
unsigned long j[NQ];
compute_neighbors(x, y, j); // These interior targets never wrap.
stream_pull_load(k, f, in, j);
if (x >= 4 && y >= 4 && x < NX - 4 && y < NY - 4)
patch_collide(f, k);
}
stream_pull_store(k, f, out);
}
extern "C" __global__ void PatchCoarseStep(
const float* in, float* out, const float* fine,
const uint64_t* targets, const uint64_t* children, unsigned int count,
unsigned int x0, unsigned int y0, unsigned int x1, unsigned int y1,
unsigned int wall_aware)
{
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x >= NX || y >= NY) return;
// Both buffers were poisoned at initialization; covered cells stay invalid.
if (x >= x0 && x < x1 && y >= y0 && y < y1) return;
if (wall_aware && (y == 0 || y == NY - 1)) return;
unsigned long j[NQ];
compute_neighbors(x, y, j);
float f[NQ];
stream_pull_load(k, f, in, j);
if (x >= x0 - 1 && x <= x1 && y >= y0 - 1 && y <= y1) {
for (unsigned int q = 1; q < NQ; ++q) {
uint64_t slot = index_f(k, q);
unsigned int lo = 0, hi = count;
while (lo < hi) {
unsigned int mid = lo + (hi - lo) / 2;
if (targets[mid] < slot) lo = mid + 1;
else hi = mid;
}
if (lo < count && targets[lo] == slot) {
const uint64_t* s = children + 4ul * lo;
f[q] = 0.25f * (fine[s[0]] + fine[s[1]] + fine[s[2]] + fine[s[3]]);
#ifdef CELERIS_REFINEMENT_AUDIT
++incoming_counts[slot];
#endif
}
}
}
// Includes the eight grazing-corner links with UNCOVERED coarse donors.
// All replacement arrivals must be in place before computing rho/u.
if (wall_aware) {
if (y == 1) reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(x));
else if (y == NY - 2) reconstruct_top_reservoir(f, RHO);
}
patch_collide(f, k);
stream_pull_store(k, f, out);
}
@@ -21,6 +21,14 @@
// i.e. the population traveling in direction i at node n
// arrived from the node at n - c_i = neighbor in direction opp(i)
// ---------------------------------------------------------------------------
// load_node: all NQ populations of node n as stored (no streaming).
__device__ __forceinline__ void load_node(const fpxx* __restrict__ fi, unsigned long n,
float* __restrict__ f)
{
for (int q = 0; q < NQ; ++q)
f[q] = load_ddf(fi, index_f(n, q));
}
__device__ __forceinline__ void stream_pull_load(unsigned long n,
float* __restrict__ f,
const fpxx* __restrict__ fi,
+6 -3
View File
@@ -362,7 +362,10 @@ class FixedPatchStepper:
raise RuntimeError("level omega/tau mismatch")
host = self._boundary_hosts[index]
n_objects = len(host.objects) if name == "PatchFineStep" and host is not None else 0
build = compiler.generate_config(config, n_objects=n_objects)
block_root_bc = 0
if index == 0 and layout is not None:
block_root_bc = 2 if self._channel else 1 if self._wall_aware else 0
build = compiler.generate_config(config, n_objects=n_objects, block_root_bc=block_root_bc)
try:
ptx = compiler.compile_kernel(self.ctx.sm_arch, source=source, build_dir=build)
module = compiler.load_module(ptx)
@@ -946,8 +949,8 @@ class FixedPatchStepper:
if exc[0] is None:
raise
# Read-only compatibility views over the level lists. Internal code must not use them;
# REMOVE AFTER the existing tests/docs migrate to hierarchy/_layouts/_level_* lists.
# Read-only compatibility views over the level lists, kept on purpose for the existing
# tests and docs. Internal code must not use them; it reads hierarchy/_layouts/_level_*.
patch = property(lambda self: self.hierarchy.pairs[0])
child_patch = property(lambda self: self.hierarchy.pairs[1] if self.hierarchy.depth > 2 else None)
base_layout = property(lambda self: self._layouts[0])
+6 -4
View File
@@ -8,7 +8,7 @@ from dataclasses import dataclass, replace
import numpy as np
from .blocks import blocked_link_tables
from ..body.jets import jet_link_shapes, surface_link_velocity
from ..body.jets import _surface_links
from .descriptors import D2Q9_EX, D2Q9_EY
from .refinement import ACTIVE, OUTER_GHOST
@@ -142,11 +142,13 @@ def _prepare_blocked_boundary(boundary, layout, cfg, jets, *, _shifted_blocked=F
tables.update(body_id=np.asarray(body_id, np.int32), rx=np.asarray(rx, np.float32),
ry=np.asarray(ry, np.float32), aperture=np.zeros(len(dirs), bool))
actuated = bool(jets) or any(o.surface_actuation is not None for o in boundary.objects)
jet_cache = None
if actuated:
cfg = replace(cfg)
cfg.validate()
tables["surface_velocity"], tables["aperture"] = surface_link_velocity(
tables["surface_velocity"], tables["aperture"], shape, jet_id = _surface_links(
boundary.objects, cfg, boundary.flags, fluid_idx, dirs, body_id, rx, ry,
jets=jets, hits=boundary.curved_hits, _shifted_blocked=_shifted_blocked)
return _BlockedBoundary(boundary, tables, jets, actuated,
jet_link_shapes(jets, boundary.curved_hits) if jets else None)
if jets:
jet_cache = (shape, jet_id)
return _BlockedBoundary(boundary, tables, jets, actuated, jet_cache)
@@ -0,0 +1,66 @@
"""Checkpoint save -> load round trip for a scene with passive bodies. Requires GPU.
A fresh Simulation with the same scene, restored from the checkpoint must hold the same populations,
flags, body count and step count, and must stay bitwise identical to the original
after both continue for the same number of steps. Checkpoints with Jets are
refused (documented limitation) and are not covered here.
"""
import os
import tempfile
import unittest
import numpy as np
from ._cases import make_sim
# Three separated cylinders (pinball-like); cut-link receivers stay in bulk fluid.
CIRCLES = (((40, 32), 6), ((56, 22), 6), ((56, 42), 6))
def scene():
sim = make_sim()
sim.initialize()
for center, radius in CIRCLES:
sim.add_body("circle", center=center, radius=radius)
sim.sync_bodies()
return sim
class CheckpointBodiesGPU(unittest.TestCase):
def test_round_trip_is_bitwise_and_continues_identically(self):
sim = scene()
sim.run(60)
with tempfile.NamedTemporaryFile(suffix=".h5", delete=False) as f:
path = f.name
restored = None
try:
saved = sim.save_checkpoint(path)
# The loader requires the same body count and geometry: rebuild the scene.
restored = scene()
restored.load_checkpoint(saved)
self.assertEqual(restored.bodies.count, len(CIRCLES))
self.assertEqual(restored.stepper.step_count, sim.stepper.step_count)
np.testing.assert_array_equal(restored.get_flags(), sim.get_flags())
a, b = sim.get_ddf(), restored.get_ddf()
self.assertEqual(a.dtype, b.dtype)
self.assertEqual(a.tobytes(), b.tobytes())
sim.run(40)
restored.run(40)
self.assertEqual(restored.stepper.step_count, sim.stepper.step_count)
a, b = sim.get_ddf(), restored.get_ddf()
self.assertTrue(np.isfinite(a).all())
self.assertEqual(a.tobytes(), b.tobytes())
for i in range(len(CIRCLES)):
self.assertTrue(np.all(np.isfinite(restored.read_force(i))))
finally:
if restored is not None:
restored.close()
os.unlink(path)
sim.close()
if __name__ == "__main__":
unittest.main()
+2 -1
View File
@@ -16,6 +16,7 @@
#define NZ 1
#define NT 1
#define COLLISION_MODEL 0
#define BLOCK_ROOT_BC 2
#define USE_LES 0
#define RHO 1.f
#define U0 .025f
@@ -40,7 +41,7 @@ inline float __fdiv_rn(float a,float b){return a/b;}
#include "boundary/inlet_outlet.cuh"
#include "step/one_step_double.cu"
#else
constexpr unsigned char KIND_OUTER=2;
#include "core/cell_kind.cuh"
#include "step/block_channel.cu"
#endif
int main(){
+46
View File
@@ -0,0 +1,46 @@
"""Check historical source manifests against their own recorded commit, not today's checkout."""
import io
from pathlib import Path
import subprocess
import tarfile
ROOT = Path(__file__).resolve().parents[2]
def _require(commit):
"""Missing recorded history is an explicit failure, never a skip or a repin."""
found = subprocess.run(
["git", "-C", str(ROOT), "cat-file", "-e", commit + "^{commit}"],
capture_output=True,
)
if found.returncode:
raise AssertionError(
f"Required history missing: commit {commit} is not in this clone; "
"fetch it instead of repinning or substituting"
)
def tracked(commit, prefix="src"):
"""Tracked paths under prefix at the recorded commit."""
_require(commit)
return set(
subprocess.check_output(
["git", "-C", str(ROOT), "ls-tree", "-r", "--name-only", commit, prefix],
text=True,
).splitlines()
)
def historical_tree(commit, dest):
"""Extract the tracked tree of the recorded commit into dest."""
_require(commit)
archive = subprocess.check_output(
["git", "-C", str(ROOT), "archive", "--format=tar", commit]
)
dest = Path(dest)
dest.mkdir(parents=True, exist_ok=True)
extra = {"filter": "data"} if hasattr(tarfile, "data_filter") else {}
with tarfile.open(fileobj=io.BytesIO(archive)) as tar:
tar.extractall(dest, **extra)
return dest
+4 -2
View File
@@ -1,6 +1,7 @@
"""Host contract of PeriodicBlockStepper.step(n): one stream wait per call (mocked driver)."""
import ast
from contextlib import nullcontext
import importlib.util
from pathlib import Path
import types
from unittest.mock import Mock
@@ -8,14 +9,15 @@ import unittest
import numpy as np
import CelerisLab.lbm.block_stepper as b
from CelerisLab.config import LBMConfig
from CelerisLab.lbm.blocks import BlockLayout2D
def host_module():
"""block_stepper with pycuda/compiler replaced by mocks; no GPU context."""
path = Path(b.__file__)
# Locate the source without importing it: importing block_stepper loads pycuda.driver,
# which CPU-only import-boundary tests (test_build_contracts) must never see.
path = Path(importlib.util.find_spec("CelerisLab.lbm.block_stepper").origin)
tree = ast.parse(path.read_text())
tree.body = [node for node in tree.body if not (
isinstance(node, ast.Import) and any(a.name == "pycuda.driver" for a in node.names)
+26
View File
@@ -83,6 +83,32 @@ assert abs(sum(compute_feq_d2q9(1.,0.,0.)) - 1.) < 1e-6
class BuildContracts(unittest.TestCase):
def test_block_root_modes_are_rendered_and_recorded(self):
with tempfile.TemporaryDirectory() as root:
for mode in (0, 1, 2):
with self.subTest(mode=mode):
build = compiler.generate_config(LBMConfig(), build_root=root, block_root_bc=mode)
try:
method = Path(build, "config/config_method.h").read_text()
self.assertIn(f"#define BLOCK_ROOT_BC {mode}\n", method)
metadata = json.loads(Path(build, ".celeris-build").read_text())
self.assertEqual(metadata["block_root_bc"], mode)
self.assertNotIn("block_root_bc", metadata["config"])
finally:
compiler.cleanup_build_dir(build)
build = compiler.generate_config(LBMConfig(), build_root=root)
try:
self.assertEqual(json.loads(Path(build, ".celeris-build").read_text())["block_root_bc"], 0)
finally:
compiler.cleanup_build_dir(build)
def test_invalid_block_root_mode_fails_before_snapshot(self):
for mode in (-1, 3, True, False, 1.0, "1", None):
with self.subTest(mode=mode), patch.object(compiler.tempfile, "mkdtemp") as allocate:
with self.assertRaisesRegex(ValueError, "block_root_bc"):
compiler.generate_config(LBMConfig(), block_root_bc=mode)
allocate.assert_not_called()
@staticmethod
def fake_compile(command, check):
source = Path(command[-3])
+2 -1
View File
@@ -50,7 +50,7 @@ def test_masks_and_minimal_footprint():
def test_actual_closure_read_conflicts():
pair, layout, targets = transfer((6, 6, 18, 14))
validate_channel_transfers(pair, layout, targets)
for target in (layout.slot(layout.owner[5, -1]-1, 2),
for target in (layout.slot(int(layout.owner[5, -1]) - 1, 2), # integer under NumPy 1 and 2
layout.slot(layout.owner[1, 4], 4),
layout.slot(layout.owner[-2, 4], 3)):
with pytest.raises(ValueError, match="immutable"):
@@ -97,5 +97,6 @@ def test_evolved_ghost_snapshot_and_physical_download():
with pytest.raises(ValueError, match="finite"):
backend.initialize(*arrays)
assert backend._levels[0].kernel == host.compiler.load_module.return_value.get_function("BlockChannelStep")
assert [call.kwargs["block_root_bc"] for call in host.compiler.generate_config.call_args_list] == [2, 0]
finally:
backend.close()
+9 -1
View File
@@ -8,6 +8,7 @@ import numpy as np
ROOT = Path(__file__).resolve().parents[2]
REPO = ROOT
sys.path.insert(0, str(ROOT / "tests/validation"))
sys.path.insert(0, str(ROOT / "tests/unit"))
sys.path.insert(0, str(REPO / "src"))
from run_jet_circle_cv import (
CX,
@@ -19,6 +20,7 @@ from run_jet_circle_cv import (
cv_from_macros,
aperture_probe,
measured_aperture_flux,
PIN,
verify,
)
@@ -166,7 +168,13 @@ class Host(unittest.TestCase):
self.assertEqual(a["steps"] * 4, b["steps"])
def test_source_pin(self):
verify(REPO)
from tempfile import TemporaryDirectory
from historical_tree import historical_tree
with self.assertRaisesRegex(RuntimeError, "Source hash mismatch"):
verify(REPO)
with TemporaryDirectory() as directory:
verify(historical_tree(PIN, directory))
if __name__ == "__main__":
+27 -16
View File
@@ -15,6 +15,10 @@ sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "validation"))
from run_uniform_les_flow import source_guard, decode_macros, HERE, HEADER
ROOT = HERE.parents[1]
sys.path.insert(0, str(ROOT / "tests/unit"))
from historical_tree import historical_tree, tracked
A135 = "a135252ae78b47b57d7e77d19001700b3bedac8b"
class GuardTests(unittest.TestCase):
@@ -42,19 +46,23 @@ class GuardTests(unittest.TestCase):
decode_macros(invalid, w, cx, cy)
def test_current_corrected_pins(self):
self.assertEqual(
source_guard(ROOT, "corrected", "a135252")["base_commit"],
"a135252ae78b47b57d7e77d19001700b3bedac8b",
)
with TemporaryDirectory() as directory:
root = historical_tree(A135, directory)
with patch("run_uniform_les_flow.subprocess.check_output", return_value="historical-tree\n"):
identity = source_guard(root, "corrected", "a135252")
# base_commit comes from the pins file; the mock only feeds checkout_commit.
self.assertEqual(identity["base_commit"], A135)
self.assertEqual(identity["checkout_commit"], "historical-tree")
def test_corrected_and_tamper_guard(self):
pins = json.loads((HERE / "data/les_flow/source_pins_a135252.json").read_text())
with TemporaryDirectory() as directory:
root = Path(directory)
origin = historical_tree(A135, Path(directory) / "origin")
root = Path(directory) / "fixture"
for name in pins["sources"]:
dest = root / name
dest.parent.mkdir(parents=True, exist_ok=True)
shutil.copyfile(ROOT / name, dest)
shutil.copyfile(origin / name, dest)
with patch(
"run_uniform_les_flow.subprocess.check_output",
return_value="isolated-fixture\n",
@@ -87,23 +95,26 @@ class GuardTests(unittest.TestCase):
def test_current_manifest_pins_actual_corrected_tree(self):
pins = json.loads((HERE / "data/les_flow/source_pins_a135252.json").read_text())
tracked = subprocess.check_output(["git", "-C", str(ROOT), "ls-files", "src"], text=True)
self.assertEqual(set(pins["sources"]), set(tracked.splitlines()))
for name, expected in pins["sources"].items():
self.assertEqual(hashlib.sha256((ROOT / name).read_bytes()).hexdigest(), expected)
self.assertEqual(pins["sources"][HEADER], pins["candidate_header_sha256"])
self.assertNotEqual(pins["sources"][HEADER], pins["baseline_header_sha256"])
with self.assertRaisesRegex(RuntimeError, "Source pin mismatch"):
source_guard(ROOT, "baseline", "a135252")
self.assertEqual(pins["base_commit"], A135)
self.assertEqual(set(pins["sources"]), tracked(A135))
with TemporaryDirectory() as directory:
root = historical_tree(A135, directory)
for name, expected in pins["sources"].items():
self.assertEqual(hashlib.sha256((root / name).read_bytes()).hexdigest(), expected)
self.assertEqual(pins["sources"][HEADER], pins["candidate_header_sha256"])
self.assertNotEqual(pins["sources"][HEADER], pins["baseline_header_sha256"])
with self.assertRaisesRegex(RuntimeError, "Source pin mismatch"):
source_guard(root, "baseline", "a135252")
def test_named_baseline_accepts_only_reversed_header(self):
pins = json.loads((HERE / "data/les_flow/source_pins_a135252.json").read_text())
with TemporaryDirectory() as directory:
root = Path(directory)
origin = historical_tree(A135, Path(directory) / "origin")
root = Path(directory) / "fixture"
for name in pins["sources"]:
dest = root / name
dest.parent.mkdir(parents=True, exist_ok=True)
shutil.copyfile(ROOT / name, dest)
shutil.copyfile(origin / name, dest)
subprocess.run(["git", "apply", "--reverse", str(HERE / "les/data/smagorinsky_normalization.patch")], cwd=root, check=True)
with patch("run_uniform_les_flow.subprocess.check_output", return_value="baseline-fixture\n"):
source_guard(root, "baseline", "a135252")
+1 -1
View File
@@ -162,7 +162,7 @@ class ProjectedKernelTests(unittest.TestCase):
targets = child.owner.ravel()[p.ghost_fine]
donors = parent.owner.ravel()[p.ghost_coarse]
source = (Path(__file__).resolve().parents[2] /
"src/CelerisLab/lbm/kernels/blocks.cu").read_text()
"src/CelerisLab/lbm/kernels/step/block_transfer.cu").read_text()
kernel = source.split('extern "C" __global__ void BlockExplodeProjected(',1)[1]
kernel = 'void BlockExplodeProjected(' + kernel.split('// Targets are receiver-local SOURCE slots',1)[0]
def array(name, values, kind):
+6 -2
View File
@@ -5,6 +5,7 @@ import pytest
ROOT = Path(__file__).resolve().parents[2]
sys.path.insert(0, str(ROOT / "tests/validation"))
sys.path.insert(0, str(ROOT / "tests/unit"))
import re100_helpers as h
import periodic as p
@@ -62,13 +63,16 @@ def test_periodic_synthetic():
assert not p.assess(rows()[:100])["periodic"]
def test_production_unmodified():
def test_production_unmodified(tmp_path):
from historical_tree import historical_tree
m = json.loads(
(ROOT / "tests/validation/data/dfg_re100/source_manifest.json").read_text()
)
assert m["commit"] == h.PIN
origin = historical_tree(m["commit"], tmp_path)
for name, digest in m["production_sha256"].items():
assert hashlib.sha256((ROOT / name).read_bytes()).hexdigest() == digest, name
assert hashlib.sha256((origin / name).read_bytes()).hexdigest() == digest, name
def test_reference():
+53 -6
View File
@@ -7,7 +7,10 @@ import numpy as np
import unittest
sys.path.insert(0, str(Path(__file__).resolve().parents[1] / "validation"))
sys.path.insert(0, str(Path(__file__).resolve().parent))
import run_refined_dfg as r
from historical_tree import historical_tree
from tempfile import TemporaryDirectory
h, periodic = r.helpers(r.DEFAULT_SOURCE)
@@ -17,11 +20,24 @@ def test_source_pin():
assert historical["commit"].startswith("cd3a657")
mrt_history = json.loads((r.HERE / "data/refined_dfg/source_manifest_mrt.json").read_text())
assert mrt_history["commit"].startswith("3ea131e")
assert r.verify_source(r.DEFAULT_SOURCE, "MRT", "window")["commit"].startswith("a233c90")
with unittest.TestCase().assertRaisesRegex(RuntimeError, "Source mismatch"):
r.verify_source(r.DEFAULT_SOURCE, "MRT", "window")
with unittest.TestCase().assertRaisesRegex(RuntimeError, "Source mismatch"):
r.verify_source(r.DEFAULT_SOURCE, "MRT")
with unittest.TestCase().assertRaisesRegex(RuntimeError, "Source mismatch"):
r.verify_source(r.DEFAULT_SOURCE, "SRT")
with TemporaryDirectory() as directory:
for collision, mode, prefix in (
("SRT", "last_base_step", "cd3a657"),
("MRT", "last_base_step", "3ea131e"),
("MRT", "window", "a233c90"),
):
name = {"SRT": "source_manifest.json", "MRT": "source_manifest_mrt.json"}[collision]
name = "source_manifest_window.json" if mode == "window" else name
commit = json.loads((r.HERE / "data/refined_dfg" / name).read_text())["commit"]
assert commit.startswith(prefix)
root = historical_tree(commit, Path(directory) / commit)
assert r.verify_source(root, collision, mode)["commit"] == commit
def test_matched_units_and_planes():
@@ -141,16 +157,47 @@ def test_axis_specific_sampling_phase():
np.testing.assert_allclose(fine_y - ((k + 43) - 0.5) * 0.00125, 0, atol=1e-15)
def frozen_manifest(observation_mode):
"""Frozen MRT manifest JSON as recorded, without hashing the working tree."""
name = "source_manifest_window.json" if observation_mode == "window" else "source_manifest_mrt.json"
return json.loads((r.HERE / "data/refined_dfg" / name).read_text())
def selected_manifest(collision, observation_mode):
"""verify_source's manifest selection with its working-tree hash comparison neutralised.
Source identity is owned by test_source_pin; this exercises only which frozen manifest
verify_source selects and returns for (collision, observation_mode).
"""
from unittest.mock import patch
class Digest(str):
def __ne__(self, other):
return False
class Hash:
def __init__(self, data=b""):
pass
def hexdigest(self):
return Digest()
with patch.object(r.hashlib, "sha256", Hash):
return r.verify_source(r.DEFAULT_SOURCE, collision, observation_mode)
def test_full_requires_completed_pilot():
with unittest.TestCase().assertRaisesRegex(ValueError, "pilot-results"):
r.verify_pilot(None, r.verify_source(r.DEFAULT_SOURCE, "MRT", "window"))
r.verify_pilot(None, frozen_manifest("window"))
def pilot_fixture(path, observation_mode="window", jet_scenario="none", case_selection="both", uniform_diameter=80, placement="wide"):
"""Synthetic contract fixture only; never a physical result or run receipt."""
# Receipt-logic fixture: the frozen manifest verify_source returns for a matching tree.
# Source identity itself is checked only by test_source_pin.
manifest = (r.verify_source(r.DEFAULT_SOURCE, "MRT", "window", jet_scenario)
if observation_mode == "window" else
json.loads((r.HERE / "data/refined_dfg/source_manifest_mrt.json").read_text()))
if observation_mode == "window" and jet_scenario != "none" else
frozen_manifest(observation_mode))
selection = dict(case_selection=case_selection, uniform_diameter=uniform_diameter, placement=placement)
(path / "plan.json").write_text(
json.dumps(
@@ -348,7 +395,7 @@ def test_record_actual_checkout_separately():
actual = checkout_provenance(r.DEFAULT_SOURCE)
assert len(actual["checkout_revision"]) == 40
assert isinstance(actual["checkout_dirty"], bool)
assert r.verify_source(r.DEFAULT_SOURCE, "MRT", "window")["production_origin_commit"].startswith("d930")
assert selected_manifest("MRT", "window")["production_origin_commit"].startswith("d930")
def test_window_plan_preserves_physics():
@@ -358,7 +405,7 @@ def test_window_plan_preserves_physics():
assert baseline.pop("observation_mode") == "last_base_step"
assert window.pop("observation_mode") == "window"
assert baseline == window
assert r.verify_source(r.DEFAULT_SOURCE, collision, "window")["observation_mode"] == "window"
assert selected_manifest(collision, "window")["observation_mode"] == "window"
with unittest.TestCase().assertRaises(ValueError):
r.verify_source(r.DEFAULT_SOURCE, "MRT", "unknown")
+14 -5
View File
@@ -22,13 +22,16 @@ class RefinementMRTHost(unittest.TestCase):
def test_srt_statements_and_guards_preserved(self):
"""MRT is an additive compile-time branch; SRT statements and other guards remain."""
blocks = (KERNELS / "blocks.cu").read_text()
blocks = (KERNELS / "blocks.cu").read_text() + (KERNELS / "step/block_step.cu").read_text()
channel = (KERNELS / "step/block_channel.cu").read_text()
refinement = (KERNELS / "refinement.cu").read_text()
block_collision = (KERNELS / "operators/block_collision.cuh").read_text()
self.assertIn('#include "operators/block_collision.cuh"', blocks)
self.assertEqual(blocks.count("block_collide(f);"), 1)
refinement = (KERNELS / "refinement.cu").read_text() + (KERNELS / "step/patch_step.cu").read_text()
for text in ("compute_feq_from_momentum(drho, jx, jy, ux, uy, feq);",
"compute_feq(rho, ux, uy, feq);",
"collide_srt_no_force(f, feq, d_params.omega);"):
self.assertIn(text, blocks)
self.assertIn(text, block_collision)
self.assertIn("collide_srt(f, feq, forcing, omega);", channel)
self.assertIn("collide_srt_no_force(f, feq, d_params.omega);", refinement)
for source in (blocks, refinement):
@@ -37,8 +40,14 @@ class RefinementMRTHost(unittest.TestCase):
self.assertIn("STORE_PRECISION != 0", source)
self.assertIn('#if COLLISION_MODEL == 2\n#include "operators/collision_mrt.cuh"\n#endif', source)
self.assertIn("USE_DDF_SHIFTING ||", refinement) # dense/mixed shifted stays rejected
for source in (blocks, channel, refinement):
self.assertEqual(source.count("collide_mrt_no_force("), 1)
# 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, refinement):
self.assertEqual(source.count("collide_mrt_unforced(f);"), 1)
self.assertNotIn("collide_mrt_no_force(", source)
self.assertNotIn("compute_rho_u(f, rho, ux, uy);\n collide_mrt", source)
self.assertNotRegex(source, r"compute_rho_u\(f, rho, ux, uy\);\s*collide_mrt")
if __name__ == "__main__": unittest.main()
+4 -2
View File
@@ -78,13 +78,15 @@ class ShiftedStorage(unittest.TestCase):
def test_shifted_density_and_tag_validation(self):
old,new = topology(A),topology(B)
a = tuple(np.where(m[None], np.float32(-.001), np.nan).repeat(9,axis=0) for m in old.leaf_masks)
a = tuple(np.where(m[None], np.float32(-.001), np.float32(np.nan)).repeat(9,axis=0)
for m in old.leaf_masks) # FP32 under NumPy 1 and 2 promotion rules
remap(old,new,a)
with self.assertRaisesRegex(ValueError,'representation'): remap(old,new,a,False)
with self.assertRaisesRegex(ValueError,'FP32'): remap(old,new,tuple(x.astype(np.float64) for x in a))
a[0][:,old.leaf_masks[0]] = -1
with self.assertRaisesRegex(ValueError,'positive density'): remap(old,new,a)
rest = tuple(np.where(m[None],np.float32(0),np.nan).repeat(9,axis=0) for m in old.leaf_masks)
rest = tuple(np.where(m[None],np.float32(0),np.float32(np.nan)).repeat(9,axis=0)
for m in old.leaf_masks)
remap(old,new,rest)
with self.assertRaisesRegex(ValueError,'positive density'):
remap_synchronized(old,new,rest,base_count=7,level_steps=(7,14,28),phase='postcollision_g')
@@ -12,13 +12,15 @@ from source_provenance import checkout_provenance
import re100_helpers
import run_jet_circle_cv
sys.path.insert(0, str(ROOT / "tests/unit"))
from historical_tree import historical_tree, tracked
ORIGIN = "47b05bc8abfc957fe876b0c63bdb386e3c6cbaa5"
def test_manifests_cover_current_source_and_preserve_origin():
names = subprocess.check_output(
["git", "-C", str(ROOT), "ls-files", "src"], text=True
).splitlines()
def test_manifests_cover_current_source_and_preserve_origin(tmp_path):
assert re100_helpers.PIN == run_jet_circle_cv.PIN == ORIGIN
manifests = {}
for group, hashes_key, revision_key in (
("dfg_re100", "production_sha256", "commit"),
("jet_circle_cv", "files", "revision"),
@@ -29,18 +31,22 @@ def test_manifests_cover_current_source_and_preserve_origin():
).read_text()
)
assert manifest[revision_key] == ORIGIN
assert set(manifest[hashes_key]) == set(names)
assert manifest["source_origin_commits"][-1] == ORIGIN
assert (
manifest["historical_checkpoints"]["production"]
== "d79e3f6f0a424f5bcc3a1a6d610313e49daf97a2"
)
assert "not inherited from d79" in manifest["experiment"]
for name, expected in manifest[hashes_key].items():
manifests[group] = manifest[hashes_key]
# History-dependent checks last, so a missing commit masks nothing above.
names = tracked(ORIGIN)
origin = historical_tree(ORIGIN, tmp_path)
for group, hashes in manifests.items():
assert set(hashes) == set(names), group
for name, expected in hashes.items():
assert (
hashlib.sha256((ROOT / name).read_bytes()).hexdigest() == expected
hashlib.sha256((origin / name).read_bytes()).hexdigest() == expected
), name
assert re100_helpers.PIN == run_jet_circle_cv.PIN == ORIGIN
def test_runtime_checkout_identity_is_separate_from_source_origin(tmp_path):
@@ -192,6 +192,7 @@ class CpuContracts(unittest.TestCase):
backend._lock = RLock()
backend._status = "ready"
backend._step_count = 0
backend._ddf_shifting = False # unshifted storage, as for the default config
backend._stream = Mock()
backend._module = Mock()
backend._copy = backend._step = Mock()
@@ -27,7 +27,7 @@ def state(ny, nx):
def compile_module(cfg, context, source):
cfg.validate()
build = compiler.generate_config(cfg)
build = compiler.generate_config(cfg, block_root_bc=2 if source == "blocks.cu" else 0)
try:
module = compiler.load_module(compiler.compile_kernel(context.sm_arch, source=source, build_dir=build))
finally:
@@ -282,7 +282,8 @@ def test_fine_surface_four_hops(speed):
@pytest.mark.parametrize("missing", ["surface_velocity_gpu", "link_exchange_gpu"])
def test_missing_surface_device_buffer_rejects_before_launch(missing):
backend = implementation.FixedPatchStepper.__new__(implementation.FixedPatchStepper)
backend.fine_boundary = SimpleNamespace(objects=[object_spec(.002)])
backend._boundary_hosts = [None, SimpleNamespace(objects=[object_spec(.002)])]
backend._level_jets = [(), ()]
backend._fine_curved = SimpleNamespace(fluid_idx_gpu=object(),
surface_velocity_gpu=object(), link_exchange_gpu=object())
setattr(backend._fine_curved, missing, None)
@@ -223,6 +223,8 @@ class BlockedBaseCPU(unittest.TestCase):
self.assertIs(collide[4], coarse.kind)
walls = [a for ptr, a in m["uploads"] if isinstance(ptr, str) and ptr == "d_block_wall"]
self.assertEqual(len(walls), 1)
self.assertEqual([call.kwargs["block_root_bc"] for call in
compiler.generate_config.call_args_list], [1, 0])
width, height, slot, speed = struct.unpack("=IIIf", walls[0])
self.assertEqual((width, height, slot), (NX, NY, 8))
self.assertAlmostEqual(speed, .001)
@@ -108,7 +108,9 @@ class FaultGuardsCPU(unittest.TestCase):
def fake_backend(self):
backend = PeriodicPatchStepper.__new__(PeriodicPatchStepper)
backend.fine_layout = None # dense fine storage
backend._layouts = [None, None] # dense storage on both levels
backend._ddf_shifting = False # unshifted on every level, as _check_storage expects
backend._level_cfgs = [LBMConfig(nx=16, ny=16, viscosity=.1)] * 2 # ddf_shifting off
backend.ctx = Mock(activate=nullcontext)
backend._lock = RLock()
backend._status = "ready"
@@ -158,10 +160,21 @@ class FaultGuardsCPU(unittest.TestCase):
with self.assertRaisesRegex(RuntimeError, "four-cell"):
PeriodicPatchStepper(cfg, context=context, refinements=[(4, 4, 12, 12)])
p = implementation.FixedPatch2D.from_config(cfg, (4, 4, 12, 12))
with patch.object(cuda, "Stream"), patch.object(implementation.FixedPatch2D, "from_config",
return_value=replace(p, tau_coarse=.9)):
with self.assertRaisesRegex(RuntimeError, "omega/tau"):
from CelerisLab.lbm import hierarchy
bad = patch.object(implementation.FixedPatch2D, "from_config",
return_value=replace(p, tau_coarse=.9))
# A pair tau inconsistent with its level config is rejected twice over, both before
# any device allocation or upload: first by the hierarchy chain check, and, if that
# is bypassed, by the per-level omega/tau check in _build_level.
with patch.object(cuda, "Stream"), bad, patch.object(cuda, "mem_alloc") as alloc, \
patch.object(cuda, "memcpy_htod") as upload:
with self.assertRaisesRegex(RuntimeError, "config chain is inconsistent"):
PeriodicPatchStepper(cfg, context=context, refinements=[p.bounds])
with patch.object(hierarchy, "_check_chain"), \
self.assertRaisesRegex(RuntimeError, "omega/tau"):
PeriodicPatchStepper(cfg, context=context, refinements=[p.bounds])
alloc.assert_not_called()
upload.assert_not_called()
module = Mock()
module.get_global.return_value = (Mock(), 16)
with patch.object(cuda, "Stream"), patch.object(compiler, "generate_config"), \
@@ -249,19 +249,41 @@ class ProjectedPairsTests(unittest.TestCase):
self.assertFalse(assess(incomplete)["comparisons"][3]["passed"])
def test_pinned_build_sources_match_host_generated_inputs(self):
from CelerisLab.cuda import compiler_v2
"""Check the pinned fingerprints using NUMERICAL_REVISION's renderer and kernels.
The imported LBMConfig is unchanged from that revision; this is a host
build-input check, not a reconstruction of the historical CUDA toolchain."""
import importlib.util
import tarfile
for size in (128, 256):
cfg = replace(run_case.config(size, 100), nx=10, ny=10)
cfg.validate()
root = compiler_v2.generate_config(cfg)
try:
actual = json.loads((Path(root) / ".celeris-build").read_text())
expected = build_metadata(size)
self.assertEqual(actual, {key: expected[key] for key in
("source_sha256", "config", "n_objects")})
finally:
compiler_v2.cleanup_build_dir(root)
kernels = "src/CelerisLab/lbm/kernels"
renderer = "src/CelerisLab/cuda/compiler_v2.py"
if subprocess.run(["git", "-C", str(REPO), "cat-file", "-e",
f"{NUMERICAL_REVISION}^{{commit}}"], capture_output=True).returncode:
self.fail(f"{NUMERICAL_REVISION} is not in this clone; the pinned build "
"sources cannot be checked")
archive = subprocess.run(["git", "-C", str(REPO), "archive", "--format=tar",
NUMERICAL_REVISION, kernels, renderer],
capture_output=True, check=True).stdout
with tempfile.TemporaryDirectory() as snapshot:
with tarfile.open(fileobj=io.BytesIO(archive)) as tar:
tar.extractall(snapshot)
spec = importlib.util.spec_from_file_location(
"CelerisLab.cuda._compiler_a135252", Path(snapshot) / renderer)
compiler = importlib.util.module_from_spec(spec)
spec.loader.exec_module(compiler)
with patch.object(compiler, "_KERNEL_DIR", str(Path(snapshot) / kernels)):
for size in (128, 256):
with self.subTest(size=size):
cfg = replace(run_case.config(size, 100), nx=10, ny=10)
cfg.validate()
root = compiler.generate_config(cfg)
try:
actual = json.loads((Path(root) / ".celeris-build").read_text())
expected = build_metadata(size)
self.assertEqual(actual, {key: expected[key] for key in
("source_sha256", "config", "n_objects")})
finally:
compiler.cleanup_build_dir(root)
def test_source_preflight_requires_clean_exact_revision_and_tree(self):
import oracle