3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
dumux/freeflow/rans/twoeq/kepsilon/problem.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_PROBLEM_HH
25#define DUMUX_KEPSILON_PROBLEM_HH
26
27#include <numeric>
28
36
37#include "model.hh"
38
39namespace Dumux {
40
47template<class TypeTag>
48class RANSProblemImpl<TypeTag, TurbulenceModel::kepsilon> : public RANSProblemBase<TypeTag>
49{
50 using ParentType = RANSProblemBase<TypeTag>;
51 using Implementation = GetPropType<TypeTag, Properties::Problem>;
53
55 using GridView = typename GridGeometry::GridView;
56 using Element = typename GridView::template Codim<0>::Entity;
57
59 using GridFaceVariables = typename GridVariables::GridFaceVariables;
60 using ElementFaceVariables = typename GridFaceVariables::LocalView;
61 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
62 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
63
64 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
65 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
66
70 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
74 using Indices = typename ModelTraits::Indices;
75
76 static constexpr bool enableEnergyBalance = ModelTraits::enableEnergyBalance();
77 static constexpr bool isCompositional = ModelTraits::numFluidComponents() > 1;
78
79 // account for the offset of the cell center privars within the PrimaryVariables container
80 static constexpr auto cellCenterOffset = ModelTraits::numEq() - CellCenterPrimaryVariables::dimension;
81 static_assert(cellCenterOffset == ModelTraits::dim(), "cellCenterOffset must equal dim for staggered NavierStokes");
82
83public:
84
86 RANSProblemImpl(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
87 : ParentType(gridGeometry, paramGroup)
88 {
89 yPlusThreshold_ = getParamFromGroup<Scalar>(this->paramGroup(), "KEpsilon.YPlusThreshold", 30);
90 useStoredEddyViscosity_ = getParamFromGroup<bool>(this->paramGroup(), "RANS.UseStoredEddyViscosity", false);
91 }
92
97 {
98 ParentType::updateStaticWallProperties();
99
100 // update size and initial values of the global vectors
101 matchingPointIdx_.resize(this->gridGeometry().elementMapper().size(), 0);
102 storedDensity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
103 storedDissipation_.resize(this->gridGeometry().elementMapper().size(), 0.0);
104 storedTurbulentKineticEnergy_.resize(this->gridGeometry().elementMapper().size(), 0.0);
105 storedDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
106 zeroEqDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
107 }
108
114 void updateDynamicWallProperties(const SolutionVector& curSol)
115 {
116 ParentType::updateDynamicWallProperties(curSol);
117
118 // update the stored eddy viscosities
119 for (const auto& element : elements(this->gridGeometry().gridView()))
120 {
121 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
122
123 auto fvGeometry = localView(this->gridGeometry());
124 fvGeometry.bindElement(element);
125 for (auto&& scv : scvs(fvGeometry))
126 {
127 const int dofIdx = scv.dofIndex();
128 const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
129 PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
130 auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
131 // NOTE: first update the turbulence quantities
132 storedDissipation_[elementIdx] = elemSol[0][Indices::dissipationEqIdx];
133 storedTurbulentKineticEnergy_[elementIdx] = elemSol[0][Indices::turbulentKineticEnergyEqIdx];
134 // NOTE: then update the volVars
135 VolumeVariables volVars;
136 volVars.update(elemSol, asImp_(), element, scv);
137 storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscosity();
138 storedDensity_[elementIdx] = volVars.density();
139 }
140 }
141
142 // get matching point for k-epsilon wall function
143 unsigned int numElementsInNearWallRegion = 0;
144 for (const auto& element : elements(this->gridGeometry().gridView()))
145 {
146 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
147 unsigned int wallNormalAxis = asImp_().wallNormalAxis_[elementIdx];
148 unsigned int neighborIdx0 = asImp_().neighborIdx_[elementIdx][wallNormalAxis][0];
149 unsigned int neighborIdx1 = asImp_().neighborIdx_[elementIdx][wallNormalAxis][1];
150 numElementsInNearWallRegion = inNearWallRegion(elementIdx)
151 ? numElementsInNearWallRegion + 1
152 : numElementsInNearWallRegion + 0;
153 if ((!inNearWallRegion(elementIdx) && (inNearWallRegion(neighborIdx0) || inNearWallRegion(neighborIdx1)))
154 || (!inNearWallRegion(elementIdx) && elementIdx == asImp_().wallElementIdx_[elementIdx])
155 || (inNearWallRegion(elementIdx) && (asImp_().wallElementIdx_[neighborIdx0] != asImp_().wallElementIdx_[neighborIdx1])))
156 {
157 matchingPointIdx_[asImp_().wallElementIdx_[elementIdx]] = elementIdx;
158 }
159 }
160 std::cout << "numElementsInNearWallRegion: " << numElementsInNearWallRegion << std::endl;
161
162 // calculate the potential zeroeq eddy viscosities for two-layer model
163 for (const auto& element : elements(this->gridGeometry().gridView()))
164 {
165 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
166 zeroEqDynamicEddyViscosity_[elementIdx] = zeroEqEddyViscosityModel(elementIdx);
167 }
168
169 // then make them match at the matching point
170 static const auto enableZeroEqScaling
171 = getParamFromGroup<bool>(this->paramGroup(), "KEpsilon.EnableZeroEqScaling", true);
172 if (enableZeroEqScaling)
173 {
174 for (const auto& element : elements(this->gridGeometry().gridView()))
175 {
176 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
177 unsigned int matchingPointIdx = matchingPointIdx_[asImp_().wallElementIdx_[elementIdx]];
178
179 Scalar scalingFactor = storedDynamicEddyViscosity_[matchingPointIdx]
180 / zeroEqDynamicEddyViscosity_[matchingPointIdx];
181 if (!isMatchingPoint(elementIdx)
182 && !std::isnan(scalingFactor) && !std::isinf(scalingFactor))
183 {
184 zeroEqDynamicEddyViscosity_[elementIdx] *= scalingFactor;
185 }
186 }
187 for (const auto& element : elements(this->gridGeometry().gridView()))
188 {
189 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
190 unsigned int matchingPointIdx = matchingPointIdx_[asImp_().wallElementIdx_[elementIdx]];
191 if (isMatchingPoint(elementIdx))
192 {
193 zeroEqDynamicEddyViscosity_[matchingPointIdx] = storedDynamicEddyViscosity_[matchingPointIdx];
194 }
195 }
196 }
197 }
198
202 const bool inNearWallRegion(unsigned int elementIdx) const
203 {
204 unsigned int wallElementIdx = asImp_().wallElementIdx_[elementIdx];
205 unsigned int matchingPointIdx = matchingPointIdx_[wallElementIdx];
206 return wallElementIdx == matchingPointIdx ? yPlusNominal(elementIdx) < yPlusThreshold_
207 : yPlus(elementIdx) < yPlusThreshold_;
208 }
209
213 const bool isMatchingPoint(unsigned int elementIdx) const
214 { return matchingPointIdx_[asImp_().wallElementIdx_[elementIdx]] == elementIdx; }
215
219 const Scalar yPlus(unsigned int elementIdx) const
220 {
221 return asImp_().wallDistance_[elementIdx] * uStar(elementIdx)
222 / asImp_().kinematicViscosity_[elementIdx];
223 }
227 const Scalar yPlusNominal(unsigned int elementIdx) const
228 {
229 return asImp_().wallDistance_[elementIdx] * uStarNominal(elementIdx)
230 / asImp_().kinematicViscosity_[elementIdx];
231 }
232
236 const Scalar zeroEqEddyViscosityModel(unsigned int elementIdx) const
237 {
238 using std::abs;
239 using std::exp;
240 using std::sqrt;
241
242 // use VanDriest's model
243 Scalar yPlusValue = yPlus(elementIdx);
244 Scalar mixingLength = 0.0;
245 if (yPlusValue > 0.0)
246 {
247 mixingLength = asImp_().karmanConstant() * asImp_().wallDistance_[elementIdx]
248 * (1.0 - exp(-yPlusValue / 26.0 ))
249 / sqrt(1.0 - exp(-0.26 * yPlusValue));
250 }
251
252 unsigned int wallNormalAxis = asImp_().wallNormalAxis_[elementIdx];
253 unsigned int flowNormalAxis = asImp_().flowNormalAxis_[elementIdx];
254 Scalar velocityGradient = asImp_().velocityGradients_[elementIdx][flowNormalAxis][wallNormalAxis];
255 return mixingLength * mixingLength * abs(velocityGradient) * storedDensity_[elementIdx];
256 }
257
259 const Scalar uStar(unsigned int elementIdx) const
260 {
261 using std::abs;
262 using std::sqrt;
263 unsigned int wallElementIdx = asImp_().wallElementIdx_[elementIdx];
264 unsigned int wallNormalAxis = asImp_().wallNormalAxis_[elementIdx];
265 unsigned int flowNormalAxis = asImp_().flowNormalAxis_[elementIdx];
266 return sqrt(asImp_().kinematicViscosity_[wallElementIdx]
267 * abs(asImp_().velocityGradients_[wallElementIdx][flowNormalAxis][wallNormalAxis]));
268 }
269
271 const Scalar uStarNominal(unsigned int elementIdx) const
272 {
273 using std::pow;
274 using std::sqrt;
275 unsigned int matchingPointIdx = matchingPointIdx_[asImp_().wallElementIdx_[elementIdx]];
276 return pow(cMu(), 0.25)
277 * sqrt(storedTurbulentKineticEnergy_[matchingPointIdx]);
278 }
279
283 const Scalar dissipationWallFunction(unsigned int elementIdx) const
284 {
285 return uStarNominal(elementIdx) * uStarNominal(elementIdx) * uStarNominal(elementIdx)
286 / asImp_().karmanConstant() / asImp_().wallDistance_[elementIdx];
287 }
288
292 const Scalar turbulentKineticEnergyWallFunction(unsigned int elementIdx) const
293 {
294 unsigned int wallElementIdx = asImp_().wallElementIdx_[elementIdx];
295 unsigned int matchingPointIdx = matchingPointIdx_[wallElementIdx];
296 return storedTurbulentKineticEnergy_[matchingPointIdx];
297 }
298
300 const Scalar tangentialMomentumWallFunction(unsigned int elementIdx, Scalar velocity) const
301 {
302 using std::log;
303 Scalar velocityNominal = uStarNominal(elementIdx) * (1.0 / asImp_().karmanConstant() * log(yPlusNominal(elementIdx)) + 5.0);
304 return uStarNominal(elementIdx) * uStarNominal(elementIdx)
305 * velocity / velocityNominal;
306 }
307
309 bool useWallFunction(const Element& element,
310 const SubControlVolumeFace& localSubFace,
311 const int& eqIdx) const
312 {
313 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
314 auto bcTypes = asImp_().boundaryTypes(element, localSubFace);
315 return asImp_().isOnWall(localSubFace)
316 && bcTypes.isDirichlet(eqIdx)
317 && isMatchingPoint(elementIdx);
318 }
319
321 FacePrimaryVariables wallFunction(const Element& element,
322 const FVElementGeometry& fvGeometry,
323 const ElementVolumeVariables& elemVolVars,
324 const ElementFaceVariables& elemFaceVars,
325 const SubControlVolumeFace& scvf,
326 const SubControlVolumeFace& localSubFace) const
327 {
328 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
329 return FacePrimaryVariables(asImp_().tangentialMomentumWallFunction(elementIdx, elemFaceVars[scvf].velocitySelf())
330 * elemVolVars[scvf.insideScvIdx()].density());
331 }
332
334 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
335 typename std::enable_if_t<eB && compositional, int> = 0>
336 CellCenterPrimaryVariables wallFunction(const Element& element,
337 const FVElementGeometry& fvGeometry,
338 const ElementVolumeVariables& elemVolVars,
339 const ElementFaceVariables& elemFaceVars,
340 const SubControlVolumeFace& scvf) const
341 {
342 return wallFunctionComponent(element, fvGeometry, elemVolVars, elemFaceVars, scvf)
343 + wallFunctionEnergy(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
344 }
345
347 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
348 typename std::enable_if_t<!eB && compositional, int> = 0>
349 CellCenterPrimaryVariables wallFunction(const Element& element,
350 const FVElementGeometry& fvGeometry,
351 const ElementVolumeVariables& elemVolVars,
352 const ElementFaceVariables& elemFaceVars,
353 const SubControlVolumeFace& scvf) const
354 { return wallFunctionComponent(element, fvGeometry, elemVolVars, elemFaceVars, scvf); }
355
357 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
358 typename std::enable_if_t<eB && !compositional, int> = 0>
359 CellCenterPrimaryVariables wallFunction(const Element& element,
360 const FVElementGeometry& fvGeometry,
361 const ElementVolumeVariables& elemVolVars,
362 const ElementFaceVariables& elemFaceVars,
363 const SubControlVolumeFace& scvf) const
364 { return wallFunctionEnergy(element, fvGeometry, elemVolVars, elemFaceVars, scvf); }
365
367 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
368 typename std::enable_if_t<!eB && !compositional, int> = 0>
369 CellCenterPrimaryVariables wallFunction(const Element& element,
370 const FVElementGeometry& fvGeometry,
371 const ElementVolumeVariables& elemVolVars,
372 const ElementFaceVariables& elemFaceVars,
373 const SubControlVolumeFace& scvf) const
374 { return CellCenterPrimaryVariables(0.0); }
375
377 CellCenterPrimaryVariables wallFunctionComponent(const Element& element,
378 const FVElementGeometry& fvGeometry,
379 const ElementVolumeVariables& elemVolVars,
380 const ElementFaceVariables& elemFaceVars,
381 const SubControlVolumeFace& scvf) const
382 {
383 using std::log;
384 auto wallFunctionFlux = CellCenterPrimaryVariables(0.0);
385 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
386
387 // component mass fluxes
388 for (int compIdx = 0; compIdx < ModelTraits::numFluidComponents(); ++compIdx)
389 {
390 if (ModelTraits::replaceCompEqIdx() == compIdx)
391 continue;
392
393 Scalar schmidtNumber = elemVolVars[scvf.insideScvIdx()].kinematicViscosity()
394 / elemVolVars[scvf.insideScvIdx()].diffusionCoefficient(0, compIdx);
395 Scalar moleToMassConversionFactor = ModelTraits::useMoles()
396 ? 1.0 : FluidSystem::molarMass(compIdx);
397 wallFunctionFlux[compIdx] +=
398 -1.0 * (asImp_().dirichlet(element, scvf)[Indices::conti0EqIdx + compIdx]
399 - elemVolVars[scvf.insideScvIdx()].moleFraction(compIdx))
400 * elemVolVars[scvf.insideScvIdx()].molarDensity()
401 * moleToMassConversionFactor
402 * uStarNominal(elementIdx)
403 / asImp_().turbulentSchmidtNumber()
404 / (1. / asImp_().karmanConstant() * log(yPlusNominal(elementIdx) * 9.793)
405 + pFunction(schmidtNumber, asImp_().turbulentSchmidtNumber()));
406 }
407
408 if (ModelTraits::replaceCompEqIdx() < ModelTraits::numFluidComponents())
409 {
410 wallFunctionFlux[ModelTraits::replaceCompEqIdx()] =
411 -std::accumulate(wallFunctionFlux.begin(), wallFunctionFlux.end(), 0.0);
412 }
413
414 return wallFunctionFlux;
415 }
416
418 CellCenterPrimaryVariables wallFunctionEnergy(const Element& element,
419 const FVElementGeometry& fvGeometry,
420 const ElementVolumeVariables& elemVolVars,
421 const ElementFaceVariables& elemFaceVars,
422 const SubControlVolumeFace& scvf) const
423 {
424 using std::log;
425 auto wallFunctionFlux = CellCenterPrimaryVariables(0.0);
426 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
427 // energy fluxes
428 Scalar prandtlNumber = elemVolVars[scvf.insideScvIdx()].kinematicViscosity()
429 * elemVolVars[scvf.insideScvIdx()].density()
430 * elemVolVars[scvf.insideScvIdx()].heatCapacity()
431 / elemVolVars[scvf.insideScvIdx()].thermalConductivity();
432 wallFunctionFlux[Indices::energyEqIdx - cellCenterOffset] +=
433 -1.0 * (asImp_().dirichlet(element, scvf)[Indices::temperatureIdx]
434 - elemVolVars[scvf.insideScvIdx()].temperature())
435 * elemVolVars[scvf.insideScvIdx()].density()
436 * elemVolVars[scvf.insideScvIdx()].heatCapacity()
437 * uStarNominal(elementIdx)
438 / asImp_().turbulentPrandtlNumber()
439 / (1. / asImp_().karmanConstant() * log(yPlusNominal(elementIdx) * 9.793)
440 + pFunction(prandtlNumber, asImp_().turbulentPrandtlNumber()));
441
442 return wallFunctionFlux;
443 }
444
446 const Scalar pFunction(Scalar molecularNumber, Scalar turbulentNumber) const
447 {
448 using std::pow;
449 using std::exp;
450 return 9.24
451 * (pow(molecularNumber / turbulentNumber, 0.75) - 1.0)
452 * (1.0 + 0.28 * exp(-0.007 * molecularNumber / turbulentNumber));
453 }
454
456 const Scalar cMu() const
457 {
458 return 0.09;
459 }
460
461public:
462 std::vector<unsigned int> matchingPointIdx_;
463 std::vector<Scalar> storedDensity_;
464 std::vector<Scalar> storedDissipation_;
465 std::vector<Scalar> storedTurbulentKineticEnergy_;
466 std::vector<Scalar> storedDynamicEddyViscosity_;
467 std::vector<Scalar> zeroEqDynamicEddyViscosity_;
470
471private:
473 Implementation &asImp_()
474 { return *static_cast<Implementation *>(this); }
475
477 const Implementation &asImp_() const
478 { return *static_cast<const Implementation *>(this); }
479};
480
481} // end namespace Dumux
482
483#endif
Base class for all staggered fv problems.
Free function to get the local view of a grid cache object.
The available discretization methods in Dumux.
The available free flow turbulence models in Dumux.
TurbulenceModel
The available free flow turbulence models in Dumux.
Definition: turbulencemodel.hh:37
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
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
Property tag Indices
Definition: porousmediumflow/sequential/properties.hh:59
forward declare
Definition: dumux/freeflow/rans/problem.hh:41
Reynolds-Averaged Navier-Stokes problem base class.
Definition: dumux/freeflow/rans/problem.hh:57
const Scalar tangentialMomentumWallFunction(unsigned int elementIdx, Scalar velocity) const
Returns the nominal wall shear stress (accounts for poor approximation of viscous sublayer)
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:300
std::vector< Scalar > storedDensity_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:463
void updateDynamicWallProperties(const SolutionVector &curSol)
Update the dynamic (solution dependent) relations to the walls.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:114
bool useStoredEddyViscosity_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:468
const Scalar yPlus(unsigned int elementIdx) const
Returns the value at an element center.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:219
std::vector< Scalar > storedDynamicEddyViscosity_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:466
const Scalar zeroEqEddyViscosityModel(unsigned int elementIdx) const
Returns the kinematic eddy viscosity of a 0-Eq. model.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:236
const bool isMatchingPoint(unsigned int elementIdx) const
Returns if an element is the matching point.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:213
const Scalar cMu() const
Returns the constant.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:456
CellCenterPrimaryVariables wallFunction(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Returns the flux for non-isothermal and compositional RANS models.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:336
const bool inNearWallRegion(unsigned int elementIdx) const
Returns if an element is located in the near-wall region.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:202
CellCenterPrimaryVariables wallFunctionEnergy(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Returns the energy wall-function flux.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:418
FacePrimaryVariables wallFunction(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const SubControlVolumeFace &localSubFace) const
Returns an additional wall function momentum flux (only needed for RANS models)
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:321
const Scalar uStarNominal(unsigned int elementIdx) const
Returns the nominal wall shear stress velocity (accounts for poor approximation of viscous sublayer)
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:271
const Scalar uStar(unsigned int elementIdx) const
Returns the wall shear stress velocity.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:259
const Scalar dissipationWallFunction(unsigned int elementIdx) const
Returns the dissipation calculated from the wall function consideration.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:283
RANSProblemImpl(std::shared_ptr< const GridGeometry > gridGeometry, const std::string &paramGroup="")
The constructor sets the gravity, if desired by the user.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:86
const Scalar yPlusNominal(unsigned int elementIdx) const
Returns the nominal value at an element center.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:227
const Scalar pFunction(Scalar molecularNumber, Scalar turbulentNumber) const
Returns the value of the P-function after Jayatilleke .
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:446
CellCenterPrimaryVariables wallFunctionComponent(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Returns the component wall-function flux.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:377
std::vector< Scalar > storedTurbulentKineticEnergy_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:465
std::vector< Scalar > zeroEqDynamicEddyViscosity_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:467
void updateStaticWallProperties()
Correct size of the static (solution independent) wall variables.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:96
bool useWallFunction(const Element &element, const SubControlVolumeFace &localSubFace, const int &eqIdx) const
Checks whether a wall function should be used.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:309
std::vector< unsigned int > matchingPointIdx_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:462
std::vector< Scalar > storedDissipation_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:464
const Scalar turbulentKineticEnergyWallFunction(unsigned int elementIdx) const
Returns the turbulentKineticEnergy calculated from the wall function consideration.
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:292
Scalar yPlusThreshold_
Definition: dumux/freeflow/rans/twoeq/kepsilon/problem.hh:469
Declares all properties used in Dumux.
Adaption of the fully implicit scheme to the tracer transport model.
The local element solution class for staggered methods.