12#ifndef DUMUX_NAVIERSTOKES_STAGGERED_PROBLEM_HH
13#define DUMUX_NAVIERSTOKES_STAGGERED_PROBLEM_HH
15#include <dune/common/exceptions.hh>
16#include <dune/common/typetraits.hh>
31template<
class TypeTag>
38 using GridView =
typename GridGeometry::GridView;
39 using Element =
typename GridView::template Codim<0>::Entity;
42 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
43 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
44 using FaceVariables =
typename GridFaceVariables::FaceVariables;
45 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
46 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
49 using FVElementGeometry =
typename GridGeometry::LocalView;
50 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
51 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
56 dim = GridView::dimension,
57 dimWorld = GridView::dimensionworld
60 using GlobalPosition =
typename SubControlVolumeFace::GlobalPosition;
61 using VelocityVector = Dune::FieldVector<Scalar, dimWorld>;
62 using GravityVector = Dune::FieldVector<Scalar, dimWorld>;
74 if (getParamFromGroup<bool>(
paramGroup,
"Problem.EnableGravity"))
75 gravity_[dim-1] = -9.81;
77 enableInertiaTerms_ = getParamFromGroup<bool>(
paramGroup,
"Problem.EnableInertiaTerms");
94 {
return enableInertiaTerms_; }
97 template <
class SolutionVector,
class G = Gr
idGeometry>
98 typename std::enable_if<G::discMethod == DiscretizationMethods::staggered, void>::type
100 const SubControlVolumeFace& scvf,
101 const PrimaryVariables& initSol)
const
103 sol[GridGeometry::faceIdx()][scvf.dofIndex()][0] = initSol[Indices::velocity(scvf.directionIndex())];
122 const Scalar factor = 8.0)
const
124 static_assert(dim == 2,
"Pseudo 3D wall friction may only be used in 2D");
125 return -factor * velocity *
viscosity / (height*height);
129 template <
class ElementVolumeVariables,
class ElementFaceVariables,
class G = Gr
idGeometry>
130 typename std::enable_if<G::discMethod == DiscretizationMethods::staggered, Scalar>::type
132 const ElementVolumeVariables& elemVolVars,
133 const ElementFaceVariables& elemFaceVars,
135 const Scalar factor = 8.0)
const
137 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
138 const Scalar
viscosity = elemVolVars[scvf.insideScvIdx()].effectiveViscosity();
147 Scalar
permeability(
const Element& element,
const SubControlVolumeFace& scvf)
const
149 DUNE_THROW(Dune::NotImplemented,
150 "When using the Beavers-Joseph-Saffman boundary condition, "
151 "the permeability must be returned in the actual problem"
160 Scalar
alphaBJ(
const SubControlVolumeFace& scvf)
const
162 DUNE_THROW(Dune::NotImplemented,
163 "When using the Beavers-Joseph-Saffman boundary condition, "
164 "the alpha value must be returned in the actual problem"
171 Scalar
betaBJ(
const Element& element,
const SubControlVolumeFace& scvf,
const GlobalPosition& tangentialVector)
const
173 const Scalar interfacePermeability = interfacePermeability_(element, scvf, tangentialVector);
175 return asImp_().alphaBJ(scvf) / sqrt(interfacePermeability);
183 return VelocityVector(0.0);
190 const SubControlVolume& scv,
191 const SubControlVolumeFace& ownScvf,
192 const SubControlVolumeFace& faceOnPorousBoundary,
193 const Scalar velocitySelf,
194 const Scalar tangentialVelocityGradient)
const
197 GlobalPosition orientation = ownScvf.unitOuterNormal();
198 orientation[ownScvf.directionIndex()] = 1.0;
202 const Scalar
betaBJ = asImp_().betaBJ(element, faceOnPorousBoundary, orientation);
203 const Scalar distanceNormalToBoundary = (faceOnPorousBoundary.center() - scv.center()).two_norm();
205 return (tangentialVelocityGradient*distanceNormalToBoundary
206 + asImp_().porousMediumVelocity(element, faceOnPorousBoundary) * orientation *
betaBJ * distanceNormalToBoundary
207 + velocitySelf) / (
betaBJ*distanceNormalToBoundary + 1.0);
212 Scalar interfacePermeability_(
const Element& element,
const SubControlVolumeFace& scvf,
const GlobalPosition& tangentialVector)
const
214 const auto& K = asImp_().permeability(element, scvf);
217 if constexpr (Dune::IsNumber<std::decay_t<
decltype(K)>>::value)
220 return vtmv(tangentialVector, K, tangentialVector);
224 Implementation &asImp_()
225 {
return *
static_cast<Implementation *
>(
this); }
228 const Implementation &asImp_()
const
229 {
return *
static_cast<const Implementation *
>(
this); }
231 GravityVector gravity_;
232 bool enableInertiaTerms_;
Base class for all finite-volume problems.
Definition: common/fvproblem.hh:43
const std::string & paramGroup() const
The parameter group in which to retrieve runtime parameters.
Definition: common/fvproblem.hh:528
const GridGeometry & gridGeometry() const
The finite volume grid geometry.
Definition: common/fvproblem.hh:524
Navier-Stokes staggered problem base class.
Definition: freeflow/navierstokes/staggered/problem.hh:33
const GravityVector & gravity() const
Returns the acceleration due to gravity.
Definition: freeflow/navierstokes/staggered/problem.hh:87
const Scalar beaversJosephVelocity(const Element &element, const SubControlVolume &scv, const SubControlVolumeFace &ownScvf, const SubControlVolumeFace &faceOnPorousBoundary, const Scalar velocitySelf, const Scalar tangentialVelocityGradient) const
Returns the slip velocity at a porous boundary based on the Beavers-Joseph(-Saffman) condition.
Definition: freeflow/navierstokes/staggered/problem.hh:189
std::enable_if< G::discMethod==DiscretizationMethods::staggered, Scalar >::type pseudo3DWallFriction(const SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const Scalar height, const Scalar factor=8.0) const
Convenience function for staggered grid implementation.
Definition: freeflow/navierstokes/staggered/problem.hh:131
NavierStokesStaggeredProblem(std::shared_ptr< const GridGeometry > gridGeometry, const std::string ¶mGroup="")
The constructor.
Definition: freeflow/navierstokes/staggered/problem.hh:70
bool enableInertiaTerms() const
Returns whether interia terms should be considered.
Definition: freeflow/navierstokes/staggered/problem.hh:93
Scalar permeability(const Element &element, const SubControlVolumeFace &scvf) const
Returns the intrinsic permeability of required as input parameter for the Beavers-Joseph-Saffman boun...
Definition: freeflow/navierstokes/staggered/problem.hh:147
Scalar pseudo3DWallFriction(const Scalar velocity, const Scalar viscosity, const Scalar height, const Scalar factor=8.0) const
An additional drag term can be included as source term for the momentum balance to mimic 3D flow beha...
Definition: freeflow/navierstokes/staggered/problem.hh:119
Scalar alphaBJ(const SubControlVolumeFace &scvf) const
Returns the alpha value required as input parameter for the Beavers-Joseph-Saffman boundary condition...
Definition: freeflow/navierstokes/staggered/problem.hh:160
VelocityVector porousMediumVelocity(const Element &element, const SubControlVolumeFace &scvf) const
Returns the velocity in the porous medium (which is 0 by default according to Saffmann).
Definition: freeflow/navierstokes/staggered/problem.hh:181
Scalar betaBJ(const Element &element, const SubControlVolumeFace &scvf, const GlobalPosition &tangentialVector) const
Returns the beta value which is the alpha value divided by the square root of the (scalar-valued) int...
Definition: freeflow/navierstokes/staggered/problem.hh:171
std::enable_if< G::discMethod==DiscretizationMethods::staggered, void >::type applyInitialFaceSolution(SolutionVector &sol, const SubControlVolumeFace &scvf, const PrimaryVariables &initSol) const
Applies the initial face solution (velocities on the faces). Specialization for staggered grid discre...
Definition: freeflow/navierstokes/staggered/problem.hh:99
Base class for all staggered finite-volume problems.
Definition: staggeredfvproblem.hh:36
Defines all properties used in Dumux.
Dune::DenseMatrix< MAT >::value_type vtmv(const Dune::DenseVector< V1 > &v1, const Dune::DenseMatrix< MAT > &M, const Dune::DenseVector< V2 > &v2)
Evaluates the scalar product of a vector v2, projected by a matrix M, with a vector v1.
Definition: math.hh:851
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The available discretization methods in Dumux.
std::string viscosity(int phaseIdx) noexcept
I/O name of viscosity for multiphase systems.
Definition: name.hh:62
Base class for all staggered fv problems.