diff --git a/dumux/discretization/staggered/freeflow/facevariables.hh b/dumux/discretization/staggered/freeflow/facevariables.hh
index 5cadb11af6f3537d9bb8a10cde4d4e59627ff248..b84abd27c42d494478df2c7a86fb5ad9d8f82809 100644
--- a/dumux/discretization/staggered/freeflow/facevariables.hh
+++ b/dumux/discretization/staggered/freeflow/facevariables.hh
@@ -86,24 +86,6 @@ public:
velocitySelf_ = faceSol[scvf.dofIndex()];
velocityOpposite_ = faceSol[scvf.dofIndexOpposingFace()];
- // lambda to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest
- auto makeGhostFace = [](const auto& pos)
- {
- return SubControlVolumeFace(pos, std::vector<unsigned int>{0,0});
- };
-
- // lambda to check whether there is a parallel face neighbor
- auto hasParallelNeighbor = [](const auto& subFaceData)
- {
- return subFaceData.outerParallelFaceDofIdx >= 0;
- };
-
- // lambda to check whether there is a normal face neighbor
- auto hasNormalNeighbor = [](const auto& subFaceData)
- {
- return subFaceData.normalPair.second >= 0;
- };
-
// handle all sub faces
for(int i = 0; i < scvf.pairData().size(); ++i)
{
@@ -112,21 +94,12 @@ public:
// treat the velocities normal to the face
velocityNormalInside_[i] = faceSol[subFaceData.normalPair.first];
- if(hasNormalNeighbor(subFaceData))
- {
+ if(scvf.hasFrontalNeighbor(i))
velocityNormalOutside_[i] = faceSol[subFaceData.normalPair.second];
- }
- else
- {
- const auto& normalFace = fvGeometry.scvf(scvf.insideScvIdx(), subFaceData.localNormalFaceIdx);
- const auto normalDirIdx = normalFace.directionIndex();
- velocityNormalOutside_[i] = problem.dirichlet(element, makeGhostFace(subFaceData.virtualOuterNormalFaceDofPos))[Indices::velocity(normalDirIdx)];
- }
// treat the velocity parallel to the face
- velocityParallel_[i] = hasParallelNeighbor(subFaceData) ?
- velocityParallel_[i] = faceSol[subFaceData.outerParallelFaceDofIdx] :
- problem.dirichlet(element, makeGhostFace(subFaceData.virtualOuterParallelFaceDofPos))[Indices::velocity(scvf.directionIndex())];
+ if(scvf.hasParallelNeighbor(i))
+ velocityParallel_[i] = faceSol[subFaceData.outerParallelFaceDofIdx];
}
}
diff --git a/dumux/discretization/staggered/freeflow/subcontrolvolumeface.hh b/dumux/discretization/staggered/freeflow/subcontrolvolumeface.hh
index be10af48e93e590ec2d76bd46bd4576e4626da73..426657563cfc9ce92f8fac15456551a3391028f9 100644
--- a/dumux/discretization/staggered/freeflow/subcontrolvolumeface.hh
+++ b/dumux/discretization/staggered/freeflow/subcontrolvolumeface.hh
@@ -31,7 +31,6 @@
#include <dumux/discretization/staggered/subcontrolvolumeface.hh>
#include <dumux/discretization/staggered/freeflow/staggeredgeometryhelper.hh>
#include <dumux/common/properties.hh>
-#include <dumux/common/optional.hh>
#include <typeinfo>
@@ -82,34 +81,32 @@ public:
unitOuterNormal_(is.centerUnitOuterNormal()),
scvfIndex_(scvfIndex),
scvIndices_(scvIndices),
- boundary_(is.boundary())
+ boundary_(is.boundary()),
+ dofIdx_(geometryHelper.dofIndex()),
+ dofIdxOpposingFace_(geometryHelper.dofIndexOpposingFace()),
+ selfToOppositeDistance_(geometryHelper.selfToOppositeDistance()),
+ pairData_(geometryHelper.pairData()),
+ localFaceIdx_(geometryHelper.localFaceIndex()),
+ dirIdx_(geometryHelper.directionIndex()),
+ normalInPosCoordDir_(unitOuterNormal_[directionIndex()] > 0.0),
+ outerNormalScalar_(unitOuterNormal_[directionIndex()]),
+ isGhostFace_(false)
{
corners_.resize(isGeometry.corners());
for (int i = 0; i < isGeometry.corners(); ++i)
corners_[i] = isGeometry.corner(i);
-
- dofIdx_ = geometryHelper.dofIndex();
- oppositeIdx_ = geometryHelper.dofIndexOpposingFace();
- selfToOppositeDistance_ = geometryHelper.selfToOppositeDistance();
-
- pairData_ = geometryHelper.pairData();
- localFaceIdx_ = geometryHelper.localFaceIndex();
- dirIdx_ = geometryHelper.directionIndex();
- normalInPosCoordDir_ = unitOuterNormal()[directionIndex()] > 0.0;
- outerNormalScalar_ = unitOuterNormal()[directionIndex()];
- isGhostFace_ = false;
}
- //! Constructor for a ghost face outside of the domain. Only needed to retrieve the center and scvIndices
- FreeFlowStaggeredSubControlVolumeFace(const GlobalPosition& dofPosition,
- const std::vector<GridIndexType>& scvIndices)
- {
- isGhostFace_ = true;
- center_ = dofPosition;
- scvIndices_ = scvIndices;
- scvfIndex_ = -1;
- dofIdx_ = -1;
- }
+ //! Constructor for a ghost face outside of the domain. Only needed to retrieve the center and scvIndices
+ FreeFlowStaggeredSubControlVolumeFace(const GlobalPosition& dofPosition,
+ const std::vector<GridIndexType>& scvIndices,
+ const int dofIdx = -1)
+ : center_(dofPosition),
+ scvfIndex_(-1),
+ scvIndices_(scvIndices),
+ dofIdx_(dofIdx),
+ isGhostFace_(true)
+ {}
//! The center of the sub control volume face
const GlobalPosition& center() const
@@ -189,7 +186,7 @@ public:
//! The global index of the dof living on the opposing face
GridIndexType dofIndexOpposingFace() const
{
- return oppositeIdx_;
+ return dofIdxOpposingFace_;
}
//! The local index of this sub control volume face
@@ -234,6 +231,21 @@ public:
return pairData_;
}
+ bool isGhostFace() const
+ {
+ return isGhostFace_;
+ }
+
+ bool hasParallelNeighbor(const int localFaceIdx) const
+ {
+ return pairData_[localFaceIdx].outerParallelFaceDofIdx >= 0;
+ }
+
+ bool hasFrontalNeighbor(const int localFaceIdx) const
+ {
+ return pairData_[localFaceIdx].normalPair.second >= 0;
+ }
+
private:
Dune::GeometryType geomType_;
std::vector<GlobalPosition> corners_;
@@ -245,7 +257,7 @@ private:
bool boundary_;
int dofIdx_;
- int oppositeIdx_;
+ int dofIdxOpposingFace_;
Scalar selfToOppositeDistance_;
std::array<PairData<Scalar, GlobalPosition>, numPairs> pairData_;
int localFaceIdx_;
diff --git a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
index 4df4442acc2c92ef19e51ce82970e715c33e0bb4..d693511e762366607c143f7e234c5637f92d2afe 100644
--- a/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
+++ b/dumux/freeflow/navierstokes/staggered/fluxvariables.hh
@@ -283,16 +283,10 @@ public:
const auto& normalFace = fvGeometry.scvf(eIdx, scvf.pairData()[localSubFaceIdx].localNormalFaceIdx);
// Check if we have a symmetry boundary condition. If yes, the tangental part of the momentum flux can be neglected.
- if(scvf.pairData()[localSubFaceIdx].outerParallelFaceDofIdx < 0)
+ if(!scvf.hasParallelNeighbor(localSubFaceIdx))
{
- // Lambda to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest.
- auto makeGhostFace = [eIdx] (const GlobalPosition& pos)
- {
- return SubControlVolumeFace(pos, std::vector<unsigned int>{eIdx,eIdx});
- };
-
// Use the ghost face to check if there is a symmetry boundary condition and skip any further steps if yes.
- const auto bcTypes = problem.boundaryTypes(element, makeGhostFace(scvf.pairData()[localSubFaceIdx].virtualOuterParallelFaceDofPos));
+ const auto bcTypes = problem.boundaryTypes(element, makeParallelGhostFace_(scvf, localSubFaceIdx));
if(bcTypes.isSymmetry())
continue;
}
@@ -346,7 +340,18 @@ private:
// Get the velocities at the current (own) scvf and at the parallel one at the neighboring scvf.
const Scalar velocitySelf = faceVars.velocitySelf();
- const Scalar velocityParallel = faceVars.velocityParallel(localSubFaceIdx);
+
+ // Lambda to conveniently get the outer parallel velocity for normal faces that are on the boundary
+ // and therefore have no neighbor. Calls the problem to retrieve a fixed value set on the boundary.
+ auto getParallelVelocityFromBoundary = [&]()
+ {
+ const auto ghostFace = makeParallelGhostFace_(scvf, localSubFaceIdx);
+ return problem.dirichlet(element, ghostFace)[Indices::velocity(scvf.directionIndex())];
+ };
+
+ const Scalar velocityParallel = scvf.hasParallelNeighbor(localSubFaceIdx) ?
+ faceVars.velocityParallel(localSubFaceIdx)
+ : getParallelVelocityFromBoundary();
// Get the volume variables of the own and the neighboring element
const auto& insideVolVars = elemVolVars[normalFace.insideScvIdx()];
@@ -420,7 +425,19 @@ private:
// the outer one at the respective staggered face of the element on the other side of the
// current scvf.
const Scalar innerNormalVelocity = faceVars.velocityNormalInside(localSubFaceIdx);
- const Scalar outerNormalVelocity = faceVars.velocityNormalOutside(localSubFaceIdx);
+
+ // Lambda to conveniently get the outer normal velocity for faces that are on the boundary
+ // and therefore have no neighbor. Calls the problem to retrieve a fixed value set on the boundary.
+ auto getNormalVelocityFromBoundary = [&]()
+ {
+ const auto ghostFace = makeNormalGhostFace_(scvf, localSubFaceIdx);
+ return problem.dirichlet(element, ghostFace)[Indices::velocity(normalFace.directionIndex())];
+ };
+
+ const Scalar outerNormalVelocity = scvf.hasFrontalNeighbor(localSubFaceIdx) ?
+ faceVars.velocityNormalOutside(localSubFaceIdx)
+ : getNormalVelocityFromBoundary();
+
// Calculate the velocity gradient in positive coordinate direction.
const Scalar normalDeltaV = scvf.normalInPosCoordDir() ?
@@ -435,7 +452,18 @@ private:
// For the parallel derivative, get the velocities at the current (own) scvf
// and at the parallel one at the neighboring scvf.
const Scalar innerParallelVelocity = faceVars.velocitySelf();
- const Scalar outerParallelVelocity = faceVars.velocityParallel(localSubFaceIdx);
+
+ // Lambda to conveniently get the outer parallel velocity for normal faces that are on the boundary
+ // and therefore have no neighbor. Calls the problem to retrieve a fixed value set on the boundary.
+ auto getParallelVelocityFromBoundary = [&]()
+ {
+ const auto ghostFace = makeParallelGhostFace_(scvf, localSubFaceIdx);
+ return problem.dirichlet(element, ghostFace)[Indices::velocity(scvf.directionIndex())];
+ };
+
+ const Scalar outerParallelVelocity = scvf.hasParallelNeighbor(localSubFaceIdx) ?
+ faceVars.velocityParallel(localSubFaceIdx)
+ : getParallelVelocityFromBoundary();
// The velocity gradient already accounts for the orientation
// of the staggered face's outer normal vector.
@@ -494,6 +522,27 @@ private:
// Account for the orientation of the face at the boundary,
return outflow * scvf.directionSign();
}
+
+private:
+
+ //! helper functiuob to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest
+ SubControlVolumeFace makeGhostFace_(const SubControlVolumeFace& ownScvf, const GlobalPosition& pos) const
+ {
+ return SubControlVolumeFace(pos, std::vector<unsigned int>{ownScvf.insideScvIdx(), ownScvf.outsideScvIdx()}, ownScvf.dofIndex());
+ };
+
+ //! helper functiuob to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest
+ SubControlVolumeFace makeNormalGhostFace_(const SubControlVolumeFace& ownScvf, const int localSubFaceIdx) const
+ {
+ return makeGhostFace_(ownScvf, ownScvf.pairData(localSubFaceIdx).virtualOuterNormalFaceDofPos);
+ };
+
+ //! helper functiuob to conveniently create a ghost face which is outside the domain, parallel to the scvf of interest
+ SubControlVolumeFace makeParallelGhostFace_(const SubControlVolumeFace& ownScvf, const int localSubFaceIdx) const
+ {
+ return makeGhostFace_(ownScvf, ownScvf.pairData(localSubFaceIdx).virtualOuterParallelFaceDofPos);
+ };
+
};
} // end namespace