24#ifndef DUMUX_DISCRETIZATION_STAGGERED_FREEFLOW_FACEVARIABLES_HH
25#define DUMUX_DISCRETIZATION_STAGGERED_FREEFLOW_FACEVARIABLES_HH
34template<
class Scalar,
int upwindSchemeOrder>
39 std::array<Scalar, upwindSchemeOrder-1>
forward{};
40 std::array<Scalar, upwindSchemeOrder-1>
backward{};
60template<
class FacePrimaryVariables,
int dim,
int upwindSchemeOrder>
63 static constexpr int numPairs = (dim == 2) ? 2 : 4;
64 static constexpr bool useHigherOrder = upwindSchemeOrder > 1;
66 using Scalar =
typename FacePrimaryVariables::block_type;
78 inAxisVelocities_.
self = priVars[0];
91 template<
class SolVector,
class Problem,
class Element,
92 class FVElementGeometry,
class SubControlVolumeFace>
93 void update(
const SolVector& faceSol,
94 const Problem& problem,
95 const Element& element,
96 const FVElementGeometry& fvGeometry,
97 const SubControlVolumeFace& scvf)
99 inAxisVelocities_.
self = faceSol[scvf.dofIndex()];
100 inAxisVelocities_.
opposite = faceSol[scvf.dofIndexOpposingFace()];
102 addHigherOrderInAxisVelocities_(faceSol, scvf, std::integral_constant<bool, useHigherOrder>{});
105 for (
int i = 0; i < velocityParallel_.size(); ++i)
106 velocityParallel_[i].fill(0.0);
108 for (
int i = 0; i < scvf.pairData().size(); ++i)
110 const auto& subFaceData = scvf.pairData(i);
113 velocityLateralInside_[i] = faceSol[subFaceData.lateralPair.first];
115 if (scvf.hasOuterLateral(i))
116 velocityLateralOutside_[i] = faceSol[subFaceData.lateralPair.second];
119 for (
int j = 0; j < upwindSchemeOrder; j++)
121 if (scvf.hasParallelNeighbor(i,j))
122 velocityParallel_[i][j] = faceSol[subFaceData.parallelDofs[j]];
132 return inAxisVelocities_.
self;
148 template<
bool enable = useHigherOrder, std::enable_if_t<enable,
int> = 0>
151 return inAxisVelocities_.
backward[backwardIdx];
159 template<
bool enable = useHigherOrder, std::enable_if_t<enable,
int> = 0>
162 return inAxisVelocities_.
forward[forwardIdx];
173 return velocityParallel_[localSubFaceIdx][parallelDegreeIdx];
183 return velocityLateralInside_[localSubFaceIdx];
193 return velocityLateralOutside_[localSubFaceIdx];
198 template<
class SolVector,
class SubControlVolumeFace>
199 void addHigherOrderInAxisVelocities_(
const SolVector& faceSol,
const SubControlVolumeFace& scvf, std::false_type) {}
201 template<
class SolVector,
class SubControlVolumeFace>
202 void addHigherOrderInAxisVelocities_(
const SolVector& faceSol,
const SubControlVolumeFace& scvf, std::true_type)
207 for (
int i = 0; i < scvf.axisData().inAxisForwardDofs.size(); i++)
209 if (scvf.hasForwardNeighbor(i))
210 inAxisVelocities_.
forward[i]= faceSol[scvf.axisData().inAxisForwardDofs[i]];
215 for (
int i = 0; i < scvf.axisData().inAxisBackwardDofs.size(); i++)
217 if (scvf.hasBackwardNeighbor(i))
218 inAxisVelocities_.
backward[i] = faceSol[scvf.axisData().inAxisBackwardDofs[i]];
222 InAxisVelocities inAxisVelocities_;
223 std::array<std::array<Scalar, upwindSchemeOrder>, numPairs> velocityParallel_;
224 std::array<Scalar, numPairs> velocityLateralInside_;
225 std::array<Scalar, numPairs> velocityLateralOutside_;
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
Definition: facevariables.hh:36
std::array< Scalar, upwindSchemeOrder-1 > forward
Definition: facevariables.hh:39
Scalar opposite
Definition: facevariables.hh:38
std::array< Scalar, upwindSchemeOrder-1 > backward
Definition: facevariables.hh:40
Scalar self
Definition: facevariables.hh:37
The face variables class for free flow staggered grid models. Contains all relevant velocities for th...
Definition: facevariables.hh:62
Scalar velocityForward(const int forwardIdx) const
Returns the velocity at a forward face.
Definition: facevariables.hh:160
void update(const SolVector &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:93
Scalar velocityOpposite() const
Returns the velocity at the opposing face.
Definition: facevariables.hh:138
Scalar velocityParallel(const int localSubFaceIdx, const int parallelDegreeIdx=0) const
Returns the velocity at a parallel face.
Definition: facevariables.hh:171
void updateOwnFaceOnly(const FacePrimaryVariables &priVars)
Partial update of the face variables. Only the face itself is considered.
Definition: facevariables.hh:76
Scalar velocitySelf() const
Returns the velocity at the face itself.
Definition: facevariables.hh:130
Scalar velocityLateralInside(const int localSubFaceIdx) const
Returns the velocity at the inner normal face.
Definition: facevariables.hh:181
Scalar velocityLateralOutside(const int localSubFaceIdx) const
Returns the velocity at the outer normal face.
Definition: facevariables.hh:191
Scalar velocityBackward(const int backwardIdx) const
Returns the velocity at a backward face.
Definition: facevariables.hh:149