3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/rans/twoeq/lowrekepsilon/staggered/localresidual.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_STAGGERED_LOWREKEPSILON_LOCAL_RESIDUAL_HH
25#define DUMUX_STAGGERED_LOWREKEPSILON_LOCAL_RESIDUAL_HH
26
27#include <dune/common/hybridutilities.hh>
31
32namespace Dumux {
33
34// forward declaration
35template<class TypeTag, class BaseLocalResidual, DiscretizationMethod discMethod>
36class LowReKEpsilonResidualImpl;
37
42template<class TypeTag, class BaseLocalResidual>
43class LowReKEpsilonResidualImpl<TypeTag, BaseLocalResidual, DiscretizationMethod::staggered>
44: public BaseLocalResidual
45{
46 using ParentType = BaseLocalResidual;
47
49
50 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
51 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
52 using VolumeVariables = typename GridVolumeVariables::VolumeVariables;
53
54 using GridFluxVariablesCache = typename GridVariables::GridFluxVariablesCache;
55 using ElementFluxVariablesCache = typename GridFluxVariablesCache::LocalView;
57
58 using GridFaceVariables = typename GridVariables::GridFaceVariables;
59 using ElementFaceVariables = typename GridFaceVariables::LocalView;
60
64 using Element = typename GridView::template Codim<0>::Entity;
65 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
66 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
67 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
71
72 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
73 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
74
75public:
76 using ParentType::ParentType;
77
79 CellCenterPrimaryVariables computeStorageForCellCenter(const Problem& problem,
80 const SubControlVolume& scv,
81 const VolumeVariables& volVars) const
82 {
83 CellCenterPrimaryVariables storage = ParentType::computeStorageForCellCenter(problem, scv, volVars);
84
85 storage[turbulentKineticEnergyEqIdx] = volVars.turbulentKineticEnergy();
86 storage[dissipationEqIdx] = volVars.dissipationTilde();
87
88 return storage;
89 }
90
91 CellCenterPrimaryVariables computeSourceForCellCenter(const Problem& problem,
92 const Element &element,
93 const FVElementGeometry& fvGeometry,
94 const ElementVolumeVariables& elemVolVars,
95 const ElementFaceVariables& elemFaceVars,
96 const SubControlVolume &scv) const
97 {
98 CellCenterPrimaryVariables source = ParentType::computeSourceForCellCenter(problem, element, fvGeometry,
99 elemVolVars, elemFaceVars, scv);
100
101 const auto& volVars = elemVolVars[scv];
102
103 // production
104 source[turbulentKineticEnergyEqIdx] += 2.0 * volVars.kinematicEddyViscosity()
105 * volVars.stressTensorScalarProduct();
106 source[dissipationEqIdx] += volVars.cOneEpsilon() * volVars.fOne()
107 * volVars.dissipationTilde() / volVars.turbulentKineticEnergy()
108 * 2.0 * volVars.kinematicEddyViscosity()
109 * volVars.stressTensorScalarProduct();
110
111 // destruction
112 source[turbulentKineticEnergyEqIdx] -= volVars.dissipationTilde();
113 source[dissipationEqIdx] -= volVars.cTwoEpsilon() * volVars.fTwo()
114 * volVars.dissipationTilde() * volVars.dissipationTilde()
115 / volVars.turbulentKineticEnergy();
116
117 // dampening functions
118 source[turbulentKineticEnergyEqIdx] -= volVars.dValue();
119 source[dissipationEqIdx] += volVars.eValue();
120
121 return source;
122 }
123};
124} // end namespace Dumux
125
126#endif
The available discretization methods in Dumux.
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
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
Definition: freeflow/rans/twoeq/lowrekepsilon/localresidual.hh:36
CellCenterPrimaryVariables computeSourceForCellCenter(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolume &scv) const
Definition: freeflow/rans/twoeq/lowrekepsilon/staggered/localresidual.hh:91
CellCenterPrimaryVariables computeStorageForCellCenter(const Problem &problem, const SubControlVolume &scv, const VolumeVariables &volVars) const
Evaluate fluxes entering or leaving the cell center control volume.
Definition: freeflow/rans/twoeq/lowrekepsilon/staggered/localresidual.hh:79
Declares all properties used in Dumux.