12#ifndef DUMUX_KEPSILON_STAGGERED_FLUXVARIABLES_HH
13#define DUMUX_KEPSILON_STAGGERED_FLUXVARIABLES_HH
31template<
class TypeTag,
class BaseFluxVariables,
class DiscretizationMethod>
32class KEpsilonFluxVariablesImpl;
34template<
class TypeTag,
class BaseFluxVariables>
36:
public BaseFluxVariables
38 using ParentType = BaseFluxVariables;
42 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
43 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
44 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
46 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
47 using FluxVariablesCache =
typename GridFluxVariablesCache::FluxVariablesCache;
49 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
50 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
51 using FaceVariables =
typename GridFaceVariables::FaceVariables;
56 using FVElementGeometry =
typename GridGeometry::LocalView;
57 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
58 using GlobalPosition =
typename SubControlVolumeFace::GlobalPosition;
60 using GridView =
typename GridGeometry::GridView;
62 using Element =
typename GridView::template Codim<0>::Entity;
67 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
68 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
76 const Element &element,
77 const FVElementGeometry& fvGeometry,
78 const ElementVolumeVariables& elemVolVars,
79 const ElementFaceVariables& elemFaceVars,
80 const SubControlVolumeFace &scvf,
81 const FluxVariablesCache& fluxVarsCache)
83 CellCenterPrimaryVariables flux = ParentType::computeMassFlux(problem, element, fvGeometry,
84 elemVolVars, elemFaceVars, scvf, fluxVarsCache);
87 auto upwindTermK = [](
const auto& volVars)
89 return volVars.turbulentKineticEnergy() * volVars.density();
91 auto upwindTermEpsilon = [](
const auto& volVars)
93 return volVars.dissipation() * volVars.density();
96 flux[turbulentKineticEnergyEqIdx]
97 = ParentType::advectiveFluxForCellCenter(problem, fvGeometry, elemVolVars, elemFaceVars, scvf, upwindTermK);
98 flux[dissipationEqIdx ]
99 = ParentType::advectiveFluxForCellCenter(problem, fvGeometry, elemVolVars, elemFaceVars, scvf, upwindTermEpsilon);
102 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
103 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
104 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
105 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
108 Scalar insideCoeff_k = (insideVolVars.dynamicEddyViscosity() / insideVolVars.sigmaK()) + insideVolVars.viscosity();
109 Scalar outsideCoeff_k = (outsideVolVars.dynamicEddyViscosity() / outsideVolVars.sigmaK()) + outsideVolVars.viscosity();
110 Scalar insideCoeff_e = (insideVolVars.dynamicEddyViscosity() / insideVolVars.sigmaEpsilon()) + insideVolVars.viscosity();
111 Scalar outsideCoeff_e = (outsideVolVars.dynamicEddyViscosity() / outsideVolVars.sigmaEpsilon()) + outsideVolVars.viscosity();
114 insideCoeff_k *= insideVolVars.extrusionFactor();
115 outsideCoeff_k *= outsideVolVars.extrusionFactor();
116 insideCoeff_e *= insideVolVars.extrusionFactor();
117 outsideCoeff_e *= outsideVolVars.extrusionFactor();
119 Scalar coeff_k = 0.0;
120 Scalar coeff_e = 0.0;
124 coeff_k = insideCoeff_k;
125 coeff_e = insideCoeff_e;
126 distance = (insideScv.dofPosition() - scvf.ipGlobal()).two_norm();
134 distance = (outsideScv.dofPosition() - insideScv.dofPosition()).two_norm();
137 const auto bcTypes = problem.boundaryTypes(element, scvf);
139 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::turbulentKineticEnergyEqIdx)
140 || bcTypes.isSymmetry()
141 || bcTypes.hasWall())))
143 if (!(insideVolVars.isMatchingPoint() && outsideVolVars.isMatchingPoint())
144 || !(insideVolVars.isMatchingPoint() && outsideVolVars.inNearWallRegion())
145 || !(insideVolVars.inNearWallRegion() && outsideVolVars.isMatchingPoint()))
147 flux[turbulentKineticEnergyEqIdx]
149 * (insideVolVars.turbulentKineticEnergy() - outsideVolVars.turbulentKineticEnergy())
150 * Extrusion::area(fvGeometry, scvf);
154 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::dissipationEqIdx)
155 || bcTypes.isSymmetry())))
157 flux[dissipationEqIdx]
159 * (insideVolVars.dissipation() - outsideVolVars.dissipation())
160 * Extrusion::area(fvGeometry, scvf);
169 const Element& element,
170 const SubControlVolumeFace& scvf,
171 const FVElementGeometry& fvGeometry,
172 const ElementVolumeVariables& elemVolVars,
173 const ElementFaceVariables& elemFaceVars,
174 const GridFluxVariablesCache& gridFluxVarsCache)
176 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
178 return ParentType::computeFrontalMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
179 + ParentType::computeLateralMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
180 + 2.0 / ModelTraits::dim() * insideVolVars.density() * insideVolVars.turbulentKineticEnergy()
181 * Extrusion::area(fvGeometry, scvf) * scvf.directionSign() * insideVolVars.extrusionFactor();
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.
Definition: freeflow/rans/twoeq/kepsilon/staggered/fluxvariables.hh:75
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/rans/twoeq/kepsilon/staggered/fluxvariables.hh:168
The flux variables class for the k-epsilon model using the staggered grid discretization.
Definition: freeflow/rans/twoeq/kepsilon/fluxvariables.hh:22
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
Base class for the flux variables living on a sub control volume face.
constexpr Scalar arithmeticMean(Scalar x, Scalar y, Scalar wx=1.0, Scalar wy=1.0) noexcept
Calculate the (weighted) arithmetic mean of two scalar values.
Definition: math.hh:38
static ctype distance(const Dune::FieldVector< ctype, dimWorld > &a, const Dune::FieldVector< ctype, dimWorld > &b)
Compute the shortest distance between two points.
Definition: distance.hh:282
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The available discretization methods in Dumux.
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:166