12#ifndef DUMUX_NAVIERSTOKES_STAGGERED_FLUXVARIABLES_HH
13#define DUMUX_NAVIERSTOKES_STAGGERED_FLUXVARIABLES_HH
40template<
class TypeTag,
class DiscretizationMethod>
41class NavierStokesFluxVariablesImpl;
43template<
class TypeTag>
46 typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView,
47 typename GetPropType<TypeTag, Properties::GridVolumeVariables>::LocalView,
48 typename GetPropType<TypeTag, Properties::GridFluxVariablesCache>::LocalView>
52 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
53 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
54 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
56 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
57 using FluxVariablesCache =
typename GridFluxVariablesCache::FluxVariablesCache;
59 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
60 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
61 using FaceVariables =
typename GridFaceVariables::FaceVariables;
65 using FVElementGeometry =
typename GridGeometry::LocalView;
66 using GridView =
typename GridGeometry::GridView;
67 using Element =
typename GridView::template Codim<0>::Entity;
70 using Indices =
typename ModelTraits::Indices;
71 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
72 using FaceFrontalSubControlVolumeFace =
typename GridGeometry::Traits::FaceFrontalSubControlVolumeFace;
73 using FaceLateralSubControlVolumeFace =
typename GridGeometry::Traits::FaceLateralSubControlVolumeFace;
80 static constexpr bool normalizePressure = getPropValue<TypeTag, Properties::NormalizePressure>();
82 using GlobalPosition =
typename Element::Geometry::GlobalCoordinate;
84 static constexpr auto cellCenterIdx = GridGeometry::cellCenterIdx();
85 static constexpr auto faceIdx = GridGeometry::faceIdx();
87 static constexpr int upwindSchemeOrder = getPropValue<TypeTag, Properties::UpwindSchemeOrder>();
88 static constexpr bool useHigherOrder = upwindSchemeOrder > 1;
104 template<
class UpwindTerm>
106 const FVElementGeometry& fvGeometry,
107 const ElementVolumeVariables& elemVolVars,
108 const ElementFaceVariables& elemFaceVars,
109 const SubControlVolumeFace &scvf,
110 UpwindTerm upwindTerm)
112 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
113 const bool insideIsUpstream = scvf.directionSign() ==
sign(velocity);
114 static const Scalar upwindWeight = getParamFromGroup<Scalar>(problem.paramGroup(),
"Flux.UpwindWeight");
116 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
117 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
119 const auto& upstreamVolVars = insideIsUpstream ? insideVolVars : outsideVolVars;
120 const auto& downstreamVolVars = insideIsUpstream ? outsideVolVars : insideVolVars;
122 const Scalar flux = (upwindWeight * upwindTerm(upstreamVolVars) +
123 (1.0 - upwindWeight) * upwindTerm(downstreamVolVars))
124 * velocity * Extrusion::area(fvGeometry, scvf) * scvf.directionSign();
126 return flux * extrusionFactor_(elemVolVars, scvf);
145 const Element& element,
146 const FVElementGeometry& fvGeometry,
147 const ElementVolumeVariables& elemVolVars,
148 const ElementFaceVariables& elemFaceVars,
149 const SubControlVolumeFace& scvf,
150 const FluxVariablesCache& fluxVarsCache)
153 auto upwindTerm = [](
const auto& volVars) {
return volVars.density(); };
156 CellCenterPrimaryVariables result(0.0);
157 result[Indices::conti0EqIdx - ModelTraits::dim()] = advectiveFluxForCellCenter(problem, fvGeometry, elemVolVars, elemFaceVars, scvf, upwindTerm);
166 const Element& element,
167 const SubControlVolumeFace& scvf,
168 const FVElementGeometry& fvGeometry,
169 const ElementVolumeVariables& elemVolVars,
170 const ElementFaceVariables& elemFaceVars,
171 const GridFluxVariablesCache& gridFluxVarsCache)
173 return computeFrontalMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache) +
174 computeLateralMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache);
195 const Element& element,
196 const SubControlVolumeFace& scvf,
197 const FVElementGeometry& fvGeometry,
198 const ElementVolumeVariables& elemVolVars,
199 const ElementFaceVariables& elemFaceVars,
200 const GridFluxVariablesCache& gridFluxVarsCache)
202 FacePrimaryVariables frontalFlux(0.0);
205 const Scalar velocitySelf = elemFaceVars[scvf].velocitySelf();
206 const Scalar velocityOpposite = elemFaceVars[scvf].velocityOpposite();
207 const auto& faceVars = elemFaceVars[scvf];
210 if (problem.enableInertiaTerms())
213 const Scalar transportingVelocity = (velocitySelf + velocityOpposite) * 0.5;
214 const bool selfIsUpstream = scvf.directionSign() !=
sign(transportingVelocity);
218 * transportingVelocity * -1.0 * scvf.directionSign();
223 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
226 const Scalar velocityGrad_ii = VelocityGradients::velocityGradII(scvf, faceVars) * scvf.directionSign();
228 static const bool enableUnsymmetrizedVelocityGradient
229 = getParamFromGroup<bool>(problem.paramGroup(),
"FreeFlow.EnableUnsymmetrizedVelocityGradient",
false);
230 const Scalar factor = enableUnsymmetrizedVelocityGradient ? 1.0 : 2.0;
231 frontalFlux -= factor * insideVolVars.effectiveViscosity() * velocityGrad_ii;
237 const Scalar
pressure = normalizePressure ?
238 insideVolVars.pressure() - problem.initial(scvf)[Indices::pressureIdx]
239 : insideVolVars.pressure();
243 frontalFlux +=
pressure * -1.0 * scvf.directionSign();
247 const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
248 FaceFrontalSubControlVolumeFace frontalFace(scv.center(), scvf.area());
249 return frontalFlux * Extrusion::area(fvGeometry, frontalFace) * insideVolVars.extrusionFactor();
270 const Element& element,
271 const SubControlVolumeFace& scvf,
272 const FVElementGeometry& fvGeometry,
273 const ElementVolumeVariables& elemVolVars,
274 const ElementFaceVariables& elemFaceVars,
275 const GridFluxVariablesCache& gridFluxVarsCache)
277 FacePrimaryVariables lateralFlux(0.0);
278 const auto& faceVars = elemFaceVars[scvf];
279 const std::size_t numSubFaces = scvf.pairData().size();
283 std::optional<BoundaryTypes> currentScvfBoundaryTypes;
285 currentScvfBoundaryTypes.emplace(problem.boundaryTypes(element, scvf));
288 for (
int localSubFaceIdx = 0; localSubFaceIdx < numSubFaces; ++localSubFaceIdx)
290 const auto eIdx = scvf.insideScvIdx();
292 const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx);
295 std::optional<BoundaryTypes> lateralFaceBoundaryTypes;
299 if (lateralFace.boundary())
312 lateralFaceBoundaryTypes.emplace(problem.boundaryTypes(element, lateralFace));
318 if (currentScvfBoundaryTypes)
321 if (currentScvfBoundaryTypes->isNeumann(Indices::velocity(lateralFace.directionIndex())))
324 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
325 const auto& lateralStaggeredFaceCenter = lateralStaggeredFaceCenter_(scvf, localSubFaceIdx);
327 lateralFlux += problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, lateralBoundaryFace)[Indices::velocity(lateralFace.directionIndex())]
328 * extrusionFactor_(elemVolVars, lateralFace) * Extrusion::area(fvGeometry, lateralScvf) * lateralFace.directionSign();
337 if (currentScvfBoundaryTypes->isBeaversJoseph(Indices::velocity(lateralFace.directionIndex())) && lateralFaceBoundaryTypes &&
338 lateralFaceBoundaryTypes->isNeumann(Indices::velocity(scvf.directionIndex())))
340 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
341 const auto& lateralStaggeredFaceCenter = lateralStaggeredFaceCenter_(scvf, localSubFaceIdx);
343 lateralFlux += problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, lateralBoundaryFace)[Indices::velocity(scvf.directionIndex())]
344 * extrusionFactor_(elemVolVars, lateralFace) * Extrusion::area(fvGeometry, lateralScvf) * lateralFace.directionSign();
352 if (lateralFace.boundary())
356 if (lateralFaceBoundaryTypes->isSymmetry())
360 if (lateralFaceBoundaryTypes->isNeumann(Indices::velocity(scvf.directionIndex())))
364 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
365 const auto& lateralStaggeredFaceCenter = lateralStaggeredFaceCenter_(scvf, localSubFaceIdx);
367 lateralFlux += problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, lateralBoundaryFace)[Indices::velocity(scvf.directionIndex())]
368 * elemVolVars[lateralFace.insideScvIdx()].extrusionFactor() * Extrusion::area(fvGeometry, lateralScvf);
374 if (incorporateWallFunction_(lateralFlux, problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars, localSubFaceIdx))
379 if (lateralFaceBoundaryTypes)
381 std::bitset<3> admittableBcTypes;
382 admittableBcTypes.set(0, lateralFaceBoundaryTypes->isDirichlet(Indices::pressureIdx));
383 admittableBcTypes.set(1, lateralFaceBoundaryTypes->isDirichlet(Indices::velocity(scvf.directionIndex())));
384 admittableBcTypes.set(2, lateralFaceBoundaryTypes->isBeaversJoseph(Indices::velocity(scvf.directionIndex())));
385 if (admittableBcTypes.count() != 1)
387 DUNE_THROW(Dune::InvalidStateException,
"Invalid boundary conditions for lateral scvf "
388 "for the momentum equations at global position " << lateralStaggeredFaceCenter_(scvf, localSubFaceIdx)
389 <<
", current scvf global position " << scvf.center());
394 if (problem.enableInertiaTerms())
395 lateralFlux += computeAdvectivePartOfLateralMomentumFlux_(problem, fvGeometry, element,
396 scvf, elemVolVars, elemFaceVars,
398 currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
401 lateralFlux += computeDiffusivePartOfLateralMomentumFlux_(problem, fvGeometry, element,
402 scvf, elemVolVars, faceVars,
403 currentScvfBoundaryTypes, lateralFaceBoundaryTypes,
427 const SubControlVolumeFace& scvf,
428 const FVElementGeometry& fvGeometry,
429 const ElementVolumeVariables& elemVolVars,
430 const ElementFaceVariables& elemFaceVars)
const
432 FacePrimaryVariables inOrOutflow(0.0);
433 const auto& element = fvGeometry.element();
434 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
437 if (problem.enableInertiaTerms())
439 const Scalar velocitySelf = elemFaceVars[scvf].velocitySelf();
440 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
441 const auto& upVolVars = (scvf.directionSign() ==
sign(velocitySelf)) ?
442 insideVolVars : outsideVolVars;
444 inOrOutflow += velocitySelf * velocitySelf * upVolVars.density();
448 const Scalar boundaryPressure = normalizePressure
449 ? (problem.dirichlet(element, scvf)[Indices::pressureIdx] -
450 problem.initial(scvf)[Indices::pressureIdx])
451 : problem.dirichlet(element, scvf)[Indices::pressureIdx];
452 inOrOutflow += boundaryPressure;
455 return inOrOutflow * scvf.directionSign() * Extrusion::area(fvGeometry, scvf) * insideVolVars.extrusionFactor();
481 FacePrimaryVariables computeAdvectivePartOfLateralMomentumFlux_(
const Problem& problem,
482 const FVElementGeometry& fvGeometry,
483 const Element& element,
484 const SubControlVolumeFace& scvf,
485 const ElementVolumeVariables& elemVolVars,
486 const ElementFaceVariables& elemFaceVars,
487 const GridFluxVariablesCache& gridFluxVarsCache,
488 const std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
489 const std::optional<BoundaryTypes>& lateralFaceBoundaryTypes,
490 const int localSubFaceIdx)
492 const auto eIdx = scvf.insideScvIdx();
493 const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx);
497 const Scalar transportingVelocity = [&]()
499 const auto& faceVars = elemFaceVars[scvf];
500 if (!scvf.boundary())
501 return faceVars.velocityLateralInside(localSubFaceIdx);
505 const auto bcTypes = problem.boundaryTypes(element, scvf);
507 if (bcTypes.isDirichlet(Indices::velocity(lateralFace.directionIndex())))
511 const auto& lateralBoundaryFacePos = lateralStaggeredFaceCenter_(scvf, localSubFaceIdx);
513 return problem.dirichlet(element, lateralBoundaryFace)[Indices::velocity(lateralFace.directionIndex())];
515 else if (bcTypes.isBeaversJoseph(Indices::velocity(lateralFace.directionIndex())))
517 return VelocityGradients::beaversJosephVelocityAtCurrentScvf(problem, element, fvGeometry, scvf, faceVars,
518 currentScvfBoundaryTypes, lateralFaceBoundaryTypes, localSubFaceIdx);
521 return faceVars.velocityLateralInside(localSubFaceIdx);
525 const bool selfIsUpstream = lateralFace.directionSign() ==
sign(transportingVelocity);
526 StaggeredUpwindHelper<TypeTag, upwindSchemeOrder> upwindHelper(element, fvGeometry, scvf, elemFaceVars, elemVolVars, gridFluxVarsCache.staggeredUpwindMethods());
527 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
528 return upwindHelper.computeUpwindLateralMomentum(selfIsUpstream, lateralFace, localSubFaceIdx, currentScvfBoundaryTypes, lateralFaceBoundaryTypes)
529 * transportingVelocity * lateralFace.directionSign() * Extrusion::area(fvGeometry, lateralScvf) * extrusionFactor_(elemVolVars, lateralFace);
554 FacePrimaryVariables computeDiffusivePartOfLateralMomentumFlux_(
const Problem& problem,
555 const FVElementGeometry& fvGeometry,
556 const Element& element,
557 const SubControlVolumeFace& scvf,
558 const ElementVolumeVariables& elemVolVars,
559 const FaceVariables& faceVars,
560 const std::optional<BoundaryTypes>& currentScvfBoundaryTypes,
561 const std::optional<BoundaryTypes>& lateralFaceBoundaryTypes,
562 const int localSubFaceIdx)
564 const auto eIdx = scvf.insideScvIdx();
565 const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx);
567 FacePrimaryVariables lateralDiffusiveFlux(0.0);
569 static const bool enableUnsymmetrizedVelocityGradient
570 = getParamFromGroup<bool>(problem.paramGroup(),
"FreeFlow.EnableUnsymmetrizedVelocityGradient",
false);
575 const auto& insideVolVars = elemVolVars[lateralFace.insideScvIdx()];
576 const auto& outsideVolVars = elemVolVars[lateralFace.outsideScvIdx()];
579 const Scalar muAvg = lateralFace.boundary()
580 ? insideVolVars.effectiveViscosity()
581 : (insideVolVars.effectiveViscosity() + outsideVolVars.effectiveViscosity()) * 0.5;
584 if (!enableUnsymmetrizedVelocityGradient)
586 if (!scvf.boundary() ||
587 currentScvfBoundaryTypes->isDirichlet(Indices::velocity(lateralFace.directionIndex())) ||
588 currentScvfBoundaryTypes->isBeaversJoseph(Indices::velocity(lateralFace.directionIndex())))
590 const Scalar velocityGrad_ji = VelocityGradients::velocityGradJI(problem, element, fvGeometry, scvf, faceVars, currentScvfBoundaryTypes, lateralFaceBoundaryTypes, localSubFaceIdx);
592 lateralDiffusiveFlux -= muAvg * velocityGrad_ji * lateralFace.directionSign();
599 if (!lateralFace.boundary() || !lateralFaceBoundaryTypes->isDirichlet(Indices::pressureIdx))
601 const Scalar velocityGrad_ij = VelocityGradients::velocityGradIJ(problem, element, fvGeometry, scvf, faceVars, currentScvfBoundaryTypes, lateralFaceBoundaryTypes, localSubFaceIdx);
602 lateralDiffusiveFlux -= muAvg * velocityGrad_ij * lateralFace.directionSign();
606 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
607 return lateralDiffusiveFlux * Extrusion::area(fvGeometry, lateralScvf) * extrusionFactor_(elemVolVars, lateralFace);
624 const GlobalPosition& lateralStaggeredFaceCenter_(
const SubControlVolumeFace& scvf,
const int localSubFaceIdx)
const
626 return scvf.pairData(localSubFaceIdx).lateralStaggeredFaceCenter;
643 GlobalPosition lateralStaggeredSCVFCenter_(
const SubControlVolumeFace& lateralFace,
644 const SubControlVolumeFace& currentFace,
645 const int localSubFaceIdx)
const
647 auto pos = lateralStaggeredFaceCenter_(currentFace, localSubFaceIdx) + lateralFace.center();
653 static Scalar extrusionFactor_(
const ElementVolumeVariables& elemVolVars,
const SubControlVolumeFace& scvf)
655 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
656 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
657 return harmonicMean(insideVolVars.extrusionFactor(), outsideVolVars.extrusionFactor());
661 template<
class ...Args,
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<!turbulenceModel, int> = 0>
662 bool incorporateWallFunction_(Args&&... args)
const
666 template<
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<turbulenceModel,
int> = 0>
667 bool incorporateWallFunction_(FacePrimaryVariables& lateralFlux,
668 const Problem& problem,
669 const Element& element,
670 const FVElementGeometry& fvGeometry,
671 const SubControlVolumeFace& scvf,
672 const ElementVolumeVariables& elemVolVars,
673 const ElementFaceVariables& elemFaceVars,
674 const std::size_t localSubFaceIdx)
const
676 const auto eIdx = scvf.insideScvIdx();
677 const auto& lateralFace = fvGeometry.scvf(eIdx, scvf.pairData(localSubFaceIdx).localLateralFaceIdx);
679 if (problem.useWallFunction(element, lateralFace, Indices::velocity(scvf.directionIndex())))
681 FaceLateralSubControlVolumeFace lateralScvf(lateralStaggeredSCVFCenter_(lateralFace, scvf, localSubFaceIdx), 0.5*lateralFace.area());
682 const auto lateralBoundaryFace =
makeStaggeredBoundaryFace(lateralFace, lateralStaggeredFaceCenter_(scvf, localSubFaceIdx));
683 lateralFlux += problem.wallFunction(element, fvGeometry, elemVolVars, elemFaceVars, scvf, lateralBoundaryFace)[Indices::velocity(scvf.directionIndex())]
684 * extrusionFactor_(elemVolVars, lateralFace) * Extrusion::area(fvGeometry, lateralScvf);
Base class for the flux variables living on a sub control volume face.
Definition: fluxvariablesbase.hh:33
static Scalar advectiveFluxForCellCenter(const Problem &problem, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, UpwindTerm upwindTerm)
Returns the advective flux over a sub control volume face.
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:105
FacePrimaryVariables computeLateralMomentumFlux(const Problem &problem, const Element &element, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const GridFluxVariablesCache &gridFluxVarsCache)
Returns the momentum flux over the staggered faces perpendicular to the scvf where the velocity dof o...
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:269
FacePrimaryVariables inflowOutflowBoundaryFlux(const Problem &problem, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars) const
Returns the momentum flux over an inflow or outflow boundary face.
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:426
FacePrimaryVariables computeFrontalMomentumFlux(const Problem &problem, const Element &element, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const GridFluxVariablesCache &gridFluxVarsCache)
Returns the frontal part of the momentum flux. This treats the flux over the staggered face at the ce...
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:194
FacePrimaryVariables computeMomentumFlux(const Problem &problem, const Element &element, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const GridFluxVariablesCache &gridFluxVarsCache)
Returns the momentum flux over all staggered faces.
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:165
GetPropType< TypeTag, Properties::HeatConductionType > HeatConductionType
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:92
CellCenterPrimaryVariables computeMassFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const FluxVariablesCache &fluxVarsCache)
Computes the flux for the cell center residual (mass balance).
Definition: freeflow/navierstokes/staggered/fluxvariables.hh:144
The flux variables class for the Navier-Stokes model using the staggered grid discretization.
Definition: freeflow/navierstokes/fluxvariables.hh:23
The upwinding variables class for the Navier-Stokes model using the staggered grid discretization.
Definition: staggeredupwindhelper.hh:36
FacePrimaryVariables computeUpwindFrontalMomentum(const bool selfIsUpstream) const
Returns the momentum in the frontal direction.
Definition: staggeredupwindhelper.hh:103
Helper class for calculating the velocity gradients for the Navier-Stokes model using the staggered g...
Definition: momentum/velocitygradients.hh:46
Defines all properties used in Dumux.
Type traits for problem classes.
Some exceptions thrown in DuMux
Helper classes to compute the integration elements.
Base class for the flux variables living on a sub control volume face.
constexpr Scalar harmonicMean(Scalar x, Scalar y, Scalar wx=1.0, Scalar wy=1.0) noexcept
Calculate the (weighted) harmonic mean of two scalar values.
Definition: math.hh:57
constexpr int sign(const ValueType &value) noexcept
Sign or signum function.
Definition: math.hh:629
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
SubControlVolumeFace makeStaggeredBoundaryFace(const SubControlVolumeFace &scvf, const typename SubControlVolumeFace::GlobalPosition &newCenter)
Helper function to turn a given cell scvface into a fake boundary face.
Definition: discretization/staggered/freeflow/subcontrolvolumeface.hh:62
Define some often used mathematical functions.
The available discretization methods in Dumux.
std::string pressure(int phaseIdx) noexcept
I/O name of pressure for multiphase systems.
Definition: name.hh:22
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:166
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
typename Detail::template ProblemTraits< Problem, typename GridGeometry::DiscretizationMethod >::BoundaryTypes BoundaryTypes
Definition: common/typetraits/problem.hh:34