12#ifndef DUMUX_STAGGERED_LOWREKEPSILON_LOCAL_RESIDUAL_HH
13#define DUMUX_STAGGERED_LOWREKEPSILON_LOCAL_RESIDUAL_HH
15#include <dune/common/hybridutilities.hh>
28template<
class TypeTag,
class BaseLocalRes
idual,
class DiscretizationMethod>
29class LowReKEpsilonResidualImpl;
31template<
class TypeTag,
class BaseLocalRes
idual>
33:
public BaseLocalResidual
35 using ParentType = BaseLocalResidual;
39 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
40 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
41 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
43 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
44 using ElementFluxVariablesCache =
typename GridFluxVariablesCache::LocalView;
47 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
48 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
53 using Element =
typename GridView::template Codim<0>::Entity;
55 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
60 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
61 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
64 using ParentType::ParentType;
68 const SubControlVolume& scv,
69 const VolumeVariables& volVars)
const
71 CellCenterPrimaryVariables storage = ParentType::computeStorageForCellCenter(problem, scv, volVars);
73 storage[turbulentKineticEnergyEqIdx] = volVars.turbulentKineticEnergy() * volVars.density();
74 storage[dissipationEqIdx] = volVars.dissipationTilde() * volVars.density();
80 const Element &element,
81 const FVElementGeometry& fvGeometry,
82 const ElementVolumeVariables& elemVolVars,
83 const ElementFaceVariables& elemFaceVars,
84 const SubControlVolume &scv)
const
86 CellCenterPrimaryVariables source = ParentType::computeSourceForCellCenter(problem, element, fvGeometry,
87 elemVolVars, elemFaceVars, scv);
89 const auto& volVars = elemVolVars[scv];
92 source[turbulentKineticEnergyEqIdx] += 2.0 * volVars.dynamicEddyViscosity()
93 * volVars.stressTensorScalarProduct();
94 source[dissipationEqIdx] += volVars.cOneEpsilon() * volVars.fOne()
95 * volVars.dissipationTilde() / volVars.turbulentKineticEnergy()
96 * 2.0 * volVars.dynamicEddyViscosity()
97 * volVars.stressTensorScalarProduct();
100 source[turbulentKineticEnergyEqIdx] -= volVars.dissipationTilde() * volVars.density();
101 source[dissipationEqIdx] -= volVars.cTwoEpsilon() * volVars.fTwo() * volVars.density()
102 * volVars.dissipationTilde() * volVars.dissipationTilde()
103 / volVars.turbulentKineticEnergy();
106 source[turbulentKineticEnergyEqIdx] -= volVars.dValue() * volVars.density();
107 source[dissipationEqIdx] += volVars.eValue() * volVars.density();
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:79
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:67
Element-wise calculation of the residual for low-Reynolds k-epsilon models using the staggered discre...
Definition: freeflow/rans/twoeq/lowrekepsilon/localresidual.hh:24
Defines all properties used in Dumux.
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The available discretization methods in Dumux.