diff --git a/docs/synchronized_leaf_remap.md b/docs/synchronized_leaf_remap.md index 9d74d00..dede90e 100644 --- a/docs/synchronized_leaf_remap.md +++ b/docs/synchronized_leaf_remap.md @@ -1,11 +1,11 @@ -# Synchronized leaf G remap (CPU only) +# Synchronized leaf G remap and manual continuation This is a CPU state-transfer component, **not a runtime adaptive solver**. It changes only a separated L2 rectangle set in a three-level, body-free, periodic hierarchy. L0/L1 topology, relaxation parameters, direction order, physical frame and completed time remain unchanged. Existing geometry checks -validate child support/nesting and separation. No CUDA, scheduler, allocator, -boundary handling, collision, ghost refill or runtime installation is added. +validate child support/nesting and separation. The CPU helper itself adds no CUDA, scheduler, allocator, boundary handling, +collision or ghost refill. A bounded runtime adapter is described below. ```python from CelerisLab.lbm.refinement import FixedPatchSet2D @@ -56,7 +56,7 @@ dtype. New nonleaves are NaN. Source arrays/metadata are not changed. and new `topology.leaf_masks`. Existing `FixedPatchStepper.initialize` resets the counter: using it is **not** a continuation adapter. This component neither publishes a candidate nor retires old state. Runtime table rebuild/halo refill -and correct first-step scheduling remain separate integration work. +and correct first-step scheduling are handled only by the explicit adapter below. Tests accumulate leaf-area integrals in FP64 for each population, mass and momentum. The declared absolute tolerance is `2e-13 * sum(abs(old population @@ -68,3 +68,44 @@ is not invertible. No kinetic accuracy or stability claim follows. ```sh PYTHONPATH=src python -m unittest discover -s tests/unit -p 'test_synchronized_remap.py' -v ``` + +## Explicit runtime continuation + +`replacement = current.remap_child(rectangles)` returns a new ready +`FixedPatchStepper` at the **same absolute step_count**. It supports only +body-free periodic three-level runs whose existing L2 uses +`child_rectangles_in_parent`. Scalar or union L1 parents are supported; scalar +`child=` is deliberately rejected. L0/L1, relaxation and block sizes stay fixed. +No tagger, scheduler, in-place replacement, or general rollback is provided. + +The method holds the stepper lock, validates geometry, synchronizes its stream, +and downloads the actual current owned leaves. It certifies postcollision G +and completed 1:2:4 counts internally, applies the CPU helper, rebuilds all +geometry/transfer/kind/halo buffers, and uploads both current and spare buffers. +No old halo is reused. The regular step schedule fills transit ghosts/halos +before use. `initialize` itself still resets the count to zero. + +The old stepper remains ready, unchanged and independently usable. The caller +owns **both** objects and must close both before closing their shared CUDA +context. Candidate construction or upload failure releases only candidate +resources; an actual old-device synchronization/read failure marks the old +stepper failed, consistent with its existing failure semantics. + +```python +replacement = current.remap_child([(20, 16, 36, 40)]) +# Explicitly choose the returned continuation only after success. +current.close() +current = replacement +current.step() # advances the preserved absolute count by one +``` + +CPU lifecycle mocks cover fresh current-buffer reads, exact transfer, separate +ownership, time preservation, ordinary initialization reset, unsupported states +and geometry, and candidate/source failure cleanup. GPU tests are predeclared +in `tests/validation/multigrid/test_runtime_child_remap_gpu.py`: identity next +step bitwise, shifted A→B→A with two intervals between changes, uniform moving +and nonequilibrium states, exact CPU remap and overlap, both buffers reset, +finite leaves, exact masks and absolute counts. They were not run on the CPU +host. The separate independent FP64 capture/replay gate is max population +error **3e-6**, with exact masks; remap mass/momentum allow FP32 roundoff. +Per-population conservation is not an advancement target under collision. diff --git a/src/CelerisLab/lbm/refinement_stepper.py b/src/CelerisLab/lbm/refinement_stepper.py index f0cc2f2..55b8260 100644 --- a/src/CelerisLab/lbm/refinement_stepper.py +++ b/src/CelerisLab/lbm/refinement_stepper.py @@ -3,7 +3,8 @@ The caller owns the CudaContext and must close this backend before that context. Only synchronized postcollision_g leaf populations are public. No checkpoint, -snapshot, forcing, runtime parameter mutation or dynamic refinement is offered. +forcing or runtime parameter mutation is offered. Manual body-free periodic +L2 set replacement is available at completed base intervals via remap_child. """ from dataclasses import dataclass, replace from numbers import Integral @@ -307,6 +308,10 @@ class FixedPatchStepper: self._child = None if self.fine_boundary is not None: self._fine_leaf &= (self.fine_boundary.flags.reshape(self.patch.fine_shape) & SOLID) == 0 + # Private constructor inputs for the bounded manual L2 continuation. + self._remap_cfg = replace(cfg) + self._remap_blocks = tuple(None if b is None else tuple(b) + for b in (base_block, fine_block, child_block)) self.ctx = context self._lock = RLock() self._status = "new" @@ -489,6 +494,60 @@ class FixedPatchStepper: self._status = "failed" raise + def remap_child(self, rectangles): + """Return a NEW ready stepper with an L2 set replaced at this base time. + + Only body-free periodic three-level set children are supported; L0/L1 + and block sizes stay fixed. The lock and stream synchronization certify + current postcollision G and the completed 1:2:4 counts internally. + Both objects remain independently owned: the caller must close each, + and chooses when to stop using the unchanged old instance. Candidate + construction/upload errors never invalidate the old instance; failure + to synchronize/read its device state does, as with step/download. + """ + from .synchronized_remap import RemapTopology, remap_synchronized + + with self._lock: + self._require("ready") + if (self._wall_aware or len(self._levels) != 3 + or not isinstance(self.child_patch, FixedPatchSet2D) + or self.fine_boundary is not None or self.child_boundary is not None + or any(level.boundary is not None for level in self._levels)): + raise ValueError("remap_child requires body-free periodic three-level L2 sets") + cfg = self._remap_cfg + old = RemapTopology.from_config(cfg, self.patch, self.child_patch.bounds) + new = RemapTopology.from_config(cfg, self.patch, rectangles) + base_block, fine_block, child_block = self._remap_blocks + # Check sparse block alignment/launch size before any device read/allocation. + layout = new.child.fine_layout(child_block) + if layout.nblocks > MAX_GRID_DIM or layout.local_shape[0] > MAX_GRID_DIM: + raise ValueError("child blocks exceed launch limits") + with self.ctx.activate(): + try: + self._stream.synchronize() + source = self.download() + except BaseException: + self._status = "failed" + raise + count = self._step_count + transferred = remap_synchronized( + old, new, source, base_count=count, level_steps=(count, 2*count, 4*count), + phase="postcollision_g") + candidate = type(self)( + cfg, self.patch.bounds, context=self.ctx, base_block=base_block, + fine_block=fine_block, child_block=child_block, + child_rectangles_in_parent=new.child.bounds) + try: + candidate.initialize(*transferred.g) + candidate._step_count = count + except BaseException: + try: + candidate.close() + except Exception: + pass # Preserve the original candidate initialization failure. + raise + return candidate + def _apply_fine_boundary(self, buffer, **kw): """Stage solid source slots; exchange buffer holds this fine substep only.""" cl = self._fine_curved diff --git a/tests/unit/test_runtime_child_remap.py b/tests/unit/test_runtime_child_remap.py new file mode 100644 index 0000000..07d46be --- /dev/null +++ b/tests/unit/test_runtime_child_remap.py @@ -0,0 +1,119 @@ +"""CPU mocks execute production remap lifecycle; no GPU numerical claim.""" +import unittest +from unittest.mock import patch + +import numpy as np + +from tests.unit.test_refinement_three_level_multipatch import mock_backend +from tests.unit.test_synchronized_remap import A, B, fields, topology, remap + + +def ready_backend(): + run, module, allocations = mock_backend(children=A) + storage = {} + module.cuda.memcpy_htod.side_effect = lambda ptr, a: storage.__setitem__(id(ptr), np.array(a, copy=True)) + module.cuda.memcpy_dtoh.side_effect = lambda a, ptr: np.copyto(a, storage[id(ptr)]) + run.initialize(*fields(topology(A), np.float32)) + run._step_count = 7 + return run, module, allocations, storage + + +class RuntimeChildRemap(unittest.TestCase): + def test_current_snapshot_identity_and_count_continuation(self): + run, module, allocations, storage = ready_backend() + # Make current state distinct from initialization and spare, as after a step. + for level in run._levels: + storage[id(level.current)] *= np.float32(1.01) + before = run.download() + saved = {key: value.copy() for key, value in storage.items()} + old_allocations = list(allocations) + module.cuda.memcpy_dtoh.reset_mock() + clone = run.remap_child(A) + self.assertEqual(module.cuda.memcpy_dtoh.call_count, 3) + for a, b in zip(clone.download(), before): + np.testing.assert_array_equal(a, b) + self.assertEqual((run.step_count, clone.step_count), (7, 7)) + for key, value in saved.items(): + np.testing.assert_array_equal(storage[key], value) + for allocation in old_allocations: + allocation.free.assert_not_called() + clone.step() + self.assertEqual((run.step_count, clone.step_count), (7, 8)) + run.step() + self.assertEqual(run.step_count, 8) + clone.close(); run.close() + + def test_shift_uses_cpu_remap_and_resets_both_buffers(self): + run, module, allocations, storage = ready_backend() + before = run.download() + expected = remap(topology(A), topology(B), before).g + clone = run.remap_child(B) + for level, want, got in zip(clone._levels, expected, clone.download()): + np.testing.assert_array_equal(got, want) + for ptr in (level.current, level.spare): + np.testing.assert_array_equal(storage[id(ptr)], level.layout.to_blocked(want)) + for a, b in zip(run.download(), before): + np.testing.assert_array_equal(a, b) + clone.initialize(*expected) + self.assertEqual(clone.step_count, 0) # Existing public initialization semantics. + clone.close(); run.close() + + def test_invalid_state_geometry_and_body_reject_before_reads_allocations(self): + run, module, allocations, storage = ready_backend() + start = len(allocations) + module.cuda.memcpy_dtoh.reset_mock() + for status in ('new', 'failed', 'closed', 'advancing'): + run._status = status + with self.assertRaises(RuntimeError): run.remap_child(B) + run._status = 'ready' + for bad in ([], [(40,16,48,40)], [(16,16,32,24),(16,27,32,40)]): + with self.assertRaises(ValueError): run.remap_child(bad) + for attribute, value in (('_wall_aware', True), ('child_boundary', object()), + ('fine_boundary', object()), ('child_patch', None)): + with patch.object(run, attribute, value): + with self.assertRaises(ValueError): run.remap_child(B) + self.assertEqual(len(allocations), start) + module.cuda.memcpy_dtoh.assert_not_called() + run.close() + + def test_candidate_construction_or_upload_failure_preserves_source(self): + for failure in ('construction', 'initialize'): + run, module, allocations, storage = ready_backend() + before = run.download() + old_allocations = list(allocations) + if failure == 'construction': + operation = patch.object(module.compiler, 'compile_kernel', side_effect=RuntimeError('candidate compile')) + else: + initialize = type(run).initialize + def fail(candidate, *arrays): + upload = module.cuda.memcpy_htod.side_effect + calls = 0 + def partial_upload(pointer, array): + nonlocal calls + calls += 1 + if calls == 2: + raise RuntimeError('candidate upload') + upload(pointer, array) + with patch.object(module.cuda, 'memcpy_htod', side_effect=partial_upload): + initialize(candidate, *arrays) + operation = patch.object(type(run), 'initialize', fail) + with operation, self.assertRaisesRegex(RuntimeError, 'candidate'): + run.remap_child(B) + self.assertEqual((run._status, run.step_count), ('ready', 7)) + for a, b in zip(run.download(), before): np.testing.assert_array_equal(a, b) + for allocation in old_allocations: allocation.free.assert_not_called() + for allocation in allocations[len(old_allocations):]: allocation.free.assert_called_once() + run.close() + + def test_old_sync_or_read_failure_marks_only_old_failed(self): + for operation in ('sync', 'read'): + run, module, allocations, storage = ready_backend() + start = len(allocations) + target = run._stream.synchronize if operation == 'sync' else module.cuda.memcpy_dtoh + target.side_effect = RuntimeError('old device failure') + with self.assertRaisesRegex(RuntimeError, 'old device'): + run.remap_child(B) + self.assertEqual(run._status, 'failed') + self.assertEqual(len(allocations), start) + target.side_effect = None + run.close() diff --git a/tests/validation/multigrid/test_runtime_child_remap_gpu.py b/tests/validation/multigrid/test_runtime_child_remap_gpu.py new file mode 100644 index 0000000..07f5a72 --- /dev/null +++ b/tests/validation/multigrid/test_runtime_child_remap_gpu.py @@ -0,0 +1,106 @@ +"""Manual L2 continuation gates; run explicitly on an approved idle GPU. + +Frozen numerical gate: max |G - independent FP64 replay| <= 3e-6 on all +ACTIVE leaves, with exact masks and exact same-level overlap at remap. The +separate independent capture/replay acceptance harness supplies that numerical +comparison; these tests cover lifecycle, exact transfer and GPU continuation. +Remap mass/momentum must agree to FP32 roundoff, not per-q conservation during +advancement (collision changes individual populations). No GPU run is implied. +""" +from contextlib import ExitStack +import unittest + +import numpy as np +import pycuda.driver as cuda + +from CelerisLab.cuda.context import CudaContext +from CelerisLab.lbm.equilibrium import compute_feq_d2q9 +from CelerisLab.lbm.refinement_stepper import FixedPatchStepper +from CelerisLab.lbm.synchronized_remap import RemapTopology, remap_synchronized +from test_refinement_blocked_base import exact +from test_refinement_three_level_multipatch_gpu import CFG, PARENTS, frames, physical_state + +A = ((16, 16, 32, 40),) +B = ((20, 16, 36, 40),) +FP64_ATOL = 3e-6 + + +def make(ctx): + return FixedPatchStepper(CFG, PARENTS, context=ctx, base_block=(4, 4), + fine_block=(4, 4), child_block=(4, 4), + child_rectangles_in_parent=A) + + +def initial(run, uniform): + equilibrium = np.asarray(compute_feq_d2q9(1., .013, -.007), np.float32)[:, None, None] + return tuple(np.broadcast_to(equilibrium, level.shape).copy() if uniform else + physical_state(level.shape[1:], *frame) + for level, frame in zip(run._levels, frames(run))) + + +class RuntimeChildRemapGPU(unittest.TestCase): + def check_public(self, run): + for array, leaf in zip(run.download(), (~run.patch.covered, run._fine_leaf, run._child_leaf)): + self.assertTrue(np.isfinite(array[:, leaf]).all()) + self.assertTrue(np.isnan(array[:, ~leaf]).all()) + + def check_reset(self, run, expected): + for level, want in zip(run._levels, expected): + blocked = level.layout.to_blocked(want) + for pointer in (level.current, level.spare): + actual = np.empty(level.layout.blocked_shape, np.float32) + cuda.memcpy_dtoh(actual, pointer) + exact(self, actual, blocked, 'fresh current and spare, including halos') + + def test_identity_next_step_bitwise_and_old_still_usable(self): + with CudaContext(0) as ctx, ExitStack() as owned: + old = owned.enter_context(make(ctx)) + old.initialize(*initial(old, False)) + old.step(2) + before = old.download() + new = owned.enter_context(old.remap_child(A)) + self.assertEqual((old.step_count, new.step_count), (2, 2)) + self.check_reset(new, before) + for got, want in zip(old.download(), before): exact(self, got, want, 'old unchanged') + old.step(); new.step() + self.assertEqual((old.step_count, new.step_count), (3, 3)) + for got, want in zip(new.download(), old.download()): + exact(self, got, want, 'identity next step') + self.check_public(new) + + def test_shift_then_return_uniform_and_nonequilibrium(self): + for uniform in (True, False): + with self.subTest(uniform=uniform), CudaContext(0) as ctx, ExitStack() as owned: + current = owned.enter_context(make(ctx)) + current.initialize(*initial(current, uniform)) + current.step(2) + for rectangles in (B, A): + before, count = current.download(), current.step_count + old_mesh = RemapTopology.from_config(CFG, current.patch, current.child_patch.bounds) + new_mesh = RemapTopology.from_config(CFG, current.patch, rectangles) + expected = remap_synchronized(old_mesh, new_mesh, before, base_count=count, + level_steps=(count, 2*count, 4*count), phase='postcollision_g').g + previous = current + current = owned.enter_context(previous.remap_child(rectangles)) + self.assertEqual(current.step_count, count) + self.check_reset(current, expected) + for got, want, old_g, old_leaf, new_leaf in zip( + current.download(), expected, before, old_mesh.leaf_masks, new_mesh.leaf_masks): + exact(self, got, want, 'CPU remap exact') + overlap = old_leaf & new_leaf + exact(self, got[:, overlap], old_g[:, overlap], 'same-level overlap') + for got, want in zip(previous.download(), before): + exact(self, got, want, 'source unchanged') + current.step(2) + self.assertEqual(current.step_count, count+2) + self.assertEqual(previous.step_count, count) + self.check_public(current) + if uniform: + equilibrium = np.asarray(compute_feq_d2q9(1., .013, -.007), np.float64) + for array in current.download(): + leaf = np.isfinite(array[0]) + self.assertLessEqual(np.max(np.abs(array[:, leaf]-equilibrium[:, None])), FP64_ATOL) + + +if __name__ == '__main__': + unittest.main()