colour_bc.hpp

Go to the documentation of this file.

 
#ifndef PECLET_FLOW_VOF_COLOUR_BC_HPP
#define PECLET_FLOW_VOF_COLOUR_BC_HPP
 
#include <Kokkos_Core.hpp>
 
#include "mac_stencils.hpp"  // peclet::flow::SExec, SField, I3, L3
#include "policy.hpp"
#include "vof/advect_wy.hpp"  // UCField
 
namespace peclet::flow::vof {
 
namespace bcdetail {
KOKKOS_INLINE_FUNCTION void axisDims(I3 e, int (&dims)[3], long (&strides)[3]) {
  dims[0] = e.x;
  dims[1] = e.y;
  dims[2] = e.z;
  strides[0] = 1;
  strides[1] = e.x;
  strides[2] = static_cast<long>(e.x) * e.y;
}
}  // namespace bcdetail
 
inline void buildOutsideMask(UCField m, I3 e, int g, I3 o, I3 gs, bool px, bool py, bool pz) {
  const bool p0 = px, p1 = py, p2 = pz;
  const int q0 = gs.x, q1 = gs.y, q2 = gs.z;
  Kokkos::parallel_for(
      "vof::buildOutsideMask", MDRange3<SExec>(SExec(), {0, 0, 0}, {e.x, e.y, e.z}),
      KOKKOS_LAMBDA(int x, int y, int z) {
        const int gx = x - g + o.x, gy = y - g + o.y, gz = z - g + o.z;
        const bool out = (!p0 && (gx < 0 || gx >= q0)) || (!p1 && (gy < 0 || gy >= q1)) ||
                         (!p2 && (gz < 0 || gz >= q2));
        m(L3(x, y, z, e)) = out ? 1u : 0u;
      });
  Kokkos::fence();
}
 
inline void bcColourConst(SField f, I3 e, int g, int a, int s, double value) {
  int dims[3];
  long strides[3];
  bcdetail::axisDims(e, dims, strides);
  const int b = (a + 1) % 3, c = (a + 2) % 3;
  const long sa = strides[a], sb = strides[b], sc = strides[c];
  const int na = dims[a];
  const int lo = (s == 0) ? 0 : (na - g);
  Kokkos::parallel_for(
      "vof::bcColourConst", MDRange2<SExec>(SExec(), {0, 0}, {dims[b], dims[c]}),
      KOKKOS_LAMBDA(int p0, int p1) {
        const long base = static_cast<long>(p0) * sb + static_cast<long>(p1) * sc;
        for (int k = 0; k < g; ++k)
          f(base + static_cast<long>(lo + k) * sa) = value;
      });
  Kokkos::fence();
}
 
inline void bcColourProfile(SField f, I3 e, int g, int a, int s, SField prof, int nc) {
  int dims[3];
  long strides[3];
  bcdetail::axisDims(e, dims, strides);
  const int b = (a + 1) % 3, c = (a + 2) % 3;
  const long sa = strides[a], sb = strides[b], sc = strides[c];
  const int na = dims[a];
  const int lo = (s == 0) ? 0 : (na - g);
  Kokkos::parallel_for(
      "vof::bcColourProfile", MDRange2<SExec>(SExec(), {0, 0}, {dims[b], dims[c]}),
      KOKKOS_LAMBDA(int p0, int p1) {
        const long base = static_cast<long>(p0) * sb + static_cast<long>(p1) * sc;
        const double v = prof(static_cast<long>(p0) * nc + p1);
        for (int k = 0; k < g; ++k)
          f(base + static_cast<long>(lo + k) * sa) = v;
      });
  Kokkos::fence();
}
 
inline void bcColourBackflow(SField f, I3 e, int g, int a, int s, SField un, double backflow) {
  int dims[3];
  long strides[3];
  bcdetail::axisDims(e, dims, strides);
  const int b = (a + 1) % 3, c = (a + 2) % 3;
  const long sa = strides[a], sb = strides[b], sc = strides[c];
  const int na = dims[a];
  const int lo = (s == 0) ? 0 : (na - g);
  const int fa = (s == 0) ? (g - 1) : (na - g - 1);  // a-index of the cell owning the boundary face
  const double sgn = (s == 0) ? 1.0 : -1.0;          // sign of an INWARD normal velocity
  Kokkos::parallel_for(
      "vof::bcColourBackflow", MDRange2<SExec>(SExec(), {0, 0}, {dims[b], dims[c]}),
      KOKKOS_LAMBDA(int p0, int p1) {
        const long base = static_cast<long>(p0) * sb + static_cast<long>(p1) * sc;
        if (sgn * un(base + static_cast<long>(fa) * sa) <= 0.0)
          return;  // leaving (or exactly zero): keep the zero-gradient band
        for (int k = 0; k < g; ++k)
          f(base + static_cast<long>(lo + k) * sa) = backflow;
      });
  Kokkos::fence();
}
 
}  // namespace peclet::flow::vof
 
#endif  // PECLET_FLOW_VOF_COLOUR_BC_HPP