3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/rans/oneeq/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_ONEEQ_PROBLEM_HH
25#define DUMUX_ONEEQ_PROBLEM_HH
26
34
35#include "model.hh"
36
37namespace Dumux {
38
45template<class TypeTag>
46class RANSProblemImpl<TypeTag, TurbulenceModel::oneeq> : public RANSProblemBase<TypeTag>
47{
48
50 using Implementation = GetPropType<TypeTag, Properties::Problem>;
52 using Grid = typename GridView::Grid;
53
55 using DimVector = Dune::FieldVector<Scalar, Grid::dimension>;
56
58 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
59 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
60
64 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
66
67public:
69 RANSProblemImpl(std::shared_ptr<const GridGeometry> gridGeometry, const std::string& paramGroup = "")
70 : ParentType(gridGeometry, paramGroup)
71 { }
72
77 {
78 ParentType::updateStaticWallProperties();
79
80 // update size and initial values of the global vectors
81 storedDynamicEddyViscosity_.resize(this->gridGeometry().elementMapper().size(), 0.0);
82 storedViscosityTilde_.resize(this->gridGeometry().elementMapper().size(), 0.0);
83 storedViscosityTildeGradient_.resize(this->gridGeometry().elementMapper().size(), DimVector(0.0));
84 }
85
91 void updateDynamicWallProperties(const SolutionVector& curSol)
92 {
93 ParentType::updateDynamicWallProperties(curSol);
94
95 auto fvGeometry = localView(this->gridGeometry());
96 for (const auto& element : elements(this->gridGeometry().gridView()))
97 {
98 unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
99 fvGeometry.bindElement(element);
100
101 for (auto&& scv : scvs(fvGeometry))
102 {
103 const int dofIdx = scv.dofIndex();
104 const auto& cellCenterPriVars = curSol[GridGeometry::cellCenterIdx()][dofIdx];
105 PrimaryVariables priVars = makePriVarsFromCellCenterPriVars<PrimaryVariables>(cellCenterPriVars);
106 auto elemSol = elementSolution<typename GridGeometry::LocalView>(std::move(priVars));
107 // NOTE: first update the turbulence quantities
108 storedViscosityTilde_[elementIdx] = elemSol[0][Indices::viscosityTildeIdx];
109 // NOTE: then update the volVars
110 VolumeVariables volVars;
111 volVars.update(elemSol, asImp_(), element, scv);
112 storedDynamicEddyViscosity_[elementIdx] = volVars.calculateEddyViscosity();
113 }
114 }
115
116 // calculate cell-center-averaged gradient
117 for (const auto& element : elements(this->gridGeometry().gridView()))
118 {
119 const unsigned int elementIdx = this->gridGeometry().elementMapper().index(element);
120 fvGeometry.bindElement(element);
121
122 for (unsigned int axisIdx = 0; axisIdx < Grid::dimension; ++axisIdx)
123 {
124 const unsigned int neighborIndex0 = ParentType::neighborIndex(elementIdx, axisIdx, 0);
125 const unsigned int neighborIndex1 = ParentType::neighborIndex(elementIdx, axisIdx, 1);
126
127 // calculate cell-centered turbulentEddyViscosity (viscosityTilde) gradient
128 storedViscosityTildeGradient_[elementIdx][axisIdx]
129 = (storedViscosityTilde(neighborIndex1) - storedViscosityTilde(neighborIndex0))
130 / (ParentType::cellCenter(neighborIndex1)[axisIdx] - ParentType::cellCenter(neighborIndex0)[axisIdx]);
131 }
132
133 // Adjust for dirichlet boundary conditions
134 for (auto&& scvf : scvfs(fvGeometry))
135 {
136 const unsigned int normDim = scvf.directionIndex();
137 if (scvf.boundary() && asImp_().boundaryTypes(element, scvf).isDirichlet(Indices::viscosityTildeIdx))
138 {
139 // face Value
140 Scalar dirichletViscosityTilde = asImp_().dirichlet(element, scvf)[Indices::viscosityTildeIdx];
141
142 unsigned int neighborIndex = ParentType::neighborIndex(elementIdx, normDim, 0);
143 if (scvf.center()[normDim] < ParentType::cellCenter(elementIdx)[normDim])
144 neighborIndex = ParentType::neighborIndex(elementIdx, normDim, 1);
145
146 storedViscosityTildeGradient_[elementIdx][normDim]
147 = (storedViscosityTilde(neighborIndex) - dirichletViscosityTilde)
148 / (ParentType::cellCenter(neighborIndex)[normDim] - scvf.center()[normDim]);
149 }
150 }
151 }
152 }
153
155 {
156 static const bool useStoredEddyViscosity = getParamFromGroup<bool>(this->paramGroup(), "RANS.UseStoredEddyViscosity", false);
157 return useStoredEddyViscosity;
158 }
159
160 Scalar storedDynamicEddyViscosity(const int elementIdx) const
161 { return storedDynamicEddyViscosity_[elementIdx]; }
162
163 Scalar storedViscosityTilde(const int elementIdx) const
164 { return storedViscosityTilde_[elementIdx]; }
165
166 DimVector storedViscosityTildeGradient(const int elementIdx) const
167 { return storedViscosityTildeGradient_[elementIdx]; }
168
169private:
170 std::vector<Scalar> storedDynamicEddyViscosity_;
171 std::vector<Scalar> storedViscosityTilde_;
172 std::vector<DimVector> storedViscosityTildeGradient_;
173
175 Implementation &asImp_()
176 { return *static_cast<Implementation *>(this); }
177
179 const Implementation &asImp_() const
180 { return *static_cast<const Implementation *>(this); }
181};
182
183} // end namespace Dumux
184
185#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:38
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
Definition: adapt.hh:29
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:150
void updateStaticWallProperties()
Correct size of the static (solution independent) wall variables.
Definition: freeflow/rans/oneeq/problem.hh:76
bool useStoredEddyViscosity() const
Definition: freeflow/rans/oneeq/problem.hh:154
void updateDynamicWallProperties(const SolutionVector &curSol)
Update the dynamic (solution dependent) relations to the walls.
Definition: freeflow/rans/oneeq/problem.hh:91
Scalar storedViscosityTilde(const int elementIdx) const
Definition: freeflow/rans/oneeq/problem.hh:163
RANSProblemImpl(std::shared_ptr< const GridGeometry > gridGeometry, const std::string &paramGroup="")
The constructor sets the gravity, if desired by the user.
Definition: freeflow/rans/oneeq/problem.hh:69
DimVector storedViscosityTildeGradient(const int elementIdx) const
Definition: freeflow/rans/oneeq/problem.hh:166
Scalar storedDynamicEddyViscosity(const int elementIdx) const
Definition: freeflow/rans/oneeq/problem.hh:160
forward declare
Definition: freeflow/rans/problem.hh:44
Reynolds-Averaged Navier-Stokes problem base class.
Definition: freeflow/rans/problem.hh:60
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.