3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/rans/twoeq/kepsilon/staggered/fluxvariables.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_KEPSILON_STAGGERED_FLUXVARIABLES_HH
25#define DUMUX_KEPSILON_STAGGERED_FLUXVARIABLES_HH
26
27#include <numeric>
33
34namespace Dumux {
35
36// forward declaration
37template<class TypeTag, class BaseFluxVariables, DiscretizationMethod discMethod>
38class KEpsilonFluxVariablesImpl;
39
44template<class TypeTag, class BaseFluxVariables>
45class KEpsilonFluxVariablesImpl<TypeTag, BaseFluxVariables, DiscretizationMethod::staggered>
46: public BaseFluxVariables
47{
48 using ParentType = BaseFluxVariables;
49
51
52 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
53 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
54 using VolumeVariables = typename GridVolumeVariables::VolumeVariables;
55
56 using GridFluxVariablesCache = typename GridVariables::GridFluxVariablesCache;
57 using FluxVariablesCache = typename GridFluxVariablesCache::FluxVariablesCache;
58
59 using GridFaceVariables = typename GridVariables::GridFaceVariables;
60 using ElementFaceVariables = typename GridFaceVariables::LocalView;
61 using FaceVariables = typename GridFaceVariables::FaceVariables;
62
67 using Element = typename GridView::template Codim<0>::Entity;
68 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
70 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
71 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
73
74 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
75 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
76
77public:
78
82 CellCenterPrimaryVariables computeMassFlux(const Problem& problem,
83 const Element &element,
84 const FVElementGeometry& fvGeometry,
85 const ElementVolumeVariables& elemVolVars,
86 const ElementFaceVariables& elemFaceVars,
87 const SubControlVolumeFace &scvf,
88 const FluxVariablesCache& fluxVarsCache)
89 {
90 CellCenterPrimaryVariables flux = ParentType::computeMassFlux(problem, element, fvGeometry,
91 elemVolVars, elemFaceVars, scvf, fluxVarsCache);
92
93 // calculate advective flux
94 auto upwindTermK = [](const auto& volVars)
95 {
96 return volVars.turbulentKineticEnergy();
97 };
98 auto upwindTermEpsilon = [](const auto& volVars)
99 {
100 return volVars.dissipation();
101 };
102
103 flux[turbulentKineticEnergyEqIdx]
104 = ParentType::advectiveFluxForCellCenter(problem, elemVolVars, elemFaceVars, scvf, upwindTermK);
105 flux[dissipationEqIdx ]
106 = ParentType::advectiveFluxForCellCenter(problem, elemVolVars, elemFaceVars, scvf, upwindTermEpsilon);
107
108 // calculate diffusive flux
109 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
110 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
111 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
112 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
113
114 // effective diffusion coefficients
115 Scalar insideCoeff_k = insideVolVars.kinematicEddyViscosity() / insideVolVars.sigmaK();
116 Scalar outsideCoeff_k = outsideVolVars.kinematicEddyViscosity() / outsideVolVars.sigmaK();
117 Scalar insideCoeff_e = insideVolVars.kinematicEddyViscosity() / insideVolVars.sigmaEpsilon();
118 Scalar outsideCoeff_e = outsideVolVars.kinematicEddyViscosity() / outsideVolVars.sigmaEpsilon();
119 static const auto kEpsilonEnableKinematicViscosity_
120 = getParamFromGroup<bool>(problem.paramGroup(), "KEpsilon.EnableKinematicViscosity", true);
121 if (kEpsilonEnableKinematicViscosity_)
122 {
123 insideCoeff_k += insideVolVars.kinematicViscosity();
124 outsideCoeff_k += outsideVolVars.kinematicViscosity();
125 insideCoeff_e += insideVolVars.kinematicViscosity();
126 outsideCoeff_e += outsideVolVars.kinematicViscosity();
127 }
128
129 // scale by extrusion factor
130 insideCoeff_k *= insideVolVars.extrusionFactor();
131 outsideCoeff_k *= outsideVolVars.extrusionFactor();
132 insideCoeff_e *= insideVolVars.extrusionFactor();
133 outsideCoeff_e *= outsideVolVars.extrusionFactor();
134
135 // average and distance
136 // is more stable with simple/unweighted arithmetic mean
137 Scalar coeff_k = arithmeticMean(insideCoeff_k, outsideCoeff_k);
138 Scalar coeff_e = arithmeticMean(insideCoeff_e, outsideCoeff_e);
139 Scalar distance = 0.0;
140
141 // adapt boundary handling
142 if (scvf.boundary())
143 {
144 distance = (insideScv.dofPosition() - scvf.ipGlobal()).two_norm();
145 coeff_k = insideCoeff_k;
146 coeff_e = insideCoeff_e;
147 }
148 else
149 {
150 distance = (outsideScv.dofPosition() - insideScv.dofPosition()).two_norm();
151 }
152
153 const auto bcTypes = problem.boundaryTypes(element, scvf);
154 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::turbulentKineticEnergyEqIdx)
155 || bcTypes.isSymmetry()
156 || problem.isOnWall(scvf))))
157 {
158 if (!(insideVolVars.isMatchingPoint() && outsideVolVars.isMatchingPoint())
159 || !(insideVolVars.isMatchingPoint() && outsideVolVars.inNearWallRegion())
160 || !(insideVolVars.inNearWallRegion() && outsideVolVars.isMatchingPoint()))
161 {
162 flux[turbulentKineticEnergyEqIdx]
163 += coeff_k / distance
164 * (insideVolVars.turbulentKineticEnergy() - outsideVolVars.turbulentKineticEnergy())
165 * scvf.area();
166 }
167 }
168 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::dissipationEqIdx)
169 || bcTypes.isSymmetry())))
170 {
171 flux[dissipationEqIdx]
172 += coeff_e / distance
173 * (insideVolVars.dissipation() - outsideVolVars.dissipation())
174 * scvf.area();
175 }
176 return flux;
177 }
178
182 FacePrimaryVariables computeMomentumFlux(const Problem& problem,
183 const Element& element,
184 const SubControlVolumeFace& scvf,
185 const FVElementGeometry& fvGeometry,
186 const ElementVolumeVariables& elemVolVars,
187 const ElementFaceVariables& elemFaceVars,
188 const GridFluxVariablesCache& gridFluxVarsCache)
189 {
190 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
191
192 return ParentType::computeFrontalMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
193 + ParentType::computeLateralMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
194 + 2.0 / ModelTraits::dim() * insideVolVars.density() * insideVolVars.turbulentKineticEnergy()
195 * scvf.area() * scvf.directionSign() * insideVolVars.extrusionFactor();
196 }
197};
198
199} // end namespace
200
201#endif
The available discretization methods in Dumux.
Base class for the flux variables living on a sub control volume face.
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
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:49
make the local view function available whenever we use the grid geometry
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:149
Definition: freeflow/rans/twoeq/kepsilon/fluxvariables.hh:34
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:82
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:182
Declares all properties used in Dumux.