3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
dumux/freeflow/rans/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_RANS_PROBLEM_HH
25#define DUMUX_RANS_PROBLEM_HH
26
27#include <dune/common/fmatrix.hh>
34
35#include "model.hh"
36
37namespace Dumux {
38
40template<class TypeTag, TurbulenceModel turbulenceModel>
42
44template<class TypeTag>
46
55template<class TypeTag>
56class RANSProblemBase : public NavierStokesProblem<TypeTag>
57{
58 using ParentType = NavierStokesProblem<TypeTag>;
59 using Implementation = GetPropType<TypeTag, Properties::Problem>;
60
62
64 using FVElementGeometry = typename GridGeometry::LocalView;
65 using GridView = typename GridGeometry::GridView;
66 using Element = typename GridView::template Codim<0>::Entity;
67 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
68 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
71 using PrimaryVariables = typename VolumeVariables::PrimaryVariables;
72 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
75
76 using GlobalPosition = typename SubControlVolumeFace::GlobalPosition;
77
78 static constexpr auto dim = GridView::dimension;
79 using DimVector = GlobalPosition;
80 using DimMatrix = Dune::FieldMatrix<Scalar, dim, dim>;
81
82public:
88 RANSProblemBase(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
89 : ParentType(gridGeometry, paramGroup)
90 { }
91
99 {
100 using std::abs;
101 std::cout << "Update static wall properties. ";
103
104 // update size and initial values of the global vectors
105 wallElementIdx_.resize(this->gridGeometry().elementMapper().size());
106 wallDistance_.resize(this->gridGeometry().elementMapper().size(), std::numeric_limits<Scalar>::max());
107 neighborIdx_.resize(this->gridGeometry().elementMapper().size());
108 cellCenter_.resize(this->gridGeometry().elementMapper().size(), GlobalPosition(0.0));
109 velocity_.resize(this->gridGeometry().elementMapper().size(), DimVector(0.0));
110 velocityMaximum_.resize(this->gridGeometry().elementMapper().size(), DimVector(0.0));
111 velocityGradients_.resize(this->gridGeometry().elementMapper().size(), DimMatrix(0.0));
112 stressTensorScalarProduct_.resize(this->gridGeometry().elementMapper().size(), 0.0);
113 vorticityTensorScalarProduct_.resize(this->gridGeometry().elementMapper().size(), 0.0);
114 flowNormalAxis_.resize(this->gridGeometry().elementMapper().size(), 0);
115 wallNormalAxis_.resize(this->gridGeometry().elementMapper().size(), 1);
116 kinematicViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
117 sandGrainRoughness_.resize(this->gridGeometry().elementMapper().size(), 0.0);
118
119 // retrieve all wall intersections and corresponding elements
120 std::vector<unsigned int> wallElements;
121 std::vector<GlobalPosition> wallPositions;
122 std::vector<unsigned int> wallNormalAxisTemp;
123
124 const auto gridView = this->gridGeometry().gridView();
125 auto fvGeometry = localView(this->gridGeometry());
126
127 for (const auto& element : elements(gridView))
128 {
129 fvGeometry.bindElement(element);
130 for (const auto& scvf : scvfs(fvGeometry))
131 {
132 // only search for walls at a global boundary
133 if (!scvf.boundary())
134 continue;
135
136 if (asImp_().isOnWall(scvf))
137 {
138 wallElements.push_back(this->gridGeometry().elementMapper().index(element));
139 wallPositions.push_back(scvf.center());
140 wallNormalAxisTemp.push_back(scvf.directionIndex());
141 }
142 }
143 }
144 std::cout << "NumWallIntersections=" << wallPositions.size() << std::endl;
145 if (wallPositions.size() == 0)
146 DUNE_THROW(Dune::InvalidStateException,
147 "No wall intersections have been found. Make sure that the isOnWall(globalPos) is working properly.");
148
149
150 // search for shortest distance to wall for each element
151 for (const auto& element : elements(gridView))
152 {
153 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
154 cellCenter_[elementIdx] = element.geometry().center();
155 for (unsigned int i = 0; i < wallPositions.size(); ++i)
156 {
157 static const int problemWallNormalAxis
158 = getParamFromGroup<int>(this->paramGroup(), "RANS.WallNormalAxis", -1);
159 int searchAxis = problemWallNormalAxis;
160
161 // search along wall normal axis of the intersection
162 if (problemWallNormalAxis < 0 || problemWallNormalAxis >= dim)
163 {
164 searchAxis = wallNormalAxisTemp[i];
165 }
166
167 GlobalPosition global = element.geometry().center();
168 global -= wallPositions[i];
169 // second and argument ensures to use only aligned elements
170 if (abs(global[searchAxis]) < wallDistance_[elementIdx]
171 && abs(global[searchAxis]) < global.two_norm() + 1e-8
172 && abs(global[searchAxis]) > global.two_norm() - 1e-8)
173 {
174 wallDistance_[elementIdx] = abs(global[searchAxis]);
175 wallElementIdx_[elementIdx] = wallElements[i];
176 wallNormalAxis_[elementIdx] = searchAxis;
177 sandGrainRoughness_[elementIdx] = asImp_().sandGrainRoughnessAtPos(wallPositions[i]);
178 }
179 }
180 }
181
182 // search for neighbor Idxs
183 for (const auto& element : elements(gridView))
184 {
185 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
186 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
187 {
188 neighborIdx_[elementIdx][dimIdx][0] = elementIdx;
189 neighborIdx_[elementIdx][dimIdx][1] = elementIdx;
190 }
191
192 for (const auto& intersection : intersections(gridView, element))
193 {
194 if (intersection.boundary())
195 continue;
196
197 unsigned int neighborIdx = this->gridGeometry().elementMapper().index(intersection.outside());
198 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
199 {
200 if (abs(cellCenter_[elementIdx][dimIdx] - cellCenter_[neighborIdx][dimIdx]) > 1e-8)
201 {
202 if (cellCenter_[elementIdx][dimIdx] > cellCenter_[neighborIdx][dimIdx])
203 {
204 neighborIdx_[elementIdx][dimIdx][0] = neighborIdx;
205 }
206 if (cellCenter_[elementIdx][dimIdx] < cellCenter_[neighborIdx][dimIdx])
207 {
208 neighborIdx_[elementIdx][dimIdx][1] = neighborIdx;
209 }
210 }
211 }
212 }
213 }
214 }
215
225 void updateDynamicWallProperties(const SolutionVector& curSol)
226 {
227 using std::abs;
228 using std::max;
229 using std::min;
230 std::cout << "Update dynamic wall properties." << std::endl;
232 DUNE_THROW(Dune::InvalidStateException,
233 "You have to call updateStaticWallProperties() once before you call updateDynamicWallProperties().");
234
235 static const int flowNormalAxis
236 = getParamFromGroup<int>(this->paramGroup(), "RANS.FlowNormalAxis", -1);
237
238 // re-initialize min and max values
239 velocityMaximum_.assign(this->gridGeometry().elementMapper().size(), DimVector(1e-16));
240 velocityMinimum_.assign(this->gridGeometry().elementMapper().size(), DimVector(std::numeric_limits<Scalar>::max()));
241
242 // calculate cell-center-averaged velocities
243 for (const auto& element : elements(this->gridGeometry().gridView()))
244 {
245 auto fvGeometry = localView(this->gridGeometry());
246 fvGeometry.bindElement(element);
247 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
248
249 // calculate velocities
250 DimVector velocityTemp(0.0);
251 for (auto&& scvf : scvfs(fvGeometry))
252 {
253 const int dofIdxFace = scvf.dofIndex();
254 const auto numericalSolutionFace = curSol[GridGeometry::faceIdx()][dofIdxFace][Indices::velocity(scvf.directionIndex())];
255 velocityTemp[scvf.directionIndex()] += numericalSolutionFace;
256 }
257 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
258 velocity_[elementIdx][dimIdx] = velocityTemp[dimIdx] * 0.5; // faces are equidistant to cell center
259 }
260
261 // calculate cell-center-averaged velocity gradients, maximum, and minimum values
262 for (const auto& element : elements(this->gridGeometry().gridView()))
263 {
264 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
265 unsigned int wallElementIdx = wallElementIdx_[elementIdx];
266
267 Scalar maxVelocity = 0.0;
268 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
269 {
270 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
271 {
272 velocityGradients_[elementIdx][velIdx][dimIdx]
273 = (velocity_[neighborIdx_[elementIdx][dimIdx][1]][velIdx]
274 - velocity_[neighborIdx_[elementIdx][dimIdx][0]][velIdx])
275 / (cellCenter_[neighborIdx_[elementIdx][dimIdx][1]][dimIdx]
276 - cellCenter_[neighborIdx_[elementIdx][dimIdx][0]][dimIdx]);
277 if (abs(cellCenter_[neighborIdx_[elementIdx][dimIdx][1]][dimIdx]
278 - cellCenter_[neighborIdx_[elementIdx][dimIdx][0]][dimIdx]) < 1e-8)
279 velocityGradients_[elementIdx][velIdx][dimIdx] = 0.0;
280 }
281
282 if (abs(velocity_[elementIdx][dimIdx]) > abs(velocityMaximum_[wallElementIdx][dimIdx]))
283 {
284 velocityMaximum_[wallElementIdx][dimIdx] = velocity_[elementIdx][dimIdx];
285 }
286 if (abs(velocity_[elementIdx][dimIdx]) < abs(velocityMinimum_[wallElementIdx][dimIdx]))
287 {
288 velocityMinimum_[wallElementIdx][dimIdx] = velocity_[elementIdx][dimIdx];
289 }
290
291 if (0 <= flowNormalAxis && flowNormalAxis < dim)
292 {
293 flowNormalAxis_[elementIdx] = flowNormalAxis;
294 }
295 else if (abs(maxVelocity) < abs(velocity_[elementIdx][dimIdx]))
296 {
297 maxVelocity = abs(velocity_[elementIdx][dimIdx]);
298 flowNormalAxis_[elementIdx] = dimIdx;
299 }
300 }
301
302 auto fvGeometry = localView(this->gridGeometry());
303 fvGeometry.bindElement(element);
304 for (auto&& scvf : scvfs(fvGeometry))
305 {
306 // adapt calculations for Dirichlet condition
307 unsigned int scvfNormDim = scvf.directionIndex();
308 if (scvf.boundary())
309 {
310 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
311 {
312 if (!asImp_().boundaryTypes(element, scvf).isDirichlet(Indices::velocity(velIdx)))
313 continue;
314
315 Scalar dirichletVelocity = asImp_().dirichlet(element, scvf)[Indices::velocity(velIdx)];
316
317 unsigned int neighborIdx = neighborIdx_[elementIdx][scvfNormDim][0];
318 if (scvf.center()[scvfNormDim] < cellCenter_[elementIdx][scvfNormDim])
319 neighborIdx = neighborIdx_[elementIdx][scvfNormDim][1];
320
321 velocityGradients_[elementIdx][velIdx][scvfNormDim]
322 = (velocity_[neighborIdx][velIdx] - dirichletVelocity)
323 / (cellCenter_[neighborIdx][scvfNormDim] - scvf.center()[scvfNormDim]);
324 }
325 }
326
327 // Calculate the BJS-velocity by accounting for all sub faces.
328 std::vector<int> bjsNumFaces(dim, 0);
329 std::vector<unsigned int> bjsNeighbor(dim, 0);
330 DimVector bjsVelocityAverage(0.0);
331 DimVector normalNormCoordinate(0.0);
332 unsigned int velIdx = Indices::velocity(scvfNormDim);
333 const int numSubFaces = scvf.pairData().size();
334 for(int localSubFaceIdx = 0; localSubFaceIdx < numSubFaces; ++localSubFaceIdx)
335 {
336 const auto& lateralFace = fvGeometry.scvf(scvf.insideScvIdx(), scvf.pairData()[localSubFaceIdx].localLateralFaceIdx);
337
338 // adapt calculations for Beavers-Joseph-Saffman condition
339 unsigned int normalNormDim = lateralFace.directionIndex();
340 if (lateralFace.boundary() && (asImp_().boundaryTypes(element, lateralFace).isBJS(Indices::velocity(velIdx))))
341 {
342 unsigned int neighborIdx = neighborIdx_[elementIdx][normalNormDim][0];
343 if (lateralFace.center()[normalNormDim] < cellCenter_[elementIdx][normalNormDim])
344 neighborIdx = neighborIdx_[elementIdx][normalNormDim][1];
345
346 const SubControlVolume& scv = fvGeometry.scv(scvf.insideScvIdx());
347 bjsVelocityAverage[normalNormDim] += ParentType::bjsVelocity(element, scv, lateralFace, velocity_[elementIdx][velIdx]);
348 if (bjsNumFaces[normalNormDim] > 0 && neighborIdx != bjsNeighbor[normalNormDim])
349 DUNE_THROW(Dune::InvalidStateException, "Two different neighborIdx should not occur");
350 bjsNeighbor[normalNormDim] = neighborIdx;
351 normalNormCoordinate[normalNormDim] = lateralFace.center()[normalNormDim];
352 bjsNumFaces[normalNormDim]++;
353 }
354 }
355 for (unsigned dirIdx = 0; dirIdx < dim; ++dirIdx)
356 {
357 if (bjsNumFaces[dirIdx] == 0)
358 continue;
359
360 unsigned int neighborIdx = bjsNeighbor[dirIdx];
361 bjsVelocityAverage[dirIdx] /= bjsNumFaces[dirIdx];
362
363 velocityGradients_[elementIdx][velIdx][dirIdx]
364 = (velocity_[neighborIdx][velIdx] - bjsVelocityAverage[dirIdx])
365 / (cellCenter_[neighborIdx][dirIdx] - normalNormCoordinate[dirIdx]);
366
367 }
368 }
369 }
370
371 // calculate or call all secondary variables
372 for (const auto& element : elements(this->gridGeometry().gridView()))
373 {
374 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
375
376 Dune::FieldMatrix<Scalar, GridView::dimension, GridView::dimension> stressTensor(0.0);
377 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
378 {
379 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
380 {
381 stressTensor[dimIdx][velIdx] = 0.5 * velocityGradients_[elementIdx][dimIdx][velIdx]
382 + 0.5 * velocityGradients_[elementIdx][velIdx][dimIdx];
383 }
384 }
385 stressTensorScalarProduct_[elementIdx] = 0.0;
386 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
387 {
388 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
389 {
390 stressTensorScalarProduct_[elementIdx] += stressTensor[dimIdx][velIdx] * stressTensor[dimIdx][velIdx];
391 }
392 }
393
394 Dune::FieldMatrix<Scalar, GridView::dimension, GridView::dimension> vorticityTensor(0.0);
395 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
396 {
397 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
398 {
399 vorticityTensor[dimIdx][velIdx] = 0.5 * velocityGradients_[elementIdx][dimIdx][velIdx]
400 - 0.5 * velocityGradients_[elementIdx][velIdx][dimIdx];
401 }
402 }
403 vorticityTensorScalarProduct_[elementIdx] = 0.0;
404 for (unsigned int dimIdx = 0; dimIdx < dim; ++dimIdx)
405 {
406 for (unsigned int velIdx = 0; velIdx < dim; ++velIdx)
407 {
408 vorticityTensorScalarProduct_[elementIdx] += vorticityTensor[dimIdx][velIdx] * vorticityTensor[dimIdx][velIdx];
409 }
410 }
411
412 auto fvGeometry = localView(this->gridGeometry());
413 fvGeometry.bindElement(element);
414 for (auto&& scv : scvs(fvGeometry))
415 {
416 const int dofIdx = scv.dofIndex();
417
418 // construct a privars object from the cell center solution vector
419 const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
420 PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
421 auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
422
423 VolumeVariables volVars;
424 volVars.update(elemSol, asImp_(), element, scv);
425 kinematicViscosity_[elementIdx] = volVars.viscosity() / volVars.density();
426 }
427 }
428 }
429
437 bool useWallFunction(const Element& element,
438 const SubControlVolumeFace& scvf,
439 const int& eqIdx) const
440 { return false; }
441
445 template<class ElementVolumeVariables, class ElementFaceVariables>
446 FacePrimaryVariables wallFunction(const Element& element,
447 const FVElementGeometry& fvGeometry,
448 const ElementVolumeVariables& elemVolVars,
449 const ElementFaceVariables& elemFaceVars,
450 const SubControlVolumeFace& scvf,
451 const SubControlVolumeFace& lateralBoundaryFace) const
452 { return FacePrimaryVariables(0.0); }
453
457 template<class ElementVolumeVariables, class ElementFaceVariables>
458 CellCenterPrimaryVariables wallFunction(const Element& element,
459 const FVElementGeometry& fvGeometry,
460 const ElementVolumeVariables& elemVolVars,
461 const ElementFaceVariables& elemFaceVars,
462 const SubControlVolumeFace& scvf) const
463 { return CellCenterPrimaryVariables(0.0); }
464
470 bool isOnWall(const SubControlVolumeFace& scvf) const
471 {
472 return asImp_().isOnWallAtPos(scvf.center());
473 }
474
480 bool isOnWallAtPos(const GlobalPosition &globalPos) const
481 {
482 // Throw an exception if no walls are implemented
483 DUNE_THROW(Dune::InvalidStateException,
484 "The problem does not provide an isOnWall() method.");
485 }
486
492 Scalar sandGrainRoughnessAtPos(const GlobalPosition &globalPos) const
493 {
494 return 0.0;
495 }
496
500 const Scalar karmanConstant() const
501 { return 0.41; }
502
504 const Scalar betaOmega() const
505 {
506 return 0.0708;
507 }
508
514 {
515 static const Scalar turbulentPrandtlNumber
516 = getParamFromGroup<Scalar>(this->paramGroup(), "RANS.TurbulentPrandtlNumber", 1.0);
518 }
519
525 {
526 static const Scalar turbulentSchmidtNumber
527 = getParamFromGroup<Scalar>(this->paramGroup(), "RANS.TurbulentSchmidtNumber", 1.0);
529 }
530
531public:
533 std::vector<unsigned int> wallElementIdx_;
534 std::vector<Scalar> wallDistance_;
535 std::vector<std::array<std::array<unsigned int, 2>, dim>> neighborIdx_;
536 std::vector<GlobalPosition> cellCenter_;
537 std::vector<DimVector> velocity_;
538 std::vector<DimVector> velocityMaximum_;
539 std::vector<DimVector> velocityMinimum_;
540 std::vector<DimMatrix> velocityGradients_;
541 std::vector<Scalar> stressTensorScalarProduct_;
542 std::vector<Scalar> vorticityTensorScalarProduct_;
543 std::vector<unsigned int> wallNormalAxis_;
544 std::vector<unsigned int> flowNormalAxis_;
545 std::vector<Scalar> kinematicViscosity_;
546 std::vector<Scalar> sandGrainRoughness_;
547
548private:
550 Implementation &asImp_()
551 { return *static_cast<Implementation *>(this); }
552
554 const Implementation &asImp_() const
555 { return *static_cast<const Implementation *>(this); }
556};
557
558} // end namespace Dumux
559
560#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.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
Dune::IteratorRange< typename MultiDomainGlue< DomainGridView, TargetGridView, DomainMapper, TargetMapper >::Intersections::const_iterator > intersections(const MultiDomainGlue< DomainGridView, TargetGridView, DomainMapper, TargetMapper > &glue)
Range generator to iterate with range-based for loops over all intersections as follows: for (const a...
Definition: glue.hh:62
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
Navier-Stokes problem base class.
Definition: dumux/freeflow/navierstokes/problem.hh:63
const Scalar bjsVelocity(const Element &element, const SubControlVolume &scv, const SubControlVolumeFace &faceOnPorousBoundary, const Scalar velocitySelf) const
helper function to evaluate the slip velocity on the boundary when the Beavers-Joseph-Saffman conditi...
Definition: dumux/freeflow/navierstokes/problem.hh:229
forward declare
Definition: dumux/freeflow/rans/problem.hh:41
Reynolds-Averaged Navier-Stokes problem base class.
Definition: dumux/freeflow/rans/problem.hh:57
std::vector< std::array< std::array< unsigned int, 2 >, dim > > neighborIdx_
Definition: dumux/freeflow/rans/problem.hh:535
void updateDynamicWallProperties(const SolutionVector &curSol)
Update the dynamic (solution dependent) relations to the walls.
Definition: dumux/freeflow/rans/problem.hh:225
bool isOnWallAtPos(const GlobalPosition &globalPos) const
Returns whether a given point is on a wall.
Definition: dumux/freeflow/rans/problem.hh:480
std::vector< DimVector > velocityMinimum_
Definition: dumux/freeflow/rans/problem.hh:539
void updateStaticWallProperties()
Update the static (solution independent) relations to the walls.
Definition: dumux/freeflow/rans/problem.hh:98
FacePrimaryVariables wallFunction(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const SubControlVolumeFace &lateralBoundaryFace) const
Returns an additional wall function momentum flux.
Definition: dumux/freeflow/rans/problem.hh:446
std::vector< DimMatrix > velocityGradients_
Definition: dumux/freeflow/rans/problem.hh:540
const Scalar betaOmega() const
Returns the constant.
Definition: dumux/freeflow/rans/problem.hh:504
std::vector< GlobalPosition > cellCenter_
Definition: dumux/freeflow/rans/problem.hh:536
std::vector< Scalar > vorticityTensorScalarProduct_
Definition: dumux/freeflow/rans/problem.hh:542
std::vector< DimVector > velocityMaximum_
Definition: dumux/freeflow/rans/problem.hh:538
Scalar turbulentSchmidtNumber() const
Return the turbulent Schmidt number which is used to convert the eddy viscosity to an eddy diffusivi...
Definition: dumux/freeflow/rans/problem.hh:524
CellCenterPrimaryVariables wallFunction(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Returns an additional wall function flux for cell-centered quantities.
Definition: dumux/freeflow/rans/problem.hh:458
std::vector< Scalar > sandGrainRoughness_
Definition: dumux/freeflow/rans/problem.hh:546
std::vector< Scalar > kinematicViscosity_
Definition: dumux/freeflow/rans/problem.hh:545
Scalar turbulentPrandtlNumber() const
Return the turbulent Prandtl number which is used to convert the eddy viscosity to an eddy thermal c...
Definition: dumux/freeflow/rans/problem.hh:513
std::vector< unsigned int > wallNormalAxis_
Definition: dumux/freeflow/rans/problem.hh:543
std::vector< Scalar > stressTensorScalarProduct_
Definition: dumux/freeflow/rans/problem.hh:541
std::vector< Scalar > wallDistance_
Definition: dumux/freeflow/rans/problem.hh:534
RANSProblemBase(std::shared_ptr< const GridGeometry > gridGeometry, const std::string &paramGroup="")
The constructor.
Definition: dumux/freeflow/rans/problem.hh:88
std::vector< DimVector > velocity_
Definition: dumux/freeflow/rans/problem.hh:537
bool calledUpdateStaticWallProperties
Definition: dumux/freeflow/rans/problem.hh:532
Scalar sandGrainRoughnessAtPos(const GlobalPosition &globalPos) const
Returns the sand-grain roughness at a given position.
Definition: dumux/freeflow/rans/problem.hh:492
bool useWallFunction(const Element &element, const SubControlVolumeFace &scvf, const int &eqIdx) const
Returns whether a wall function should be used at a given face.
Definition: dumux/freeflow/rans/problem.hh:437
const Scalar karmanConstant() const
Returns the Karman constant.
Definition: dumux/freeflow/rans/problem.hh:500
std::vector< unsigned int > wallElementIdx_
Definition: dumux/freeflow/rans/problem.hh:533
std::vector< unsigned int > flowNormalAxis_
Definition: dumux/freeflow/rans/problem.hh:544
bool isOnWall(const SubControlVolumeFace &scvf) const
Returns whether a given sub control volume face is on a wall.
Definition: dumux/freeflow/rans/problem.hh:470
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.