gauge_exact_gradient.hpp

Go to the documentation of this file.

 
#ifndef PECLET_FLOW_GAUGE_EXACT_GRADIENT_HPP
#define PECLET_FLOW_GAUGE_EXACT_GRADIENT_HPP
 
#include <Kokkos_Core.hpp>
 
#include "mac_cutcell.hpp"  // CCField/CCConst, C3, CCExec
#include "policy.hpp"
 
namespace peclet::flow {
 
inline void gpCenterGrad(CCField out, CCConst p, CCConst sdf, int axis, C3 e, int g,
                         double wa = 1.0) {
  CCExec space;
  using MD = MDRange3<CCExec>;
  Kokkos::parallel_for(
      "peclet::flow::gp_center_grad", MD(space, {g, g, g}, {e.x - g, e.y - g, e.z - g}),
      KOKKOS_LAMBDA(int x, int y, int z) {
        const long sy = e.x, sz = (long)e.x * e.y;
        const long i = (long)x + (long)y * sy + (long)z * sz;
        const long sa = (axis == 0) ? 1 : (axis == 1) ? sy : sz;
        if (sdf(i) < 0.0) {
          out(i) = 0.0;
          return;
        }
        const bool am = sdf(i - sa) >= 0.0, ap = sdf(i + sa) >= 0.0;
        double gr;
        if (am && ap)
          gr = 0.5 * (p(i + sa) - p(i - sa));
        else if (ap)
          gr = (sdf(i + 2 * sa) >= 0.0) ? 0.5 * (-3.0 * p(i) + 4.0 * p(i + sa) - p(i + 2 * sa))
                                        : (p(i + sa) - p(i));
        else if (am)
          gr = (sdf(i - 2 * sa) >= 0.0) ? 0.5 * (3.0 * p(i) - 4.0 * p(i - sa) + p(i - 2 * sa))
                                        : (p(i) - p(i - sa));
        else
          gr = 0.0;
        out(i) = wa * gr;
      });
}
 
}  // namespace peclet::flow
 
#endif  // PECLET_FLOW_GAUGE_EXACT_GRADIENT_HPP