Clarify raw exchange metadata and streamline host tests

docs/development/raw_exchange.md

@@ -20,6 +20,14 @@ their geometric classification when their prescribed velocity is zero.
 All results are base-lattice units per unit out-of-plane depth. For local
 spacing and timestep h=dt, raw `(jx,jy,tz,dm)` converts to impulse
 `(h*dt*jx,h*dt*jy,h*h*dt*tz)` and mass `h*h*dm`.
+The accumulator intentionally accepts h and dt separately to make the spatial
+and temporal conversion explicit; the current acoustic hierarchy supplies dt=h.
+The stored rx/ry arms run from geometry.center to the actual wall hit, not to
+the fluid node. For this raw exchange j=c_i(a+b), moving the arm from the fluid
+node to the hit adds q*c_i, whose cross product with j is zero. Thus either arm
+gives the same raw torque (up to floating-point rounding); this equivalence
+must not be assumed for a future nonparallel stress/advective correction.
+
 Force/torque are impulse/angular impulse divided by the reported window length.
 `dm=b-a` is positive into the exterior fluid, while j is the opposite of the
 fluid's momentum change across that cut link. Interpolated BFL can introduce

src/CelerisLab/lbm/refinement_stepper.py

@@ -21,7 +21,8 @@ from .block_stepper import MAX_GRID_DIM
 from .hierarchy import build_hierarchy, is_rectangle
 from .refinement import ACTIVE, HALO, FixedPatch2D, FixedPatchSet2D, prepare_fine_boundary
 from .refinement_tables import (blocked_base_tables, blocked_fine_tables,  # noqa: F401 (re-export)
-                                blocked_kind, _prepare_blocked_boundary)
+                                blocked_kind, _prepare_blocked_boundary,
+                                BLOCKED_BOUNDARY_DEVICE_COLUMNS)
 from .descriptors import SOLID
 from .raw_exchange import summarize_exchange
 
@@ -335,8 +336,8 @@ class FixedPatchStepper:
                     b = level.boundary
                     if b is not None:
                         b.flags = self._upload(b.tables["flags"])
-                        columns = ("block", "cell", "dir", "q", "fallback", "surface_velocity")
+                        b.device = {k: self._upload(v) for k, v in b.tables.items()
-                        b.device = {k: self._upload(v) for k, v in b.tables.items() if k in columns}
+                                    if k in BLOCKED_BOUNDARY_DEVICE_COLUMNS}
                         if self._observe_raw_exchange:
                             b.device.update(rx=self._upload(b.tables["rx"]),
                                             ry=self._upload(b.tables["ry"]))

src/CelerisLab/lbm/refinement_tables.py

@@ -66,6 +66,11 @@ def blocked_base_tables(patch, layout, order, wall_aware=False, *, composed_kind
             blocked_kind(layout, dense_kind))
 
 
+# Required device columns plus optional actuation; observation arms are uploaded
+# separately only when enabled. body_id/aperture remain host reduction metadata.
+BLOCKED_BOUNDARY_DEVICE_COLUMNS = ("block", "cell", "dir", "q", "fallback", "surface_velocity")
+
+
 @dataclass
 class _BlockedBoundary:
     """One level's stationary BFL data and mutable Jet velocity column."""
@@ -87,10 +92,8 @@ def _prepare_blocked_boundary(boundary, layout, cfg, jets):
         return None
     fluid_idx, dirs, q, body_id, rx, ry, _, fallback = boundary.curved_host
     tables, _ = blocked_link_tables(layout, boundary.flags, fluid_idx, dirs, q, fallback)
-    # Physical row identity excludes body_id: duplicate owners are always wrong.
+    # blocked_link_tables/map_cut_links already rejects repeated (fluid_idx, dir),
-    keys = np.column_stack((fluid_idx, dirs))
+    # independently of body_id, and repeated local writers.
-    if len(np.unique(keys, axis=0)) != len(keys):
-        raise ValueError("duplicate physical boundary cut-link")
     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)

tests/__init__.py

@@ -0,0 +1 @@
+"""CelerisLab test package; keeps shared test helpers independent of discovery paths."""

tests/unit/test_raw_exchange.py

@@ -1,11 +1,9 @@
 """CPU reduction and real host scheduling with CUDA doubles, not GPU physics."""
 from contextlib import nullcontext
 from dataclasses import replace
-from unittest.mock import Mock, patch
+from unittest.mock import Mock
-
 import numpy as np
 import unittest
-
 from CelerisLab.body.geometry.circle import CircleGeometry
 from CelerisLab.body.objects import SimObject, SurfaceActuation
 from CelerisLab.config import LBMConfig
@@ -13,7 +11,7 @@ from CelerisLab.lbm.blocks import BlockLayout2D
 from CelerisLab.lbm.refinement import FixedPatch2D, prepare_fine_boundary
 from CelerisLab.lbm.refinement_tables import _prepare_blocked_boundary
 from CelerisLab.lbm.raw_exchange import summarize_exchange
-from test_refinement_three_level_multipatch import host_module
+from tests.unit.test_refinement_three_level_multipatch import host_module
 
 
 def body(center=(16., 16.), radius=.5, speed=0.):
@@ -21,49 +19,6 @@ def body(center=(16., 16.), radius=.5, speed=0.):
                      surface_actuation=SurfaceActuation(.7, .85, speed))
 
 
-def check_weights_split_and_window(level):
-    h = 2. ** -level
-    rows = np.array([[2, -3, 5, 7], [-1, 4, -2, -6], [99, 99, 99, 99]], np.float32)
-    kw = dict(body_id=0, level=level, h=h, time_end=3, reference_center=(7, 8))
-    last = summarize_exchange(rows, [0, 0, 1], [False, True, False], window='last_substep', **kw)
-    integral = rows.astype(float)*[h*h, h*h, h**3, h*h]/h
-    full = summarize_exchange(integral, [0, 0, 1], [False, True, False], window='last_base_step', **kw)
-    np.testing.assert_array_equal(last['total']['force'], np.array([1, 1])*h)
-    assert last['total']['torque'] == 3*h*h
-    assert last['total']['mass_into_fluid'] == h*h
-    for name in ('force', 'torque', 'mass_rate'):
-        np.testing.assert_array_equal(last['total'][name], full['total'][name])
-    for name in ('force', 'impulse', 'torque', 'angular_impulse', 'mass_into_fluid'):
-        np.testing.assert_array_equal(last['total'][name],
-                                     np.asarray(last['wall'][name])+last['aperture'][name])
-    assert last['time_start'] == 3-h and full['time_start'] == 2
-    assert last['sample_count'] == 1 and full['sample_count'] == 2**level
-    assert full['reference_center'] == (7, 8)
-
-
-def test_final_substep_is_not_interval_mean():
-    kw = dict(body_id=0, level=1, h=.5, time_end=1, reference_center=(0, 0))
-    last = summarize_exchange([[4, 0, 0, 0]], [0], [False], window='last_substep', **kw)
-    full = summarize_exchange([[1.5, 0, 0, 0]], [0], [False], window='last_base_step', **kw)
-    assert last['total']['force'][0] == 2 and full['total']['force'][0] == 1.5
-
-
-def test_host_metadata_zero_jet_and_rows():
-    cfg = LBMConfig(nx=32, ny=32, viscosity=.1)
-    pair = FixedPatch2D.from_config(cfg, (4, 4, 28, 28))
-    boundary = prepare_fine_boundary(cfg, pair, [body(radius=2.15)])
-    local = replace(cfg, nx=pair.fine_shape[1], ny=pair.fine_shape[0], viscosity=.2)
-    local.validate()
-    layout = BlockLayout2D.build(*pair.fine_shape, 8, 8, periodic=(False, False))
-    b = _prepare_blocked_boundary(boundary, layout, local, ())
-    for key, idx in [('body_id', 3), ('rx', 4), ('ry', 5)]:
-        np.testing.assert_array_equal(b.tables[key], boundary.curved_host[idx])
-    assert b.tables['aperture'].any() and (~b.tables['aperture']).any()
-    assert not b.tables['surface_velocity'].any()
-    # Classification does not disappear when prescribed velocity is zero.
-    assert len(b.tables['cell']) == len(boundary.curved_host[0])
-
-
 def backend(observe=True, obj=None):
     host = host_module()
     calls = []
@@ -83,89 +38,133 @@ def backend(observe=True, obj=None):
     return run, host, calls, arrays
 
 
-def test_four_level_schedule_publication_and_fixed_storage():
+class RawExchange(unittest.TestCase):
+    def _check_weights_split_and_window(self, level):
+        h = 2. ** -level
+        rows = np.array([[2, -3, 5, 7], [-1, 4, -2, -6], [99, 99, 99, 99]], np.float32)
+        kw = dict(body_id=0, level=level, h=h, time_end=3, reference_center=(7, 8))
+        last = summarize_exchange(rows, [0, 0, 1], [False, True, False], window='last_substep', **kw)
+        integral = rows.astype(float)*[h*h, h*h, h**3, h*h]/h
+        full = summarize_exchange(integral, [0, 0, 1], [False, True, False], window='last_base_step', **kw)
+        np.testing.assert_array_equal(last['total']['force'], np.array([1, 1])*h)
+        self.assertEqual(last['total']['torque'], 3 * h * h)
+        self.assertEqual(last['total']['mass_into_fluid'], h * h)
+        for name in ('force', 'torque', 'mass_rate'):
+            np.testing.assert_array_equal(last['total'][name], full['total'][name])
+        for name in ('force', 'impulse', 'torque', 'angular_impulse', 'mass_into_fluid'):
+            np.testing.assert_array_equal(last['total'][name],
+                                         np.asarray(last['wall'][name])+last['aperture'][name])
+        self.assertEqual(last['time_start'], 3 - h)
+        self.assertEqual(full['time_start'], 2)
+        self.assertEqual(last['sample_count'], 1)
+        self.assertEqual(full['sample_count'], 2 ** level)
+        self.assertEqual(full['reference_center'], (7, 8))
+
+    def test_final_substep_is_not_interval_mean(self):
+        kw = dict(body_id=0, level=1, h=.5, time_end=1, reference_center=(0, 0))
+        last = summarize_exchange([[4, 0, 0, 0]], [0], [False], window='last_substep', **kw)
+        full = summarize_exchange([[1.5, 0, 0, 0]], [0], [False], window='last_base_step', **kw)
+        self.assertEqual(last['total']['force'][0], 2)
+        self.assertEqual(full['total']['force'][0], 1.5)
+
+    def test_host_metadata_zero_jet_and_rows(self):
+        cfg = LBMConfig(nx=32, ny=32, viscosity=.1)
+        pair = FixedPatch2D.from_config(cfg, (4, 4, 28, 28))
+        boundary = prepare_fine_boundary(cfg, pair, [body(radius=2.15)])
+        local = replace(cfg, nx=pair.fine_shape[1], ny=pair.fine_shape[0], viscosity=.2)
+        local.validate()
+        layout = BlockLayout2D.build(*pair.fine_shape, 8, 8, periodic=(False, False))
+        b = _prepare_blocked_boundary(boundary, layout, local, ())
+        for key, idx in [('body_id', 3), ('rx', 4), ('ry', 5)]:
+            np.testing.assert_array_equal(b.tables[key], boundary.curved_host[idx])
+        self.assertTrue(b.tables['aperture'].any())
+        self.assertTrue((~b.tables['aperture']).any())
+        self.assertFalse(b.tables['surface_velocity'].any())
+        # Classification does not disappear when prescribed velocity is zero.
+        self.assertEqual(len(b.tables['cell']), len(boundary.curved_host[0]))
+
+    def test_four_level_schedule_publication_and_fixed_storage(self):
         run, host, calls, arrays = backend()
         b = run._level_boundaries[3]
         rows = len(b.tables['cell'])
-    assert b.last_raw.size == rows*16 and b.interval_integral.size == rows*32
+        self.assertEqual(b.last_raw.size, rows * 16)
-    with unittest.TestCase().assertRaisesRegex(RuntimeError, 'completed'):
+        self.assertEqual(b.interval_integral.size, rows * 32)
+        with self.assertRaisesRegex(RuntimeError, 'completed'):
             run.read_raw_exchange(0, level=3)
         run.step()
         names = [name for name, _, _ in calls]
-    assert names.count('BlockAccumulateExchange') == 8
+        self.assertEqual(names.count('BlockAccumulateExchange'), 8)
         for i, (name, args, kw) in enumerate(calls):
             if name == 'BlockAccumulateExchange':
-            assert names[i-1] == 'BlockCurvedBoundary' and names[i+1] == 'BlockStep'
+                self.assertEqual(names[i - 1], 'BlockCurvedBoundary')
-            assert args[0] is b.last_raw and args[1] is b.interval_integral
+                self.assertEqual(names[i + 1], 'BlockStep')
-            assert args[3:] == (.125, .125)
+                self.assertIs(args[0], b.last_raw)
-            assert kw['stream'] is run._stream
+                self.assertIs(args[1], b.interval_integral)
-    assert host.cuda.memset_d8_async.call_count == 1
+                self.assertEqual(args[3:], (0.125, 0.125))
+                self.assertIs(kw['stream'], run._stream)
+        self.assertEqual(host.cuda.memset_d8_async.call_count, 1)
         value = run.read_raw_exchange(0, level=3)
-    assert value['sample_count'] == 8 and value['time_start'] == 0
+        self.assertEqual(value['sample_count'], 8)
+        self.assertEqual(value['time_start'], 0)
         run.step(0)
-    assert run.read_raw_exchange(0, level=3)['time_end'] == 1
+        self.assertEqual(run.read_raw_exchange(0, level=3)['time_end'], 1)
-    assert run.read_raw_exchange(0, level=3, window='last_substep')['time_start'] == .875
+        self.assertEqual(run.read_raw_exchange(0, level=3, window='last_substep')['time_start'], 0.875)
         run.step(2)
-    assert run.read_raw_exchange(0, level=3)['time_start'] == 2
+        self.assertEqual(run.read_raw_exchange(0, level=3)['time_start'], 2)
         run.initialize(*arrays)
-    with unittest.TestCase().assertRaisesRegex(RuntimeError, 'completed'):
+        with self.assertRaisesRegex(RuntimeError, 'completed'):
             run.read_raw_exchange(0, level=3)
         run.close()
 
-
+    def test_torque_reference_uses_geometry_center_not_optional_object_center(self):
-
-def test_torque_reference_uses_geometry_center_not_optional_object_center():
         for center in (None, (3., 4.)):
             obj = SimObject(91, geometry=CircleGeometry(16., 16., .5), center=center, radius=.5)
             run, _, _, _ = backend(obj=obj)
             run.step()
             result = run.read_raw_exchange(0, level=3)
-        assert result['reference_center'] == (16., 16.)
+            self.assertEqual(result['reference_center'], (16.0, 16.0))
-        assert result['reference_center'] != tuple(obj.center)
+            self.assertNotEqual(result['reference_center'], tuple(obj.center))
             run.close()
 
-def test_disabled_has_nulls_and_no_observation_allocation_or_launch():
+    def test_disabled_has_nulls_and_no_observation_allocation_or_launch(self):
         run, host, calls, _ = backend(False)
         b = run._level_boundaries[3]
-    assert b.last_raw is b.interval_integral is None
+        self.assertIs(b.last_raw, b.interval_integral)
+        self.assertIs(b.interval_integral, None)
         run.step()
         for name, args, _ in calls:
             if name == 'BlockCurvedBoundary':
-            assert tuple(args[8:]) == (0, 0, 0)
+                self.assertEqual(tuple(args[8:]), (0, 0, 0))
-    assert all(n != 'BlockAccumulateExchange' for n, _, _ in calls)
+        self.assertTrue(all((n != 'BlockAccumulateExchange' for n, _, _ in calls)))
         host.cuda.memset_d8_async.assert_not_called()
-    with unittest.TestCase().assertRaisesRegex(RuntimeError, 'disabled'):
+        with self.assertRaisesRegex(RuntimeError, 'disabled'):
             run.read_raw_exchange(0, level=3)
         run.close()
 
-
+    def test_failed_sync_never_publishes_partial_interval(self):
-def test_failed_sync_never_publishes_partial_interval():
         run, host, calls, _ = backend()
         run.step()
         run._stream.synchronize.side_effect = RuntimeError('device failure')
-    with unittest.TestCase().assertRaisesRegex(RuntimeError, 'device failure'):
+        with self.assertRaisesRegex(RuntimeError, 'device failure'):
             run.step()
-    assert not run._exchange_valid
+        self.assertFalse(run._exchange_valid)
-    with unittest.TestCase().assertRaisesRegex(RuntimeError, 'failed'):
+        with self.assertRaisesRegex(RuntimeError, 'failed'):
             run.read_raw_exchange(0, level=3)
         run._stream.synchronize.side_effect = None
         run.close()
 
-
+    def test_dense_observation_rejected_before_device_access(self):
-def test_dense_observation_rejected_before_device_access():
         host = host_module()
         ctx = Mock()
-    with unittest.TestCase().assertRaisesRegex(ValueError, 'blocked body'):
+        with self.assertRaisesRegex(ValueError, 'blocked body'):
             host.FixedPatchStepper(LBMConfig(nx=32, ny=32, viscosity=.1), (4, 4, 28, 28),
                                    context=ctx, fine_objects=[body(radius=2.15)], observe_raw_exchange=True)
         ctx.activate.assert_not_called()
 
-
+    def test_all_depth_weights(self):
-def test_all_depth_weights():
         for level in (1, 2, 3, 7):
-        check_weights_split_and_window(level)
+            self._check_weights_split_and_window(level)
 
 
-def load_tests(loader, tests, pattern):
+if __name__ == "__main__":
-    return unittest.TestSuite(unittest.FunctionTestCase(value) for name, value in globals().items()
+    unittest.main()
-                              if name.startswith("test_") and callable(value))