face_props.hpp

Go to the documentation of this file.

 
#ifndef PECLET_FLOW_FACE_PROPS_HPP
#define PECLET_FLOW_FACE_PROPS_HPP
 
#include <Kokkos_Core.hpp>
 
#include "mac_cutcell.hpp"
 
namespace peclet::flow {
 
// Constant properties — the Var kernel with this accessor is the constant-mu operator (equivalence
// check for the variable path).
struct UniformFaceProps {
  double idiag_, beta_;
  KOKKOS_INLINE_FUNCTION double idiag(long) const { return idiag_; }
  KOKKOS_INLINE_FUNCTION double beta(long, long) const { return beta_; }
};
 
// Per-cell viscosity field (+ constant density for the time diagonal). harmonic=true uses the
// harmonic face mean (continuous shear stress across a viscosity jump); false = arithmetic.
struct FieldFaceProps {
  CCConst mu;
  double rhoIdt;  // rho / dt (constant-density momentum time term)
  bool harmonic;
  KOKKOS_INLINE_FUNCTION double idiag(long) const { return rhoIdt; }
  KOKKOS_INLINE_FUNCTION double beta(long i, long j) const {
    const double a = mu(i), b = mu(j);
    if (harmonic) {
      const double s = a + b;
      return (s > 0.0) ? (2.0 * a * b / s) : 0.0;
    }
    return 0.5 * (a + b);
  }
};
 
// General accessor: viscosity constant-or-field AND density constant-or-field. The time diagonal
// for a variable density is the FACE density of the staggered velocity unknown (component stride
// sc): idiag(i) = 0.5*(rho(i)+rho(i-sc))/dt — the arithmetic mean (mass is volume-additive), and
// the SAME face density the variable-density projection uses, which is what makes discrete
// hydrostatic balance exact. Flags select the constant fallbacks so a single-variable case (only mu
// or only rho) composes.
struct VarFaceProps {
  CCConst mu;
  bool haveMu = false;
  double muC = 0.0;
  bool harmMu = false;
  CCConst rho;
  bool haveRho = false;
  double rhoIdtC = 0.0;  // rho_/dt fallback
  double idt = 0.0;      // 1/dt (variable-rho path)
  long sc = 0;           // component face stride (staggered velocity placement)
  KOKKOS_INLINE_FUNCTION double idiag(long i) const {
    return haveRho ? 0.5 * (rho(i) + rho(i - sc)) * idt : rhoIdtC;
  }
  KOKKOS_INLINE_FUNCTION double beta(long i, long j) const {
    if (!haveMu)
      return muC;
    const double a = mu(i), b = mu(j);
    if (harmMu) {
      const double s = a + b;
      return (s > 0.0) ? (2.0 * a * b / s) : 0.0;
    }
    return 0.5 * (a + b);
  }
};
 
}  // namespace peclet::flow
 
#endif  // PECLET_FLOW_FACE_PROPS_HPP