Files
CelerisLab/src/CelerisLab/lbm/kernels/step/one_step_esopull.cu
T
ClaudeandClaude Opus 5.5 2b3c2a7529 kernels: name compilation entries by storage route (uniform.cu, blocked.cu)
kernel_v2.cu -> uniform.cu (whole-domain NX/NY, Simulation) and
blocks.cu -> blocked.cu (one block window, every blocked level), with
InitTubeFlow_v2 -> InitTubeFlow. All live callers, include comments,
CONFIG.md, lifecycle doc and current-tree tests move together; no alias.
refinement.cu keeps its name until the dense-route decision.
Historical-revision tooling that names the old files (LES source pins,
timing/J campaign runners, TG run_case, supported_settings citations at
888c8038, development notes) is unchanged on purpose.
No arithmetic change; the PTX symbol InitTubeFlow_v2 is renamed.

Co-Authored-By: Claude Opus 5.5 <noreply@anthropic.com>
2026-10-26 14:07:34 +00:00

192 lines
6.5 KiB
Plaintext
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
// CelerisLab – step/one_step_esopull.cu
// Esoteric-Pull single-buffer step kernel.
// Included by uniform.cu. No standalone compilation.
//
// CURVED BC: supported via PrepareEsoPullCurvedKernel + ApplyEsoPullCurvedKernel
// in aux_kernels.cu, launched by stepper.py before EsoPullStep.
// ============================================================================
#ifndef CELERIS_STEP_ONE_STEP_ESOPULL_CU
#define CELERIS_STEP_ONE_STEP_ESOPULL_CU
// ---- Helper: apply inlet/outlet BC in esopull mode ----
#if DIM == 2
__device__ __forceinline__ void apply_boundary_esopull(
float* f, uint16_t fl, unsigned int x, unsigned int y,
const float* boundary_values)
{
// Curved-BC nodes have been pre-processed by PrepareEsoPullCurvedKernel
// + ApplyEsoPullCurvedKernel, so is_curved() cells here have valid backed
// populations. Bounce-back on obstacle cells not covered by curved
// (is_obstacle && !is_curved) is handled in EsoPullStep.
if (is_curved(fl)) return;
bool interior_y = (y > 0u) && (y < (unsigned int)(NY - 1));
if (is_inlet(fl) && interior_y) {
apply_inlet_esopull_d2q9(f, x, y, boundary_values + 6u * NX + y * NQ);
}
else if (is_outlet(fl) && interior_y) {
apply_outlet_esopull_d2q9(f, boundary_values + 6u * NX + (NY + y) * NQ);
}
else {
bounce_back_swap(f);
}
}
#endif // DIM == 2
#if DIM == 3
__device__ __forceinline__ void apply_boundary_esopull_3d(
float* f, uint16_t fl, unsigned int x, unsigned int y, unsigned int z,
const fpxx* fi, unsigned long t)
{
if (is_curved(fl)) return;
bool interior_y = (y > 0u) && (y < (unsigned int)(NY - 1));
if (is_inlet(fl) && interior_y) {
apply_inlet_esopull_d3q19(f, x, y, z, fi, t);
}
else if (is_outlet(fl) && interior_y) {
apply_outlet_esopull_d3q19(f, x, y, z, fi, t);
}
else {
bounce_back_swap(f);
}
}
#endif
// One read-only preparation for all D2 planar boundary donors, before curved
// or in-place writes. [6*NX wall values][2*NY*NQ x-boundary donor values].
// t=0 uses the initialized physical state; later t uses spatial postcollision
// values, matching the double-buffer donor reads. Slots retain DDF shifting.
__global__ void PrepareEsoPullBoundaries(const fpxx* fi, float* boundary_values,
unsigned long t)
{
#if DIM == 2
const unsigned int n = threadIdx.x + blockIdx.x * blockDim.x;
const unsigned int side = blockIdx.y;
if (side >= 2u) return;
if (n < NX) {
const unsigned int y = side == 0u ? 1u : NY - 2u;
const unsigned long k = linear_index(n, y);
unsigned long j[NQ]; compute_neighbors(n, y, j);
const unsigned int lower[3] = {4u, 6u, 7u};
const unsigned int upper[3] = {3u, 5u, 8u};
for (unsigned int a = 0; a < 3u; ++a) {
unsigned int donor = side == 0u ? lower[a] : upper[a];
#if Y_WALL_BC == 1
if (a > 0u) donor = side == 0u ? lower[3u-a] : upper[3u-a];
#endif
boundary_values[(side * NX + n) * 3u + a] =
load_physical_dir_pre_streaming(k, donor, fi, j, t);
}
}
if (n < NY) {
const unsigned int x = side == 0u ? 1u : NX - 2u;
const unsigned long k = linear_index(x, n);
unsigned long j[NQ]; compute_neighbors(x, n, j);
float* donor = boundary_values + 6u * NX + (side * NY + n) * NQ;
for (unsigned int i = 0; i < NQ; ++i)
donor[i] = load_physical_dir_pre_streaming(k, i, fi, j, t);
}
#endif
}
// ============================================================================
// EsoPullStep — Esoteric-Pull single-buffer step kernel
// ============================================================================
__global__ __launch_bounds__(NT, 2)
void EsoPullStep(
fpxx* fi, uint16_t* flag,
float* action, float* obs,
const float* boundary_values, unsigned long t)
{
#if DIM == 2
unsigned int x, y; unsigned long k;
index_from_thread(x, y, k);
if (x >= (unsigned int)NX || y >= (unsigned int)NY) return;
uint16_t fl = flag[k];
unsigned long j[NQ];
const bool interior_y = (y > 0u) && (y < (unsigned int)(NY - 1));
compute_neighbors(x, y, j);
float f[NQ];
load_f_esopull(k, f, fi, j, t);
// Pure solid early exit
if (is_solid(fl) && !has_bc(fl)) {
bounce_back_swap(f);
store_f_esopull(k, f, fi, j, t);
return;
}
// Boundary nodes
if (is_solid(fl))
apply_boundary_esopull(f, fl, x, y, boundary_values);
if (is_obstacle(fl) && !is_curved(fl))
bounce_back_swap(f);
// ---- y-wall adjacent row correction ----
if (is_fluid(fl) && (y == 1u || y == (unsigned int)(NY - 2))) {
const unsigned int side = y == 1u ? 0u : 1u;
const float* staged = boundary_values + (side * NX + x) * 3u;
if (side == 0u) { f[3] = staged[0]; f[5] = staged[1]; f[8] = staged[2]; }
else { f[4] = staged[0]; f[6] = staged[1]; f[7] = staged[2]; }
}
// Collision
const bool collide_inlet_ghost = is_inlet(fl) && interior_y
&& inlet_scheme_uses_post_collision_ghost();
if (is_fluid(fl) || collide_inlet_ghost) {
float rho_n, ux, uy;
compute_rho_u(f, rho_n, ux, uy);
collide_dispatch(f, rho_n, ux, uy);
}
store_f_esopull(k, f, fi, j, t);
#elif DIM == 3
unsigned int x, y, z; unsigned long k;
index_from_thread(x, y, z, k);
if (x >= (unsigned int)NX || y >= (unsigned int)NY || z >= (unsigned int)NZ) return;
uint16_t fl = flag[k];
unsigned long j[NQ];
const bool interior_y = (y > 0u) && (y < (unsigned int)(NY - 1));
compute_neighbors(k, j);
float f[NQ];
load_f_esopull(k, f, fi, j, t);
if (is_solid(fl) && !has_bc(fl)) {
bounce_back_swap(f);
store_f_esopull(k, f, fi, j, t);
return;
}
if (is_solid(fl))
apply_boundary_esopull_3d(f, fl, x, y, z, fi, t);
if (is_obstacle(fl) && !is_curved(fl))
bounce_back_swap(f);
if (is_fluid(fl) && (y == 1u || y == (unsigned int)(NY - 2)))
apply_wall_bb_d3q19_y_pull(y, f, fi, k);
const bool collide_inlet_ghost = is_inlet(fl) && interior_y
&& inlet_scheme_uses_post_collision_ghost();
if (is_fluid(fl) || collide_inlet_ghost) {
float rho_n, ux, uy, uz;
compute_rho_u(f, rho_n, ux, uy, uz);
collide_dispatch(f, rho_n, ux, uy, uz);
}
store_f_esopull(k, f, fi, j, t);
#endif
}
#endif // CELERIS_STEP_ONE_STEP_ESOPULL_CU