12#ifndef DUMUX_DISCRETIZATION_STAGGERED_FREEFLOW_FACEVARIABLES_HH 
   13#define DUMUX_DISCRETIZATION_STAGGERED_FREEFLOW_FACEVARIABLES_HH 
   23template<
class Scalar, 
int upwindSchemeOrder>
 
   28    std::array<Scalar, upwindSchemeOrder-1> 
forward{};
 
   29    std::array<Scalar, upwindSchemeOrder-1> 
backward{};
 
 
   49template<
class FacePrimaryVariables, 
int dim, 
int upwindSchemeOrder>
 
   52    static constexpr int numPairs = (dim == 2) ? 2 : 4;
 
   53    static constexpr bool useHigherOrder = upwindSchemeOrder > 1;
 
   55    using Scalar = 
typename FacePrimaryVariables::block_type;
 
   67        inAxisVelocities_.self = priVars[0];
 
 
   80    template<
class FaceSolution, 
class Problem, 
class Element,
 
   81             class FVElementGeometry, 
class SubControlVolumeFace>
 
   82    void update(
const FaceSolution& faceSol,
 
   83                const Problem& problem,
 
   84                const Element& element,
 
   85                const FVElementGeometry& fvGeometry,
 
   86                const SubControlVolumeFace& scvf)
 
   88        static_assert(std::decay_t<
decltype(faceSol[0])>::dimension == 1,
 
   89                      "\n\n\nVelocity primary variable must be a scalar value. \n\n Make sure to use\n\n ffSol = partial(sol, ffFaceIdx, ffCellCenterIdx);\n\n");
 
   91        inAxisVelocities_.self = faceSol[scvf.dofIndex()];
 
   92        inAxisVelocities_.opposite = faceSol[scvf.dofIndexOpposingFace()];
 
   94        addHigherOrderInAxisVelocities_(faceSol, scvf, std::integral_constant<bool, useHigherOrder>{});
 
   97        for (
int i = 0; i < velocityParallel_.size(); ++i)
 
   98            velocityParallel_[i].fill(0.0);
 
  100        for (
int i = 0; i < scvf.pairData().size(); ++i)
 
  102            const auto& subFaceData = scvf.pairData(i);
 
  105            velocityLateralInside_[i] = faceSol[subFaceData.lateralPair.first];
 
  107            if (scvf.hasOuterLateral(i))
 
  108                velocityLateralOutside_[i] = faceSol[subFaceData.lateralPair.second];
 
  111            for (
int j = 0; j < upwindSchemeOrder; j++)
 
  113                if (scvf.hasParallelNeighbor(i,j))
 
  114                    velocityParallel_[i][j] = faceSol[subFaceData.parallelDofs[j]];
 
 
  124        return inAxisVelocities_.self;
 
 
  132        return inAxisVelocities_.opposite;
 
 
  140    template<
bool enable = useHigherOrder, std::enable_if_t<enable, 
int> = 0>
 
  143        return inAxisVelocities_.backward[backwardIdx];
 
 
  151    template<
bool enable = useHigherOrder, std::enable_if_t<enable, 
int> = 0>
 
  154        return inAxisVelocities_.forward[forwardIdx];
 
 
  165        return velocityParallel_[localSubFaceIdx][parallelDegreeIdx];
 
 
  175        return velocityLateralInside_[localSubFaceIdx];
 
 
  185        return velocityLateralOutside_[localSubFaceIdx];
 
 
  190    template<
class SolVector, 
class SubControlVolumeFace>
 
  191    void addHigherOrderInAxisVelocities_(
const SolVector& faceSol, 
const SubControlVolumeFace& scvf, std::false_type) {}
 
  193    template<
class SolVector, 
class SubControlVolumeFace>
 
  194    void addHigherOrderInAxisVelocities_(
const SolVector& faceSol, 
const SubControlVolumeFace& scvf, std::true_type)
 
  199        for (
int i = 0; i < scvf.axisData().inAxisForwardDofs.size(); i++)
 
  201             if (scvf.hasForwardNeighbor(i))
 
  202                 inAxisVelocities_.
forward[i]= faceSol[scvf.axisData().inAxisForwardDofs[i]];
 
  207        for (
int i = 0; i < scvf.axisData().inAxisBackwardDofs.size(); i++)
 
  209             if (scvf.hasBackwardNeighbor(i))
 
  210                 inAxisVelocities_.backward[i] = faceSol[scvf.axisData().inAxisBackwardDofs[i]];
 
  214    InAxisVelocities inAxisVelocities_;
 
  215    std::array<std::array<Scalar, upwindSchemeOrder>, numPairs> velocityParallel_;
 
  216    std::array<Scalar, numPairs>  velocityLateralInside_;
 
  217    std::array<Scalar, numPairs>  velocityLateralOutside_;
 
 
The face variables class for free flow staggered grid models. Contains all relevant velocities for th...
Definition facevariables.hh:51
Scalar velocityForward(const int forwardIdx) const
Returns the velocity at a forward face.
Definition facevariables.hh:152
void update(const FaceSolution &faceSol, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf)
Complete update of the face variables (i.e. velocities for free flow) for a given face.
Definition facevariables.hh:82
Scalar velocityOpposite() const
Returns the velocity at the opposing face.
Definition facevariables.hh:130
Scalar velocityParallel(const int localSubFaceIdx, const int parallelDegreeIdx=0) const
Returns the velocity at a parallel face.
Definition facevariables.hh:163
void updateOwnFaceOnly(const FacePrimaryVariables &priVars)
Partial update of the face variables. Only the face itself is considered.
Definition facevariables.hh:65
Scalar velocitySelf() const
Returns the velocity at the face itself.
Definition facevariables.hh:122
Scalar velocityLateralInside(const int localSubFaceIdx) const
Returns the velocity at the inner normal face.
Definition facevariables.hh:173
Scalar velocityLateralOutside(const int localSubFaceIdx) const
Returns the velocity at the outer normal face.
Definition facevariables.hh:183
Scalar velocityBackward(const int backwardIdx) const
Returns the velocity at a backward face.
Definition facevariables.hh:141
Distance implementation details.
Definition cvfelocalresidual.hh:25
Scalar opposite
Definition facevariables.hh:36
Scalar self
Definition facevariables.hh:35
Definition facevariables.hh:25
std::array< Scalar, upwindSchemeOrder-1 > forward
Definition facevariables.hh:28
Scalar opposite
Definition facevariables.hh:27
std::array< Scalar, upwindSchemeOrder-1 > backward
Definition facevariables.hh:29
Scalar self
Definition facevariables.hh:26