Compare commits
9
Commits
| Author | SHA1 | Date | |
|---|---|---|---|
|
|
c95a7baae4 | ||
|
|
f4372a84c0 | ||
|
|
e48498d089 | ||
|
|
531739454e | ||
|
|
b93de2ce24 | ||
|
|
a454f02355 | ||
|
|
1cbdc479cf | ||
|
|
1fa498c7d7 | ||
|
|
781333407e |
@@ -207,19 +207,30 @@ def surface_link_velocity(objects, cfg, flags, fluid_idx, dir, body_id, rx, ry,
|
|||||||
jets=(), hits=None, *, _shifted_blocked=False):
|
jets=(), hits=None, *, _shifted_blocked=False):
|
||||||
"""Per-link prescribed surface velocity (n, 2) FP32 and aperture mask.
|
"""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
|
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
|
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.
|
not share a link with each other or with a legacy aperture.
|
||||||
Shifted storage is privately enabled only by the all-blocked stepper.
|
Shifted storage is privately enabled only by the all-blocked stepper.
|
||||||
"""
|
"""
|
||||||
enabled = [obj for obj in objects if obj.surface_actuation is not None]
|
enabled = [obj for obj in objects if obj.surface_actuation is not None]
|
||||||
if not enabled and not jets:
|
if not enabled and not jets:
|
||||||
return None, None
|
return None, None, None, None
|
||||||
from .preprocess.flag_overlay import validate_cut_link_ownership
|
from .preprocess.flag_overlay import validate_cut_link_ownership
|
||||||
validate_cut_link_ownership(flags, fluid_idx, dir, cfg.nq)
|
validate_cut_link_ownership(flags, fluid_idx, dir, cfg.nq)
|
||||||
surface_velocity = np.zeros((len(dir), 2), np.float32)
|
surface_velocity = np.zeros((len(dir), 2), np.float32)
|
||||||
surface_aperture = np.zeros(len(dir), bool)
|
surface_aperture = np.zeros(len(dir), bool)
|
||||||
|
shape = jet_id = None
|
||||||
for obj in enabled:
|
for obj in enabled:
|
||||||
obj.validate_surface_actuation(cfg)
|
obj.validate_surface_actuation(cfg)
|
||||||
selected = body_id == obj.obj_id
|
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")
|
raise ValueError("Jet overlaps a SurfaceActuation aperture")
|
||||||
surface_velocity[driven] = jet_link_velocity(jets, shape, jet_id)[driven]
|
surface_velocity[driven] = jet_link_velocity(jets, shape, jet_id)[driven]
|
||||||
surface_aperture |= driven
|
surface_aperture |= driven
|
||||||
return surface_velocity, surface_aperture
|
return surface_velocity, surface_aperture, shape, jet_id
|
||||||
|
|||||||
@@ -94,12 +94,10 @@ class CurvedLinkSoA:
|
|||||||
"""
|
"""
|
||||||
self.surface_velocity = self.aperture = None
|
self.surface_velocity = self.aperture = None
|
||||||
self.jets, self.jet_shape, self.jet_id = tuple(jets), None, 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,
|
objects, cfg, flags, self.fluid_idx, self.dir, self.body_id, self.rx, self.ry,
|
||||||
jets=self.jets, hits=hits)
|
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):
|
def set_jet_amplitude(self, index, amplitude, stream=None):
|
||||||
"""Rewrite jet rows from cached shapes; re-upload only the velocity column.
|
"""Rewrite jet rows from cached shapes; re-upload only the velocity column.
|
||||||
|
|||||||
@@ -168,21 +168,15 @@ extern "C" __global__ void BlockStep(const fpxx* in, fpxx* out, unsigned int hal
|
|||||||
stream_pull_load(k, f, in + base, j);
|
stream_pull_load(k, f, in + base, j);
|
||||||
if (cell_kind == KIND_ACTIVE) {
|
if (cell_kind == KIND_ACTIVE) {
|
||||||
if (d_block_wall.nx != 0) {
|
if (d_block_wall.nx != 0) {
|
||||||
const unsigned int bx = NX - 2 * halo, by = NY - 2 * halo;
|
unsigned int gx, gy;
|
||||||
const unsigned int nbx = d_block_wall.nx / bx;
|
block_global_xy(x, y, halo, d_block_wall.nx, gx, gy);
|
||||||
const unsigned int gx = (blockIdx.z % nbx) * bx + x - halo;
|
|
||||||
const unsigned int gy = (blockIdx.z / nbx) * by + y - halo;
|
|
||||||
if (gy == 1)
|
if (gy == 1)
|
||||||
reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(
|
reconstruct_bottom_velocity(f, 0.0f, bottom_wall_normal_velocity(
|
||||||
gx, d_block_wall.nx, d_block_wall.slot_width, d_block_wall.amplitude));
|
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);
|
else if (gy == d_block_wall.ny - 2) reconstruct_top_reservoir(f, RHO);
|
||||||
}
|
}
|
||||||
#if COLLISION_MODEL == 2
|
#if COLLISION_MODEL == 2
|
||||||
// Uniform MRT sequence (collide_dispatch, zero forcing); SRT path below unchanged.
|
collide_mrt_unforced(f); // 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
|
#else
|
||||||
float rho, ux, uy, feq[NQ];
|
float rho, ux, uy, feq[NQ];
|
||||||
compute_rho_u(f, rho, ux, uy);
|
compute_rho_u(f, rho, ux, uy);
|
||||||
|
|||||||
@@ -53,6 +53,18 @@ __device__ __forceinline__ void index_from_thread(unsigned int& x,
|
|||||||
k = (unsigned long)y * NX + 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)
|
// Neighbor list (periodic BCs via modular arithmetic)
|
||||||
// ---------------------------------------------------------------------------
|
// ---------------------------------------------------------------------------
|
||||||
|
|||||||
@@ -114,6 +114,16 @@ __device__ __forceinline__ void collide_mrt_no_force(float* __restrict__ g,
|
|||||||
collide_mrt(g, rho, ux, uy, Fin, omega);
|
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
|
#elif LATTICE_MODEL == LATTICE_D3Q19
|
||||||
|
|
||||||
__device__ __forceinline__ void collide_mrt(float* __restrict__ g,
|
__device__ __forceinline__ void collide_mrt(float* __restrict__ g,
|
||||||
|
|||||||
@@ -37,10 +37,7 @@ __device__ unsigned int incoming_counts[NQ * TOTAL_CELLS];
|
|||||||
|
|
||||||
__device__ __forceinline__ void patch_collide(float* f, unsigned long k) {
|
__device__ __forceinline__ void patch_collide(float* f, unsigned long k) {
|
||||||
#if COLLISION_MODEL == 2
|
#if COLLISION_MODEL == 2
|
||||||
float rho, ux, uy;
|
collide_mrt_unforced(f);
|
||||||
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
|
#else
|
||||||
float rho, ux, uy, feq[NQ];
|
float rho, ux, uy, feq[NQ];
|
||||||
compute_rho_u(f, rho, ux, uy);
|
compute_rho_u(f, rho, ux, uy);
|
||||||
|
|||||||
@@ -15,10 +15,8 @@ extern "C" __global__ void BlockChannelStep(
|
|||||||
unsigned int x, y; unsigned long k;
|
unsigned int x, y; unsigned long k;
|
||||||
index_from_thread(x, y, k);
|
index_from_thread(x, y, k);
|
||||||
if (x < halo || y < halo || x >= NX - halo || y >= NY - halo) return;
|
if (x < halo || y < halo || x >= NX - halo || y >= NY - halo) return;
|
||||||
const unsigned int bx = NX - 2 * halo, by = NY - 2 * halo;
|
unsigned int gx, gy;
|
||||||
const unsigned int nbx = d_block_channel.nx / bx;
|
block_global_xy(x, y, halo, d_block_channel.nx, gx, gy);
|
||||||
const unsigned int gx = (blockIdx.z % nbx) * bx + x - halo;
|
|
||||||
const unsigned int gy = (blockIdx.z / nbx) * by + y - halo;
|
|
||||||
const unsigned long base = (unsigned long)blockIdx.z * NQ * TOTAL_CELLS;
|
const unsigned long base = (unsigned long)blockIdx.z * NQ * TOTAL_CELLS;
|
||||||
const unsigned long cell = (unsigned long)blockIdx.z * TOTAL_CELLS + k;
|
const unsigned long cell = (unsigned long)blockIdx.z * TOTAL_CELLS + k;
|
||||||
float f[NQ];
|
float f[NQ];
|
||||||
@@ -55,11 +53,7 @@ extern "C" __global__ void BlockChannelStep(
|
|||||||
}
|
}
|
||||||
if (!wall && !outlet) {
|
if (!wall && !outlet) {
|
||||||
#if COLLISION_MODEL == 2
|
#if COLLISION_MODEL == 2
|
||||||
// Match uniform MRT arithmetic (collide_dispatch, zero forcing).
|
collide_mrt_unforced(f);
|
||||||
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
|
#else
|
||||||
// Match uniform SRT arithmetic, including its shifted equilibrium path.
|
// Match uniform SRT arithmetic, including its shifted equilibrium path.
|
||||||
float rho, ux, uy, feq[NQ], forcing[NQ] = {};
|
float rho, ux, uy, feq[NQ], forcing[NQ] = {};
|
||||||
|
|||||||
@@ -8,7 +8,7 @@ from dataclasses import dataclass, replace
|
|||||||
import numpy as np
|
import numpy as np
|
||||||
|
|
||||||
from .blocks import blocked_link_tables
|
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 .descriptors import D2Q9_EX, D2Q9_EY
|
||||||
from .refinement import ACTIVE, OUTER_GHOST
|
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),
|
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))
|
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)
|
actuated = bool(jets) or any(o.surface_actuation is not None for o in boundary.objects)
|
||||||
|
jet_cache = None
|
||||||
if actuated:
|
if actuated:
|
||||||
cfg = replace(cfg)
|
cfg = replace(cfg)
|
||||||
cfg.validate()
|
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,
|
boundary.objects, cfg, boundary.flags, fluid_idx, dirs, body_id, rx, ry,
|
||||||
jets=jets, hits=boundary.curved_hits, _shifted_blocked=_shifted_blocked)
|
jets=jets, hits=boundary.curved_hits, _shifted_blocked=_shifted_blocked)
|
||||||
return _BlockedBoundary(boundary, tables, jets, actuated,
|
if jets:
|
||||||
jet_link_shapes(jets, boundary.curved_hits) if jets else None)
|
jet_cache = (shape, jet_id)
|
||||||
|
return _BlockedBoundary(boundary, tables, jets, actuated, jet_cache)
|
||||||
|
|||||||
@@ -1,6 +1,7 @@
|
|||||||
"""Host contract of PeriodicBlockStepper.step(n): one stream wait per call (mocked driver)."""
|
"""Host contract of PeriodicBlockStepper.step(n): one stream wait per call (mocked driver)."""
|
||||||
import ast
|
import ast
|
||||||
from contextlib import nullcontext
|
from contextlib import nullcontext
|
||||||
|
import importlib.util
|
||||||
from pathlib import Path
|
from pathlib import Path
|
||||||
import types
|
import types
|
||||||
from unittest.mock import Mock
|
from unittest.mock import Mock
|
||||||
@@ -8,14 +9,15 @@ import unittest
|
|||||||
|
|
||||||
import numpy as np
|
import numpy as np
|
||||||
|
|
||||||
import CelerisLab.lbm.block_stepper as b
|
|
||||||
from CelerisLab.config import LBMConfig
|
from CelerisLab.config import LBMConfig
|
||||||
from CelerisLab.lbm.blocks import BlockLayout2D
|
from CelerisLab.lbm.blocks import BlockLayout2D
|
||||||
|
|
||||||
|
|
||||||
def host_module():
|
def host_module():
|
||||||
"""block_stepper with pycuda/compiler replaced by mocks; no GPU context."""
|
"""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 = ast.parse(path.read_text())
|
||||||
tree.body = [node for node in tree.body if not (
|
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)
|
isinstance(node, ast.Import) and any(a.name == "pycuda.driver" for a in node.names)
|
||||||
|
|||||||
@@ -37,8 +37,13 @@ class RefinementMRTHost(unittest.TestCase):
|
|||||||
self.assertIn("STORE_PRECISION != 0", source)
|
self.assertIn("STORE_PRECISION != 0", source)
|
||||||
self.assertIn('#if COLLISION_MODEL == 2\n#include "operators/collision_mrt.cuh"\n#endif', 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
|
self.assertIn("USE_DDF_SHIFTING ||", refinement) # dense/mixed shifted stays rejected
|
||||||
|
# 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, channel, refinement):
|
for source in (blocks, channel, refinement):
|
||||||
self.assertEqual(source.count("collide_mrt_no_force("), 1)
|
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)
|
||||||
|
|
||||||
|
|
||||||
if __name__ == "__main__": unittest.main()
|
if __name__ == "__main__": unittest.main()
|
||||||
|
|||||||
@@ -249,19 +249,35 @@ class ProjectedPairsTests(unittest.TestCase):
|
|||||||
self.assertFalse(assess(incomplete)["comparisons"][3]["passed"])
|
self.assertFalse(assess(incomplete)["comparisons"][3]["passed"])
|
||||||
|
|
||||||
def test_pinned_build_sources_match_host_generated_inputs(self):
|
def test_pinned_build_sources_match_host_generated_inputs(self):
|
||||||
|
"""The pinned fingerprints are what today's compiler generates from the qualified
|
||||||
|
NUMERICAL_REVISION kernels; the current kernels may legitimately differ."""
|
||||||
|
import tarfile
|
||||||
from CelerisLab.cuda import compiler_v2
|
from CelerisLab.cuda import compiler_v2
|
||||||
|
|
||||||
for size in (128, 256):
|
kernels = "src/CelerisLab/lbm/kernels"
|
||||||
cfg = replace(run_case.config(size, 100), nx=10, ny=10)
|
if subprocess.run(["git", "-C", str(REPO), "cat-file", "-e",
|
||||||
cfg.validate()
|
f"{NUMERICAL_REVISION}^{{commit}}"], capture_output=True).returncode:
|
||||||
root = compiler_v2.generate_config(cfg)
|
self.fail(f"{NUMERICAL_REVISION} is not in this clone; the pinned build "
|
||||||
try:
|
"sources cannot be checked")
|
||||||
actual = json.loads((Path(root) / ".celeris-build").read_text())
|
archive = subprocess.run(["git", "-C", str(REPO), "archive", "--format=tar",
|
||||||
expected = build_metadata(size)
|
NUMERICAL_REVISION, kernels],
|
||||||
self.assertEqual(actual, {key: expected[key] for key in
|
capture_output=True, check=True).stdout
|
||||||
("source_sha256", "config", "n_objects")})
|
with tempfile.TemporaryDirectory() as snapshot:
|
||||||
finally:
|
with tarfile.open(fileobj=io.BytesIO(archive)) as tar:
|
||||||
compiler_v2.cleanup_build_dir(root)
|
tar.extractall(snapshot)
|
||||||
|
with patch.object(compiler_v2, "_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_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)
|
||||||
|
|
||||||
def test_source_preflight_requires_clean_exact_revision_and_tree(self):
|
def test_source_preflight_requires_clean_exact_revision_and_tree(self):
|
||||||
import oracle
|
import oracle
|
||||||
|
|||||||
Reference in New Issue
Block a user