12#ifndef DUMUX_SST_STAGGERED_KOMEGA_LOCAL_RESIDUAL_HH
13#define DUMUX_SST_STAGGERED_KOMEGA_LOCAL_RESIDUAL_HH
15#include <dune/common/hybridutilities.hh>
16#include <dune/common/exceptions.hh>
30template<
class TypeTag,
class BaseLocalRes
idual,
class DiscretizationMethod>
33template<
class TypeTag,
class BaseLocalRes
idual>
34class SSTResidualImpl<TypeTag, BaseLocalResidual, DiscretizationMethods::Staggered>
35:
public BaseLocalResidual
37 using ParentType = BaseLocalResidual;
41 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
42 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
43 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
45 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
46 using ElementFluxVariablesCache =
typename GridFluxVariablesCache::LocalView;
49 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
50 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
55 using Element =
typename GridView::template Codim<0>::Entity;
57 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
59 using CellCenterResidual = CellCenterPrimaryVariables;
63 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
64 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
67 using ParentType::ParentType;
71 const SubControlVolume& scv,
72 const VolumeVariables& volVars)
const
74 CellCenterPrimaryVariables storage = ParentType::computeStorageForCellCenter(problem, scv, volVars);
76 storage[turbulentKineticEnergyEqIdx] = volVars.turbulentKineticEnergy()*volVars.density();
77 storage[dissipationEqIdx] = volVars.dissipation()*volVars.density();
83 const Element &element,
84 const FVElementGeometry& fvGeometry,
85 const ElementVolumeVariables& elemVolVars,
86 const ElementFaceVariables& elemFaceVars,
87 const SubControlVolume &scv)
const
89 CellCenterPrimaryVariables source = ParentType::computeSourceForCellCenter(problem, element, fvGeometry,
90 elemVolVars, elemFaceVars, scv);
93 const auto& volVars = elemVolVars[scv];
98 Scalar productionTerm = 2.0 * volVars.dynamicEddyViscosity() * volVars.stressTensorScalarProduct();
99 source[turbulentKineticEnergyEqIdx] += productionTerm;
100 source[dissipationEqIdx] += volVars.gammaBSL() / volVars.kinematicEddyViscosity() * productionTerm;
103 source[turbulentKineticEnergyEqIdx] -= volVars.betaStarBSL() * volVars.density() * volVars.turbulentKineticEnergy() * volVars.dissipation();
104 source[dissipationEqIdx] -= volVars.betaBSL() * volVars.density() * volVars.dissipation() * volVars.dissipation();
107 Scalar gradientProduct = 0.0;
108 for (
unsigned int i = 0; i < ModelTraits::dim(); ++i)
109 gradientProduct += volVars.storedTurbulentKineticEnergyGradient()[i]
110 * volVars.storedDissipationGradient()[i];
112 source[dissipationEqIdx] += 2.0 * volVars.density() * (1.0-volVars.F1()) * volVars.sigmaOmega2() / volVars.dissipation() * gradientProduct;
117 Scalar productionTerm = 2.0 * volVars.dynamicEddyViscosity() * volVars.stressTensorScalarProduct();
118 source[turbulentKineticEnergyEqIdx] += productionTerm;
119 source[dissipationEqIdx] += volVars.gammaSST() / volVars.kinematicEddyViscosity() * productionTerm;
122 source[turbulentKineticEnergyEqIdx] -= volVars.betaStarSST() * volVars.density() * volVars.turbulentKineticEnergy() * volVars.dissipation();
123 source[dissipationEqIdx] -= volVars.betaSST() * volVars.density() * volVars.dissipation() * volVars.dissipation();
126 Scalar gradientProduct = 0.0;
127 for (
unsigned int i = 0; i < ModelTraits::dim(); ++i)
128 gradientProduct += volVars.storedTurbulentKineticEnergyGradient()[i]
129 * volVars.storedDissipationGradient()[i];
130 source[dissipationEqIdx] += 2.0 * volVars.density() * (1.0-volVars.F1()) * volVars.sigmaOmega2() / volVars.dissipation() * gradientProduct;
133 DUNE_THROW(Dune::NotImplemented,
"\nThis SST Model is not implemented.\n");
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/sst/staggered/localresidual.hh:82
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/sst/staggered/localresidual.hh:70
Element-wise calculation of the residual for SST models using the staggered discretization.
Definition: freeflow/rans/twoeq/sst/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.
The available free flow turbulence models in Dumux.