version 3.9.0
Loading...
Searching...
No Matches
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// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
12#ifndef DUMUX_KEPSILON_PROBLEM_HH
13#define DUMUX_KEPSILON_PROBLEM_HH
14
15#include <numeric>
16
24
25#include "model.hh"
26
27namespace Dumux {
28
35template<class TypeTag>
36class RANSProblemImpl<TypeTag, TurbulenceModel::kepsilon> : public RANSProblemBase<TypeTag>
37{
39 using Implementation = GetPropType<TypeTag, Properties::Problem>;
41
43 using GridView = typename GridGeometry::GridView;
44 using Element = typename GridView::template Codim<0>::Entity;
45
47 using GridFaceVariables = typename GridVariables::GridFaceVariables;
48 using ElementFaceVariables = typename GridFaceVariables::LocalView;
49 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
50 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
51
52 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
53 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
54 using GlobalPosition = typename SubControlVolumeFace::GlobalPosition;
55
58 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
62 using Indices = typename ModelTraits::Indices;
63
64 static constexpr bool enableEnergyBalance = ModelTraits::enableEnergyBalance();
65 static constexpr bool isCompositional = ModelTraits::numFluidComponents() > 1;
66
67 // account for the offset of the cell center privars within the PrimaryVariables container
68 static constexpr auto cellCenterOffset = ModelTraits::numEq() - CellCenterPrimaryVariables::dimension;
69 static_assert(cellCenterOffset == ModelTraits::dim(), "cellCenterOffset must equal dim for staggered NavierStokes");
70
71public:
72
74 RANSProblemImpl(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
75 : ParentType(gridGeometry, paramGroup)
76 { }
77
82 {
83 if (!ParentType::isFlatWallBounded())
84 {
85 DUNE_THROW(Dune::NotImplemented, "\n Due to grid/geometric concerns, k-epsilon models should only be used for "
86 << " wall bounded flows with flat channel geometries. "
87 << "\n If your geometry is a flat channel, please set the runtime parameter RANS.IsFlatWallBounded to true. \n");
88 }
89
90 ParentType::updateStaticWallProperties();
91 // update size and initial values of the global vectors
92 matchingPointIdx_.resize(this->gridGeometry().elementMapper().size(), 0);
93 storedDissipation_.resize(this->gridGeometry().elementMapper().size(), 0.0);
94 storedTurbulentKineticEnergy_.resize(this->gridGeometry().elementMapper().size(), 0.0);
95 storedDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
96 zeroEqDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
97 }
98
104 template<class SolutionVector>
105 void updateDynamicWallProperties(const SolutionVector& curSol)
106 {
107 ParentType::updateDynamicWallProperties(curSol);
108
109 // update the stored eddy viscosities
110 auto fvGeometry = localView(this->gridGeometry());
111 for (const auto& element : elements(this->gridGeometry().gridView()))
112 {
113 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
114 fvGeometry.bindElement(element);
115 for (auto&& scv : scvs(fvGeometry))
116 {
117 const int dofIdx = scv.dofIndex();
118 const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
119 PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
120 auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
121 // NOTE: first update the turbulence quantities
122 storedDissipation_[elementIdx] = elemSol[0][Indices::dissipationEqIdx];
123 storedTurbulentKineticEnergy_[elementIdx] = elemSol[0][Indices::turbulentKineticEnergyEqIdx];
124 // NOTE: then update the volVars
125 VolumeVariables volVars;
126 volVars.update(elemSol, asImp_(), element, scv);
127 storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscosity();
128 }
129 }
130
131 // get matching point for k-epsilon wall function
132 unsigned int numElementsInNearWallRegion = 0;
133 for (const auto& element : elements(this->gridGeometry().gridView()))
134 {
135 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
136 unsigned int wallNormalAxis = asImp_().wallNormalAxis(elementIdx);
137 unsigned int neighborIndex0 = asImp_().neighborIndex(elementIdx, wallNormalAxis, 0);
138 unsigned int neighborIndex1 = asImp_().neighborIndex(elementIdx, wallNormalAxis, 1);
139 numElementsInNearWallRegion = inNearWallRegion(elementIdx)
140 ? numElementsInNearWallRegion + 1
141 : numElementsInNearWallRegion + 0;
142 if ((!inNearWallRegion(elementIdx) && (inNearWallRegion(neighborIndex0) || inNearWallRegion(neighborIndex1)))
143 || (!inNearWallRegion(elementIdx) && elementIdx == asImp_().wallElementIndex(elementIdx))
144 || (inNearWallRegion(elementIdx) && (asImp_().wallElementIndex(neighborIndex0) != asImp_().wallElementIndex(neighborIndex1))))
145 {
146 matchingPointIdx_[asImp_().wallElementIndex(elementIdx)] = elementIdx;
147 }
148 }
149 std::cout << "numElementsInNearWallRegion: " << numElementsInNearWallRegion << std::endl;
150
151 // calculate the potential zeroeq eddy viscosities for two-layer model
152 for (const auto& element : elements(this->gridGeometry().gridView()))
153 {
154 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
155 zeroEqDynamicEddyViscosity_[elementIdx] = zeroEqEddyViscosityModel(elementIdx);
156 }
157
158 // then make them match at the matching point
159 static const auto enableZeroEqScaling
160 = getParamFromGroup<bool>(this->paramGroup(), "KEpsilon.EnableZeroEqScaling", true);
161 if (enableZeroEqScaling)
162 {
163 for (const auto& element : elements(this->gridGeometry().gridView()))
164 {
165 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
166 unsigned int matchingPointIndex = matchingPointIdx(asImp_().wallElementIndex(elementIdx));
167
168 Scalar scalingFactor = storedDynamicEddyViscosity(matchingPointIndex)
169 / zeroEqDynamicEddyViscosity_[matchingPointIndex];
170 if (!isMatchingPoint(elementIdx)
171 && !std::isnan(scalingFactor) && !std::isinf(scalingFactor))
172 {
173 zeroEqDynamicEddyViscosity_[elementIdx] *= scalingFactor;
174 }
175 }
176 for (const auto& element : elements(this->gridGeometry().gridView()))
177 {
178 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
179 unsigned int matchingPointIndex = matchingPointIdx(asImp_().wallElementIndex(elementIdx));
180 if (isMatchingPoint(elementIdx))
181 {
182 zeroEqDynamicEddyViscosity_[matchingPointIndex] = storedDynamicEddyViscosity(matchingPointIndex);
183 }
184 }
185 }
186 }
187
191 bool inNearWallRegion(unsigned int elementIdx) const
192 {
193 unsigned int wallElementIdx = asImp_().wallElementIndex(elementIdx);
194 unsigned int matchingPointIndex = matchingPointIdx(wallElementIdx);
195 return (wallElementIdx == matchingPointIndex) ? yPlusNominal(elementIdx) < yPlusThreshold()
196 : yPlus(elementIdx) < yPlusThreshold();
197 }
198
202 bool isMatchingPoint(unsigned int elementIdx) const
203 { return matchingPointIdx(asImp_().wallElementIndex(elementIdx)) == elementIdx; }
204
208 const Scalar yPlus(unsigned int elementIdx) const
209 {
210 return asImp_().wallDistance(elementIdx) * uStar(elementIdx)
211 / asImp_().kinematicViscosity(elementIdx);
212 }
216 const Scalar yPlusNominal(unsigned int elementIdx) const
217 {
218 return asImp_().wallDistance(elementIdx) * uStarNominal(elementIdx)
219 / asImp_().kinematicViscosity(elementIdx);
220 }
221
225 const Scalar zeroEqEddyViscosityModel(unsigned int elementIdx) const
226 {
227 using std::abs;
228 using std::exp;
229 using std::sqrt;
230
231 // use VanDriest's model
232 Scalar yPlusValue = yPlus(elementIdx);
233 Scalar mixingLength = 0.0;
234 if (yPlusValue > 0.0)
235 {
236 mixingLength = asImp_().karmanConstant() * asImp_().wallDistance(elementIdx)
237 * (1.0 - exp(-yPlusValue / 26.0 ))
238 / sqrt(1.0 - exp(-0.26 * yPlusValue));
239 }
240
241 unsigned int wallNormalAxis = asImp_().wallNormalAxis(elementIdx);
242 unsigned int flowDirectionAxis = asImp_().flowDirectionAxis(elementIdx);
243 Scalar velocityGradient = asImp_().velocityGradient(elementIdx, flowDirectionAxis, wallNormalAxis);
244 return mixingLength * mixingLength * abs(velocityGradient) * asImp_().storedDensity(elementIdx);
245 }
246
248 const Scalar uStar(unsigned int elementIdx) const
249 {
250 using std::abs;
251 using std::sqrt;
252 unsigned int wallElementIdx = asImp_().wallElementIndex(elementIdx);
253 unsigned int wallNormalAxis = asImp_().wallNormalAxis(elementIdx);
254 unsigned int flowDirectionAxis = asImp_().flowDirectionAxis(elementIdx);
255 return sqrt(asImp_().kinematicViscosity(wallElementIdx)
256 * abs(asImp_().velocityGradient(wallElementIdx, flowDirectionAxis, wallNormalAxis)));
257 }
258
260 const Scalar uStarNominal(unsigned int elementIdx) const
261 {
262 using std::pow;
263 using std::sqrt;
264 unsigned int matchingPointIndex = matchingPointIdx(asImp_().wallElementIndex(elementIdx));
265 return pow(cMu(), 0.25) * sqrt(storedTurbulentKineticEnergy(matchingPointIndex));
266 }
267
271 const Scalar dissipationWallFunction(unsigned int elementIdx) const
272 {
273 return uStarNominal(elementIdx) * uStarNominal(elementIdx) * uStarNominal(elementIdx)
274 / asImp_().karmanConstant() / asImp_().wallDistance(elementIdx);
275 }
276
280 const Scalar turbulentKineticEnergyWallFunction(unsigned int elementIdx) const
281 {
282 unsigned int matchingPointIndex = matchingPointIdx(asImp_().wallElementIndex(elementIdx));
283 return storedTurbulentKineticEnergy(matchingPointIndex);
284 }
285
287 const Scalar tangentialMomentumWallFunction(unsigned int elementIdx, Scalar velocity) const
288 {
289 using std::log;
290 Scalar velocityNominal = uStarNominal(elementIdx) * (1.0 / asImp_().karmanConstant() * log(yPlusNominal(elementIdx)) + 5.0);
291 return uStarNominal(elementIdx) * uStarNominal(elementIdx) * velocity / velocityNominal;
292 }
293
295 bool useWallFunction(const Element& element,
296 const SubControlVolumeFace& localSubFace,
297 const int& eqIdx) const
298 {
299 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
300 auto bcTypes = asImp_().boundaryTypes(element, localSubFace);
301
302 return bcTypes.hasWall()
303 && bcTypes.isDirichlet(eqIdx)
304 && isMatchingPoint(elementIdx);
305 }
306
308 FacePrimaryVariables wallFunction(const Element& element,
309 const FVElementGeometry& fvGeometry,
310 const ElementVolumeVariables& elemVolVars,
311 const ElementFaceVariables& elemFaceVars,
312 const SubControlVolumeFace& scvf,
313 const SubControlVolumeFace& localSubFace) const
314 {
315 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
316 return FacePrimaryVariables(asImp_().tangentialMomentumWallFunction(elementIdx, elemFaceVars[scvf].velocitySelf())
317 * asImp_().storedDensity(elementIdx) );
318 }
319
321 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
322 typename std::enable_if_t<eB && compositional, int> = 0>
323 CellCenterPrimaryVariables wallFunction(const Element& element,
324 const FVElementGeometry& fvGeometry,
325 const ElementVolumeVariables& elemVolVars,
326 const ElementFaceVariables& elemFaceVars,
327 const SubControlVolumeFace& scvf) const
328 {
329 return wallFunctionComponent(element, fvGeometry, elemVolVars, elemFaceVars, scvf)
330 + wallFunctionEnergy(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
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 { return wallFunctionComponent(element, fvGeometry, elemVolVars, elemFaceVars, scvf); }
342
344 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
345 typename std::enable_if_t<eB && !compositional, int> = 0>
346 CellCenterPrimaryVariables wallFunction(const Element& element,
347 const FVElementGeometry& fvGeometry,
348 const ElementVolumeVariables& elemVolVars,
349 const ElementFaceVariables& elemFaceVars,
350 const SubControlVolumeFace& scvf) const
351 { return wallFunctionEnergy(element, fvGeometry, elemVolVars, elemFaceVars, scvf); }
352
354 template<bool eB = enableEnergyBalance, bool compositional = isCompositional,
355 typename std::enable_if_t<!eB && !compositional, int> = 0>
356 CellCenterPrimaryVariables wallFunction(const Element& element,
357 const FVElementGeometry& fvGeometry,
358 const ElementVolumeVariables& elemVolVars,
359 const ElementFaceVariables& elemFaceVars,
360 const SubControlVolumeFace& scvf) const
361 { return CellCenterPrimaryVariables(0.0); }
362
364 CellCenterPrimaryVariables wallFunctionComponent(const Element& element,
365 const FVElementGeometry& fvGeometry,
366 const ElementVolumeVariables& elemVolVars,
367 const ElementFaceVariables& elemFaceVars,
368 const SubControlVolumeFace& scvf) const
369 {
370 using std::log;
371 auto wallFunctionFlux = CellCenterPrimaryVariables(0.0);
372 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
373
374 // component mass fluxes
375 for (int compIdx = 0; compIdx < ModelTraits::numFluidComponents(); ++compIdx)
376 {
377 if (ModelTraits::replaceCompEqIdx() == compIdx)
378 continue;
379
380 Scalar schmidtNumber = elemVolVars[scvf.insideScvIdx()].kinematicViscosity()
381 / elemVolVars[scvf.insideScvIdx()].diffusionCoefficient(0, 0, compIdx);
382 Scalar moleToMassConversionFactor = ModelTraits::useMoles()
383 ? 1.0 : FluidSystem::molarMass(compIdx);
384 wallFunctionFlux[compIdx] +=
385 -1.0 * (asImp_().dirichlet(element, scvf)[Indices::conti0EqIdx + compIdx]
386 - elemVolVars[scvf.insideScvIdx()].moleFraction(compIdx))
387 * elemVolVars[scvf.insideScvIdx()].molarDensity()
388 * moleToMassConversionFactor
389 * uStarNominal(elementIdx)
390 / asImp_().turbulentSchmidtNumber()
391 / (1. / asImp_().karmanConstant() * log(yPlusNominal(elementIdx) * 9.793)
392 + pFunction(schmidtNumber, asImp_().turbulentSchmidtNumber()));
393 }
394
395 if (ModelTraits::replaceCompEqIdx() < ModelTraits::numFluidComponents())
396 {
397 wallFunctionFlux[ModelTraits::replaceCompEqIdx()] =
398 -std::accumulate(wallFunctionFlux.begin(), wallFunctionFlux.end(), 0.0);
399 }
400
401 return wallFunctionFlux;
402 }
403
405 CellCenterPrimaryVariables wallFunctionEnergy(const Element& element,
406 const FVElementGeometry& fvGeometry,
407 const ElementVolumeVariables& elemVolVars,
408 const ElementFaceVariables& elemFaceVars,
409 const SubControlVolumeFace& scvf) const
410 {
411 using std::log;
412 auto wallFunctionFlux = CellCenterPrimaryVariables(0.0);
413 unsigned int elementIdx = asImp_().gridGeometry().elementMapper().index(element);
414 // energy fluxes
415 Scalar prandtlNumber = elemVolVars[scvf.insideScvIdx()].kinematicViscosity()
416 * elemVolVars[scvf.insideScvIdx()].density()
417 * elemVolVars[scvf.insideScvIdx()].heatCapacity()
418 / elemVolVars[scvf.insideScvIdx()].thermalConductivity();
419 wallFunctionFlux[Indices::energyEqIdx - cellCenterOffset] +=
420 -1.0 * (asImp_().dirichlet(element, scvf)[Indices::temperatureIdx]
421 - elemVolVars[scvf.insideScvIdx()].temperature())
422 * elemVolVars[scvf.insideScvIdx()].density()
423 * elemVolVars[scvf.insideScvIdx()].heatCapacity()
424 * uStarNominal(elementIdx)
425 / asImp_().turbulentPrandtlNumber()
426 / (1. / asImp_().karmanConstant() * log(yPlusNominal(elementIdx) * 9.793)
427 + pFunction(prandtlNumber, asImp_().turbulentPrandtlNumber()));
428
429 return wallFunctionFlux;
430 }
431
433 const Scalar pFunction(Scalar molecularNumber, Scalar turbulentNumber) const
434 {
435 using std::pow;
436 using std::exp;
437 return 9.24
438 * (pow(molecularNumber / turbulentNumber, 0.75) - 1.0)
439 * (1.0 + 0.28 * exp(-0.007 * molecularNumber / turbulentNumber));
440 }
441
443 const Scalar cMu() const
444 { return 0.09; }
445
446 Scalar yPlusThreshold() const
447 {
448 static const Scalar yPlusThreshold = getParamFromGroup<Scalar>(this->paramGroup(), "KEpsilon.YPlusThreshold", 30.0);
449 return yPlusThreshold;
450 }
451
453 {
454 static const bool useStoredEddyViscosity = getParamFromGroup<bool>(this->paramGroup(), "RANS.UseStoredEddyViscosity", false);
455 return useStoredEddyViscosity;
456 }
457
458 Scalar storedDissipation(const int elementIdx) const
459 { return storedDissipation_[elementIdx]; }
460
461 Scalar storedTurbulentKineticEnergy(const int elementIdx) const
462 { return storedTurbulentKineticEnergy_[elementIdx]; }
463
464 Scalar storedDynamicEddyViscosity(const int elementIdx) const
465 { return storedDynamicEddyViscosity_[elementIdx]; }
466
467 Scalar zeroEqDynamicEddyViscosity(const int elementIdx) const
468 { return zeroEqDynamicEddyViscosity_[elementIdx]; }
469
470 unsigned int matchingPointIdx(const int elementIdx) const
471 { return matchingPointIdx_[elementIdx]; }
472
473private:
474 std::vector<unsigned int> matchingPointIdx_;
475 std::vector<Scalar> storedDissipation_;
476 std::vector<Scalar> storedTurbulentKineticEnergy_;
477 std::vector<Scalar> storedDynamicEddyViscosity_;
478 std::vector<Scalar> zeroEqDynamicEddyViscosity_;
479
481 Implementation &asImp_()
482 { return *static_cast<Implementation *>(this); }
483
485 const Implementation &asImp_() const
486 { return *static_cast<const Implementation *>(this); }
487};
488
489} // end namespace Dumux
490
491#endif
Reynolds-Averaged Navier-Stokes problem base class.
Definition freeflow/rans/problem.hh:47
const Scalar tangentialMomentumWallFunction(unsigned int elementIdx, Scalar velocity) const
Returns the nominal wall shear stress (accounts for poor approximation of viscous sublayer)
Definition freeflow/rans/twoeq/kepsilon/problem.hh:287
Scalar yPlusThreshold() const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:446
void updateDynamicWallProperties(const SolutionVector &curSol)
Update the dynamic (solution dependent) relations to the walls.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:105
const Scalar yPlus(unsigned int elementIdx) const
Returns the value at an element center.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:208
bool useStoredEddyViscosity() const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:452
unsigned int matchingPointIdx(const int elementIdx) const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:470
const Scalar zeroEqEddyViscosityModel(unsigned int elementIdx) const
Returns the kinematic eddy viscosity of a 0-Eq. model.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:225
const Scalar cMu() const
Returns the constant.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:443
bool inNearWallRegion(unsigned int elementIdx) const
Returns if an element is located in the near-wall region.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:191
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 freeflow/rans/twoeq/kepsilon/problem.hh:323
Scalar storedTurbulentKineticEnergy(const int elementIdx) const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:461
Scalar storedDissipation(const int elementIdx) const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:458
Scalar zeroEqDynamicEddyViscosity(const int elementIdx) const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:467
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 freeflow/rans/twoeq/kepsilon/problem.hh:405
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 freeflow/rans/twoeq/kepsilon/problem.hh:308
const Scalar uStarNominal(unsigned int elementIdx) const
Returns the nominal wall shear stress velocity (accounts for poor approximation of viscous sublayer)
Definition freeflow/rans/twoeq/kepsilon/problem.hh:260
const Scalar uStar(unsigned int elementIdx) const
Returns the wall shear stress velocity.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:248
const Scalar dissipationWallFunction(unsigned int elementIdx) const
Returns the dissipation calculated from the wall function consideration.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:271
RANSProblemImpl(std::shared_ptr< const GridGeometry > gridGeometry, const std::string &paramGroup="")
The constructor sets the gravity, if desired by the user.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:74
const Scalar yPlusNominal(unsigned int elementIdx) const
Returns the nominal value at an element center.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:216
Scalar storedDynamicEddyViscosity(const int elementIdx) const
Definition freeflow/rans/twoeq/kepsilon/problem.hh:464
const Scalar pFunction(Scalar molecularNumber, Scalar turbulentNumber) const
Returns the value of the P-function after Jayatilleke Versteeg2009a.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:433
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 freeflow/rans/twoeq/kepsilon/problem.hh:364
void updateStaticWallProperties()
Correct size of the static (solution independent) wall variables.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:81
bool useWallFunction(const Element &element, const SubControlVolumeFace &localSubFace, const int &eqIdx) const
Checks whether a wall function should be used.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:295
bool isMatchingPoint(unsigned int elementIdx) const
Returns if an element is the matching point.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:202
const Scalar turbulentKineticEnergyWallFunction(unsigned int elementIdx) const
Returns the turbulentKineticEnergy calculated from the wall function consideration.
Definition freeflow/rans/twoeq/kepsilon/problem.hh:280
forward declare
Definition freeflow/rans/problem.hh:31
Defines all properties used in Dumux.
the turbulence-model-specfic RANS problem
A single-phase, isothermal k-epsilon model.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethods::cctpfa||GridGeometry::discMethod==DiscretizationMethods::ccmpfa, CCElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for cell-centered schemes.
Definition cellcentered/elementsolution.hh:101
TurbulenceModel
The available free flow turbulence models in Dumux.
Definition turbulencemodel.hh:26
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
PrimaryVariables makePriVarsFromCellCenterPriVars(const CellCenterPrimaryVariables &cellCenterPriVars)
Helper function to create a PrimaryVariables object from CellCenterPrimaryVariables.
Definition staggered/elementsolution.hh:29
Free function to get the local view of a grid cache object.
The available discretization methods in Dumux.
Definition adapt.hh:17
The local element solution class for staggered methods.
Base class for all staggered fv problems.
The available free flow turbulence models in Dumux.