diff --git a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/exchange.rs b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/exchange.rs index 7a0f6b5..c0c2ee3 100644 --- a/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/exchange.rs +++ b/crates/specialized/rtx-cfd/src/solvers/incompressible/embedded3/exchange.rs @@ -72,6 +72,7 @@ impl Mask { } let cv = self.cv_geometry(c, p); let u0 = vals[c][idx]; + let shift0 = self.face_shift(c, p); let cell_minus = add(p, ec, -1); let cell_plus = p; for d in 0..3 { @@ -99,6 +100,17 @@ impl Mask { 0.5 * (flux(cell_minus) + flux(cell_plus)) * a_d, ) }; + // The spacing the predictor uses toward a solid + // neighbour (its shift is zero): the centroid + // spacing in a cross direction (S2-5), else the + // exchange distance. + let solid_spacing = |sign: f64| -> f64 { + if self.diffusion_centroid && d != c && !self.wall_exchange_axis { + (h[d] - sign * shift0[d]).clamp(0.25 * h[d], 2.0 * h[d]) + } else { + self.exchange_delta(&cv, d) + } + }; // Plus side: a prescribed neighbour face. if let Some(fp) = lat.face(c, add(p, ed, 1)) { if kind(c, fp) != FaceKind::Fluid { @@ -112,8 +124,7 @@ impl Mask { }; let u_face = upwind(m_plus, u0, un) + delta; let on_fluid = -rho * m_plus * (u_face - u0) - + mu * cv.ap[d][1] * a_d * (un - u0) - / self.exchange_delta(&cv, d); + + mu * cv.ap[d][1] * a_d * (un - u0) / solid_spacing(1.0); force[c] -= on_fluid; } } @@ -130,8 +141,7 @@ impl Mask { }; let u_face = upwind(m_minus, ud, u0) + delta; let on_fluid = rho * m_minus * (u_face - u0) - + mu * cv.ap[d][0] * a_d * (ud - u0) - / self.exchange_delta(&cv, d); + + mu * cv.ap[d][0] * a_d * (ud - u0) / solid_spacing(-1.0); force[c] -= on_fluid; } }