ClaudeandClaude Opus 5.5
64d54ab23e
surface_link_velocity (the shared host check used by CurvedLinkSoA, ObjectManager sync, prepare_fine_boundary and FixedPatchStepper) now rejects a jet whose selected cut links carry more than one body id, at any amplitude, before uploads. Identity is the scene-wide body id of the physical cut-link list, so one body spanning storage blocks stays valid. No body selector; empty-selection semantics unchanged. Tests: two bodies under one jet rejected (amplitude 2e-3 and 0) on the Simulation/manager path and on the fine path; separate jets on separate bodies accepted on both; one fine body crossing 2x2 block seams accepted. Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
178 lines
9.1 KiB
Python
178 lines
9.1 KiB
Python
"""Jets on the FixedPatchStepper fine level: CPU wiring checks plus GPU gates.
|
|
|
|
Jets use the same coarse cell-edge frame as fine_objects. The equivalent jet
|
|
of Astra's actuated fixture is SurfaceActuation(0.7, 0.85, U).as_jet((16, 12),
|
|
2.15). Frozen GPU gates (8 coarse steps, both coarse boundary modes):
|
|
G1 jet run: dense vs 1x1/2x2/4x2 fine blocks bit-identical (values + NaN masks);
|
|
G2 jet vs legacy actuation (dense): max |dG| over finite fine ACTIVE fluid
|
|
populations <= 1e-6 (velocity rows differ only by FP32 rounding; reported);
|
|
G3 amplitude change after 4 steps (+U -> -U), 4 more steps: dense vs blocked
|
|
bit-identical and the run differs from the unchanged run.
|
|
"""
|
|
from contextlib import nullcontext
|
|
import unittest
|
|
from unittest.mock import Mock, patch
|
|
|
|
import numpy as np
|
|
import pycuda.driver as cuda
|
|
|
|
from CelerisLab.cuda import compiler_v2 as compiler
|
|
from CelerisLab.lbm.descriptors import SOLID
|
|
from CelerisLab.lbm.refinement import ACTIVE, FixedPatch2D
|
|
from CelerisLab.lbm.refinement_stepper import FixedPatchStepper
|
|
from test_refinement_blocked_fine import (BOUNDS, FINE_BLOCKS, NX, NY, STEPS, actuation,
|
|
circle, config, fine_circle, raw, state)
|
|
|
|
SPEED = .002
|
|
|
|
|
|
def jet(speed=SPEED):
|
|
return actuation(speed).as_jet((16.0, 12.0), 2.15)
|
|
|
|
|
|
def mocked(objects, jets, fine_block):
|
|
"""Real constructor with driver doubles; returns backend and the memcpy mock."""
|
|
context = Mock(activate=nullcontext, sm_arch="sm_70")
|
|
module = Mock(get_function=Mock(return_value=Mock()))
|
|
module.get_global.return_value = (Mock(), 20)
|
|
htod = Mock()
|
|
with patch.object(cuda, "Stream"), patch.object(compiler, "generate_config"), \
|
|
patch.object(compiler, "compile_kernel"), patch.object(compiler, "cleanup_build_dir"), \
|
|
patch.object(compiler, "load_module", return_value=module), \
|
|
patch.object(cuda, "mem_alloc", side_effect=lambda size: Mock()), \
|
|
patch.object(cuda, "memcpy_htod", htod):
|
|
backend = FixedPatchStepper(config("periodic"), BOUNDS, context=context,
|
|
fine_objects=objects, fine_block=fine_block, jets=jets)
|
|
return backend, htod
|
|
|
|
|
|
class JetStepperCPU(unittest.TestCase):
|
|
def test_one_body_across_blocks_accepted(self):
|
|
# Fine circle (19.5, 15.5) R4.3 crosses the 2x2 block seams; one body, one jet.
|
|
import math
|
|
from CelerisLab.body.jets import Jet, Radial, Sector
|
|
C = (16.0, 12.0)
|
|
whole = Jet(Sector(C, .1, math.pi, 3.), Radial(C, inward=False), amplitude=0.)
|
|
backend, _ = mocked([circle()], [whole], FINE_BLOCKS[1])
|
|
_, jet_id = backend._jet_cache
|
|
self.assertGreater(len(np.unique(backend._fine_tables["block"][jet_id == 0])), 1)
|
|
backend.close()
|
|
|
|
def test_blocked_jet_velocity_matches_legacy_rows(self):
|
|
legacy, _ = mocked([circle(actuation(SPEED))], (), FINE_BLOCKS[1])
|
|
driven, _ = mocked([circle()], [jet()], FINE_BLOCKS[1])
|
|
a = legacy._fine_tables["surface_velocity"]
|
|
b = driven._fine_tables["surface_velocity"]
|
|
diff = np.abs(a - b).max()
|
|
print(f"fine jet vs legacy velocity rows: max |dU| = {diff:.3e} (exact={diff == 0})")
|
|
self.assertLessEqual(diff, 4 * np.finfo(np.float32).eps * SPEED)
|
|
legacy.close(), driven.close()
|
|
|
|
def test_blocked_amplitude_update_rewrites_only_velocity_column(self):
|
|
backend, htod = mocked([circle()], [jet()], FINE_BLOCKS[1])
|
|
before = backend._fine_tables["surface_velocity"].copy()
|
|
device = backend._link_tables["surface_velocity"]
|
|
with patch.object(cuda, "memcpy_htod", htod):
|
|
htod.reset_mock()
|
|
backend.set_jet_amplitude(0, -SPEED)
|
|
(ptr, sent), _ = htod.call_args
|
|
self.assertEqual(htod.call_count, 1)
|
|
self.assertIs(ptr, device)
|
|
np.testing.assert_allclose(sent, -before, rtol=1e-6, atol=0)
|
|
np.testing.assert_array_equal(backend._fine_tables["surface_velocity"], sent)
|
|
with self.assertRaises(ValueError):
|
|
backend.set_jet_amplitude(0, 1.5)
|
|
self.assertEqual(htod.call_count, 1)
|
|
self.assertEqual(backend._jets[0].amplitude, -SPEED)
|
|
backend.close()
|
|
|
|
def test_copy_failure_marks_failed_and_keeps_host_amplitude(self):
|
|
for fine_block in (None, FINE_BLOCKS[1]): # dense SoA and blocked tables
|
|
backend, _ = mocked([circle()], [jet()], fine_block)
|
|
status = backend._status
|
|
for bad, error in ((1.5, ValueError), (float("nan"), ValueError)):
|
|
with self.assertRaises(error):
|
|
backend.set_jet_amplitude(0, bad)
|
|
self.assertEqual(backend._status, status) # invalid input stays usable
|
|
with self.assertRaises(IndexError):
|
|
backend.set_jet_amplitude(1, -SPEED)
|
|
self.assertEqual(backend._status, status)
|
|
dense = backend._fine_curved
|
|
before = (dense.surface_velocity if dense is not None
|
|
else backend._fine_tables["surface_velocity"]).copy()
|
|
boom = Mock(side_effect=RuntimeError("copy failed"))
|
|
with patch.object(cuda, "memcpy_htod", boom), patch.object(cuda, "memcpy_htod_async", boom), \
|
|
patch.object(cuda, "pagelocked_empty", side_effect=lambda shape, dtype: np.empty(shape, dtype)):
|
|
with self.assertRaisesRegex(RuntimeError, "copy failed"):
|
|
backend.set_jet_amplitude(0, -SPEED)
|
|
self.assertEqual(backend._status, "failed", fine_block)
|
|
self.assertEqual(backend._jets[0].amplitude, SPEED)
|
|
if dense is not None:
|
|
self.assertEqual(dense.jets[0].amplitude, SPEED)
|
|
np.testing.assert_array_equal(dense.surface_velocity, before)
|
|
else:
|
|
np.testing.assert_array_equal(backend._fine_tables["surface_velocity"], before)
|
|
with self.assertRaisesRegex(RuntimeError, "failed"):
|
|
backend.step(1)
|
|
with self.assertRaisesRegex(RuntimeError, "failed"):
|
|
backend.set_jet_amplitude(0, -SPEED)
|
|
backend.close()
|
|
|
|
def test_jet_overlapping_legacy_aperture_rejected_before_cuda(self):
|
|
with patch.object(cuda, "Stream", side_effect=AssertionError("CUDA touched")):
|
|
with self.assertRaisesRegex(ValueError, "SurfaceActuation"):
|
|
FixedPatchStepper(config("periodic"), BOUNDS, context=Mock(),
|
|
fine_objects=[circle(actuation(SPEED))], jets=[jet()])
|
|
|
|
|
|
def run(mode, fine_block, *, objects, jets, switch=None):
|
|
from CelerisLab.cuda.context import CudaContext
|
|
cfg = config(mode)
|
|
patch_ = FixedPatch2D.from_config(cfg, BOUNDS)
|
|
with CudaContext(0) as context, FixedPatchStepper(
|
|
cfg, BOUNDS, context=context, coarse_boundary=mode, fine_objects=objects,
|
|
fine_block=fine_block, jets=jets) as backend:
|
|
backend.initialize(state(NY, NX, 1), state(*patch_.fine_shape, 2))
|
|
if switch is None:
|
|
backend.step(STEPS)
|
|
else:
|
|
backend.step(STEPS // 2)
|
|
backend.set_jet_amplitude(0, switch)
|
|
backend.step(STEPS - STEPS // 2)
|
|
return backend.download(), raw(backend)
|
|
|
|
|
|
class JetStepperGPU(unittest.TestCase):
|
|
def test_g1_g2_jet_layouts_identical_and_close_to_legacy(self):
|
|
for mode in ("periodic", "slot_reservoir"):
|
|
patch_, boundary = fine_circle(mode)
|
|
fluid = (patch_.fine_kind == ACTIVE) & (
|
|
(boundary.flags.reshape(patch_.fine_shape) & SOLID) == 0)
|
|
dense = run(mode, None, objects=[circle()], jets=[jet()])
|
|
self.assertTrue(np.isfinite(dense[1][1][:, fluid]).all(), mode)
|
|
for fine_block in FINE_BLOCKS:
|
|
public, raw_arrays = run(mode, fine_block, objects=[circle()], jets=[jet()])
|
|
for a, b in zip(public + tuple(raw_arrays), dense[0] + tuple(dense[1])):
|
|
np.testing.assert_array_equal(a, b, err_msg=f"G1 {mode} {fine_block}")
|
|
legacy = run(mode, None, objects=[circle(actuation(SPEED))], jets=())
|
|
diff = np.abs(dense[1][1][:, fluid] - legacy[1][1][:, fluid]).max()
|
|
print(f"G2 {mode}: jet vs legacy fine max |dG| = {diff:.3e}", flush=True)
|
|
self.assertLessEqual(diff, 1e-6, mode)
|
|
|
|
def test_g3_amplitude_switch_layouts_identical(self):
|
|
for mode in ("periodic", "slot_reservoir"):
|
|
dense = run(mode, None, objects=[circle()], jets=[jet()], switch=-SPEED)
|
|
steady = run(mode, None, objects=[circle()], jets=[jet()])
|
|
self.assertFalse(np.array_equal(dense[1][1], steady[1][1], equal_nan=True), mode)
|
|
for fine_block in FINE_BLOCKS:
|
|
public, raw_arrays = run(mode, fine_block, objects=[circle()], jets=[jet()],
|
|
switch=-SPEED)
|
|
for a, b in zip(public + tuple(raw_arrays), dense[0] + tuple(dense[1])):
|
|
np.testing.assert_array_equal(a, b, err_msg=f"G3 {mode} {fine_block}")
|
|
print(f"G3 {mode}: amplitude switch dense vs {len(FINE_BLOCKS)} layouts identical",
|
|
flush=True)
|
|
|
|
|
|
if __name__ == "__main__":
|
|
unittest.main()
|