version 3.8.0
Loading...
Searching...
No Matches
freeflow/rans/twoeq/lowrekepsilon/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_LOWREKEPSILON_PROBLEM_HH
13#define DUMUX_LOWREKEPSILON_PROBLEM_HH
14
22
23#include "model.hh"
24
25namespace Dumux {
26
33template<class TypeTag>
35{
37 using Implementation = GetPropType<TypeTag, Properties::Problem>;
39
41 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
42
45 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
47
48public:
50 RANSProblemImpl(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
51 : ParentType(gridGeometry, paramGroup)
52 { }
53
58 {
59 ParentType::updateStaticWallProperties();
60
61 // update size and initial values of the global vectors
62 storedDissipationTilde_.resize(this->gridGeometry().elementMapper().size(), 0.0);
63 storedDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
64 storedTurbulentKineticEnergy_.resize(this->gridGeometry().elementMapper().size(), 0.0);
65 }
66
72 template<class SolutionVector>
73 void updateDynamicWallProperties(const SolutionVector& curSol)
74 {
75 ParentType::updateDynamicWallProperties(curSol);
76
77 auto fvGeometry = localView(this->gridGeometry());
78 for (const auto& element : elements(this->gridGeometry().gridView()))
79 {
80 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
81 fvGeometry.bindElement(element);
82
83 for (auto&& scv : scvs(fvGeometry))
84 {
85 const int dofIdx = scv.dofIndex();
86 const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
87 PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
88 auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
89 // NOTE: first update the turbulence quantities
90 storedDissipationTilde_[elementIdx] = elemSol[0][Indices::dissipationEqIdx];
91 storedTurbulentKineticEnergy_[elementIdx] = elemSol[0][Indices::turbulentKineticEnergyEqIdx];
92 // NOTE: then update the volVars
93 VolumeVariables volVars;
94 volVars.update(elemSol, asImp_(), element, scv);
95 storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscosity();
96 }
97 }
98 }
99
101 {
102 static const bool useStoredEddyViscosity = getParamFromGroup<bool>(this->paramGroup(), "RANS.UseStoredEddyViscosity", true);
103 return useStoredEddyViscosity;
104 }
105
106 Scalar storedDissipationTilde(const int elementIdx) const
107 { return storedDissipationTilde_[elementIdx]; }
108
109 Scalar storedDynamicEddyViscosity(const int elementIdx) const
110 { return storedDynamicEddyViscosity_[elementIdx]; }
111
112 Scalar storedTurbulentKineticEnergy(const int elementIdx) const
113 { return storedTurbulentKineticEnergy_[elementIdx]; }
114
115private:
116 std::vector<Scalar> storedDissipationTilde_;
117 std::vector<Scalar> storedDynamicEddyViscosity_;
118 std::vector<Scalar> storedTurbulentKineticEnergy_;
119
121 Implementation &asImp_()
122 { return *static_cast<Implementation *>(this); }
123
125 const Implementation &asImp_() const
126 { return *static_cast<const Implementation *>(this); }
127};
128
129}
130
131#endif
Reynolds-Averaged Navier-Stokes problem base class.
Definition freeflow/rans/problem.hh:47
Scalar storedDissipationTilde(const int elementIdx) const
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:106
void updateDynamicWallProperties(const SolutionVector &curSol)
Update the dynamic (solution dependent) relations to the walls.
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:73
void updateStaticWallProperties()
Correct size of the static (solution independent) wall variables.
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:57
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/lowrekepsilon/problem.hh:50
Scalar storedTurbulentKineticEnergy(const int elementIdx) const
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:112
Scalar storedDynamicEddyViscosity(const int elementIdx) const
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:109
bool useStoredEddyViscosity() const
Definition freeflow/rans/twoeq/lowrekepsilon/problem.hh:100
forward declare
Definition freeflow/rans/problem.hh:31
Defines all properties used in Dumux.
the turbulence-model-specfic RANS problem
A single-phase, isothermal low-Reynolds 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
TurbulenceModel
The available free flow turbulence models in Dumux.
Definition turbulencemodel.hh:26
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition propertysystem.hh:296
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.