From f6b34f924da1e41ce2b7d9dca05752c24831d89a Mon Sep 17 00:00:00 2001 From: Omar Sobh Date: Sun, 20 Sep 2026 08:59:30 -0500 Subject: [PATCH] =?UTF-8?q?embedded3=20S2-7b:=20the=20moving-wall=20forms?= =?UTF-8?q?=20on=20the=20device=20=E2=80=94=20the=20ub=20table=20at=20the?= =?UTF-8?q?=20open=20part's=20centroid=20foot=20(RTX=5FE3=5FWALL=5FFOOT=3D?= =?UTF-8?q?centroid,=20table=20only),=20the=20solid=20exchange's=20wall=20?= =?UTF-8?q?velocity=20at=20the=20axis=20feet=20(RTX=5FE3=5FWALL=5FEXCHANGE?= =?UTF-8?q?=3Daxisfoot:=20per-face=20six-entry=20foot=20table,=20bit=2012?= =?UTF-8?q?=20in=20e3=5Fcut.cu),=20the=20true-normal=20wall=20flux=20throu?= =?UTF-8?q?gh=20the=20shared=20host=20table=20(RTX=5FE3=5FWALL=5FFLUX=3Dtr?= =?UTF-8?q?ue);=20refusals=20lifted=20(conv=20sides=20stay=20host-only);?= =?UTF-8?q?=20host=20=3D=20device=20green=20on=20the=20sphere=20(varying?= =?UTF-8?q?=20surface=20velocity),=20the=20cylinder=20and=20the=20moving?= =?UTF-8?q?=20circle=20with=20each=20knob=20and=20all=20three?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-Authored-By: Claude Fable 5.1 --- .../rtx-cfd/src/kernels/cuda/e3_cut.cu | 8 ++- .../incompressible/embedded3/step/device.rs | 12 ++-- .../embedded3/step/device/cut.rs | 59 ++++++++++++++++++- 3 files changed, 68 insertions(+), 11 deletions(-) diff --git a/crates/specialized/rtx-cfd/src/kernels/cuda/e3_cut.cu b/crates/specialized/rtx-cfd/src/kernels/cuda/e3_cut.cu index ac61d02..3c0a21a 100644 --- a/crates/specialized/rtx-cfd/src/kernels/cuda/e3_cut.cu +++ b/crates/specialized/rtx-cfd/src/kernels/cuda/e3_cut.cu @@ -30,6 +30,7 @@ struct E3Cut { double *cell_flux; /* scratch per cell */ const double *s_u, *s_v, *s_w; /* open-part centroid shifts per face, 3 interleaved (S2-5) */ const double *vn_u, *vn_v, *vn_w; /* the wall's normal velocity into the fluid per face (A3-i) */ + const double *foot_u, *foot_v, *foot_w; /* the wall velocity at the axis feet per face, 6 interleaved (S2-7b: [2d + side]) */ }; /* f(xi) = xi / (1 - exp(-xi)), f(0) = 1 exactly (closure.rs advancing_factor). */ @@ -206,6 +207,9 @@ __device__ double cut_face_update(const E3Params& g, const E3Ptrs& f, const E3Cu int cp[3] = { i, j, k }; /* cell plus */ const double* ubt = c == 0 ? m.ub_u : (c == 1 ? m.ub_v : m.ub_w); double ub = fast ? 0.0 : ubt[fidx]; + /* bit 12 of wall_order: the solid exchange takes the wall velocity at the axis foot (S2-7b) */ + const double* ft = c == 0 ? m.foot_u : (c == 1 ? m.foot_v : m.foot_w); + const int axis_foot = (g.wall_order & 4096) != 0; /* face position */ double x[3]; #pragma unroll @@ -318,7 +322,7 @@ __device__ double cut_face_update(const E3Params& g, const E3Ptrs& f, const E3Cu diff += mu * g_minus * a_d * cut_transverse(g, m, c, d, fidx, f_dn1, qn, -1.0, u0, ub, alpha, wall, distance, dn1, sh, ubt) / dl; } } - if (solid_up) { double kx = mu * g_plus * a_d / delta_x; wall_implicit += kx; wall_rhs += kx * up1; } + if (solid_up) { double kx = mu * g_plus * a_d / delta_x; wall_implicit += kx; wall_rhs += kx * (axis_foot ? ft[6 * fidx + 2 * d + 1] : up1); } /* fast: dl == h exactly, kx == +0.0 exactly, so the else branch below is (X + 0.0) * (up1 - u0) == X * (up1 - u0) */ else if (fast && f_up1 >= 0 && centroid) diff += (mu * g_plus * a_d / h[d]) * (up1 - u0); else if (f_up1 >= 0 && centroid) { @@ -330,7 +334,7 @@ __device__ double cut_face_update(const E3Params& g, const E3Ptrs& f, const E3Cu } else if (f_up1 >= 0) diff += mu * g_plus * a_d * (up1 - u0) / h[d]; else if (sides[d][1] == SIDE_VELOCITY) diff += mu * g_plus * a_d * (beyond_p - u0) / (0.5 * h[d]); - if (solid_dn) { double kx = mu * g_minus * a_d / delta_x; wall_implicit += kx; wall_rhs += kx * dn1; } + if (solid_dn) { double kx = mu * g_minus * a_d / delta_x; wall_implicit += kx; wall_rhs += kx * (axis_foot ? ft[6 * fidx + 2 * d] : dn1); } else if (fast && f_dn1 >= 0 && centroid) diff -= (mu * g_minus * a_d / h[d]) * (u0 - dn1); else if (f_dn1 >= 0 && centroid) { double dl = h[d] - (sh[3 * f_dn1 + d] - sh[3 * fidx + d]); diff --git a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device.rs b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device.rs index 492fa0b..2beb8f2 100644 --- a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device.rs +++ b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device.rs @@ -178,11 +178,8 @@ impl DeviceStep { "the order-2 second point at the neighbour's centroid (S2-7) is a host prototype: the device kernels do not carry it" ); assert!( - !solver.params.wall_exchange_foot - && !solver.params.conv_sides_exact - && !solver.params.wall_flux_true_normal - && !solver.params.wall_foot_centroid, - "the axis-foot exchange / exact convective sides (S2-7b) are host prototypes: the device kernels do not carry them" + !solver.params.conv_sides_exact, + "the exact convective sides (S2-7b) are a host prototype: the device kernels do not carry it" ); assert!( !solver.params.cv_sides_exact, @@ -329,7 +326,10 @@ impl DeviceStep { // bit 11: the predictor's interior fast path (PERF-3 P3-1; bit-identical // by construction and gated so; default ON, `RTX_E3_PREDICT_FAST=0` restores // the uniform path). - + 2048 * i32::from(!std::env::var("RTX_E3_PREDICT_FAST").is_ok_and(|v| v == "0")), + + 2048 * i32::from(!std::env::var("RTX_E3_PREDICT_FAST").is_ok_and(|v| v == "0")) + // bit 12: the solid exchange's wall velocity at the axis foot (S2-7b + // `axisfoot`; the centroid foot needs no bit — it is in the ub table). + + 4096 * i32::from(self.solver.params.wall_exchange_foot), dx: g.dx, dy: g.dy, dz: g.dz, diff --git a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device/cut.rs b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device/cut.rs index b8c04b4..b64c49f 100644 --- a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device/cut.rs +++ b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/step/device/cut.rs @@ -52,11 +52,11 @@ fn cut_kernels() -> &'static CutKernels { }) } -/// `struct E3Cut` in e3_cut.cu: 21 device pointers. +/// `struct E3Cut` in e3_cut.cu: 27 device pointers. #[repr(C)] #[derive(Clone, Copy)] struct E3CutPtrs { - ptrs: [u64; 24], + ptrs: [u64; 27], } unsafe impl DeviceRepr for E3CutPtrs {} unsafe impl ValidAsZeroBits for E3CutPtrs {} @@ -79,6 +79,9 @@ pub(super) struct DeviceCut { shift: [CudaSlice; 3], /// The wall's normal velocity into the fluid per face (A3-i; one dummy entry when off). vn: [CudaSlice; 3], + /// The wall velocity at the two axis feet along each direction per face, + /// six interleaved (S2-7b `axisfoot`; one dummy entry when off). + foot: [CudaSlice; 3], pub(super) merged: usize, } @@ -136,6 +139,12 @@ impl DeviceCut { let mut open: [Vec; 3] = [Vec::new(), Vec::new(), Vec::new()]; let advancing = mask.wall_advancing; let mut vn: [Vec; 3] = [vec![0.0], vec![0.0], vec![0.0]]; + // S2-7b: the wall velocity at the open part's centroid's foot (the + // host predictor's `foot_of`), and the axis feet of the solid exchange. + let centroid_foot = mask.wall_foot_centroid; + let shifts = mask.face_shift_tables().cloned(); + let axis_foot = mask.wall_exchange_foot; + let mut foot: [Vec; 3] = [vec![0.0], vec![0.0], vec![0.0]]; for c in 0..3 { let (ni, nj, nk) = match c { 0 => (nx + 1, ny, nz), @@ -161,7 +170,15 @@ impl DeviceCut { *ub_out = match shared { Some(sh) => sh[c][idx], None if dists[c][idx].abs() <= band => { - mask.surface_velocity_at(body, x, c, t) + let xf = match (&shifts, centroid_foot) { + (Some(sh), true) => [ + x[0] + sh[c][3 * idx], + x[1] + sh[c][3 * idx + 1], + x[2] + sh[c][3 * idx + 2], + ], + _ => x, + }; + mask.surface_velocity_at(body, xf, c, t) } None => 0.0, }; @@ -182,6 +199,38 @@ impl DeviceCut { }); ub[c] = ubc; open[c] = opc; + if axis_foot { + let lat_face_solid = |cc: usize, q: [i64; 3]| mask.aperture(cc, q) == Some(0.0); + foot[c] = (0..counts[c]) + .into_par_iter() + .flat_map_iter(|idx| { + let mut out = [0.0f64; 6]; + if dists[c][idx].abs() <= band { + let (k, j, i) = (idx / (nj * ni), (idx / ni) % nj, idx % ni); + let p = [i as i64, j as i64, k as i64]; + let x = [ + (i as f64 + if c == 0 { 0.0 } else { 0.5 }) * h[0], + (j as f64 + if c == 1 { 0.0 } else { 0.5 }) * h[1], + (k as f64 + if c == 2 { 0.0 } else { 0.5 }) * h[2], + ]; + let cv = mask.cv_geometry(c, p); + for d in 0..3 { + for (side, sign) in [(0usize, -1.0f64), (1, 1.0)] { + let mut q = p; + q[d] += sign as i64; + if lat_face_solid(c, q) { + let delta = mask.exchange_delta(&cv, d); + let mut xf = x; + xf[d] += sign * delta; + out[2 * d + side] = mask.surface_velocity_at(body, xf, c, t); + } + } + } + } + out + }) + .collect(); + } // A3-i: the wall's normal velocity into the fluid per face within the band. if advancing { vn[c] = (0..counts[c]) @@ -274,6 +323,7 @@ impl DeviceCut { None => [up_f(&[0.0]), up_f(&[0.0]), up_f(&[0.0])], }, vn: [up_f(&vn[0]), up_f(&vn[1]), up_f(&vn[2])], + foot: [up_f(&foot[0]), up_f(&foot[1]), up_f(&foot[2])], merged, }; if profile { @@ -366,6 +416,9 @@ impl DeviceCut { pf(&self.vn[0]), pf(&self.vn[1]), pf(&self.vn[2]), + pf(&self.foot[0]), + pf(&self.foot[1]), + pf(&self.foot[2]), ], } }