12#ifndef DUMUX_SST_STAGGERED_FLUXVARIABLES_HH
13#define DUMUX_SST_STAGGERED_FLUXVARIABLES_HH
16#include <dune/common/exceptions.hh>
33template<
class TypeTag,
class BaseFluxVariables,
class DiscretizationMethod>
34class SSTFluxVariablesImpl;
36template<
class TypeTag,
class BaseFluxVariables>
38:
public BaseFluxVariables
40 using ParentType = BaseFluxVariables;
44 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
45 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
46 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
48 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
49 using FluxVariablesCache =
typename GridFluxVariablesCache::FluxVariablesCache;
51 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
52 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
53 using FaceVariables =
typename GridFaceVariables::FaceVariables;
58 using FVElementGeometry =
typename GridGeometry::LocalView;
59 using SubControlVolumeFace =
typename GridGeometry::SubControlVolumeFace;
61 using GridView =
typename GridGeometry::GridView;
63 using Element =
typename GridView::template Codim<0>::Entity;
68 static constexpr int turbulentKineticEnergyEqIdx = Indices::turbulentKineticEnergyEqIdx - ModelTraits::dim();
69 static constexpr int dissipationEqIdx = Indices::dissipationEqIdx - ModelTraits::dim();
77 const Element &element,
78 const FVElementGeometry& fvGeometry,
79 const ElementVolumeVariables& elemVolVars,
80 const ElementFaceVariables& elemFaceVars,
81 const SubControlVolumeFace &scvf,
82 const FluxVariablesCache& fluxVarsCache)
84 CellCenterPrimaryVariables flux = ParentType::computeMassFlux(problem, element, fvGeometry,
85 elemVolVars, elemFaceVars, scvf, fluxVarsCache);
88 auto upwindTermK = [](
const auto& volVars)
89 {
return volVars.turbulentKineticEnergy() * volVars.density(); };
90 auto upwindTermOmega = [](
const auto& volVars)
91 {
return volVars.dissipation() * volVars.density(); };
93 flux[turbulentKineticEnergyEqIdx]
94 = ParentType::advectiveFluxForCellCenter(problem, fvGeometry, elemVolVars, elemFaceVars, scvf, upwindTermK);
95 flux[dissipationEqIdx]
96 = ParentType::advectiveFluxForCellCenter(problem, fvGeometry, elemVolVars, elemFaceVars, scvf, upwindTermOmega);
99 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
100 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
101 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
102 const auto& outsideVolVars = elemVolVars[scvf.outsideScvIdx()];
103 Scalar insideCoeff_k = 0.0, insideCoeff_w = 0.0, outsideCoeff_k = 0.0, outsideCoeff_w = 0.0;
107 insideCoeff_k = insideVolVars.viscosity()
108 + ( insideVolVars.sigmaKBSL() * insideVolVars.dynamicEddyViscosity() );
109 outsideCoeff_k = outsideVolVars.viscosity()
110 + ( outsideVolVars.sigmaKBSL() * outsideVolVars.dynamicEddyViscosity() );
111 insideCoeff_w = insideVolVars.viscosity()
112 + ( insideVolVars.sigmaOmegaBSL() * insideVolVars.dynamicEddyViscosity() );
113 outsideCoeff_w = outsideVolVars.viscosity()
114 + ( outsideVolVars.sigmaOmegaBSL() * outsideVolVars.dynamicEddyViscosity() );
118 insideCoeff_k = insideVolVars.viscosity()
119 + ( insideVolVars.sigmaKSST() * insideVolVars.dynamicEddyViscosity() );
120 outsideCoeff_k = outsideVolVars.viscosity()
121 + ( outsideVolVars.sigmaKSST() * outsideVolVars.dynamicEddyViscosity() );
122 insideCoeff_w = insideVolVars.viscosity()
123 + ( insideVolVars.sigmaOmegaSST() * insideVolVars.dynamicEddyViscosity() );
124 outsideCoeff_w = outsideVolVars.viscosity()
125 + ( outsideVolVars.sigmaOmegaSST() * outsideVolVars.dynamicEddyViscosity() );
128 DUNE_THROW(Dune::NotImplemented,
"\nThis SST Model is not implemented.\n");
132 insideCoeff_k *= insideVolVars.extrusionFactor();
133 outsideCoeff_k *= outsideVolVars.extrusionFactor();
134 insideCoeff_w *= insideVolVars.extrusionFactor();
135 outsideCoeff_w *= outsideVolVars.extrusionFactor();
138 Scalar coeff_k = 0.0;
139 Scalar coeff_w = 0.0;
142 distance = (insideScv.dofPosition() - scvf.ipGlobal()).two_norm();
143 coeff_k = insideCoeff_k;
144 coeff_w = insideCoeff_w;
150 (outsideScv.dofPosition() - scvf.ipGlobal()).two_norm(),
151 (insideScv.dofPosition() - scvf.ipGlobal()).two_norm());
153 (outsideScv.dofPosition() - scvf.ipGlobal()).two_norm(),
154 (insideScv.dofPosition() - scvf.ipGlobal()).two_norm());
155 distance = (outsideScv.dofPosition() - insideScv.dofPosition()).two_norm();
158 const auto bcTypes = problem.boundaryTypes(element, scvf);
159 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::turbulentKineticEnergyEqIdx)
160 || bcTypes.isSymmetry())))
162 flux[turbulentKineticEnergyEqIdx]
164 * (insideVolVars.turbulentKineticEnergy() - outsideVolVars.turbulentKineticEnergy())
165 * Extrusion::area(fvGeometry, scvf);
167 if (!(scvf.boundary() && (bcTypes.isOutflow(Indices::dissipationEqIdx)
168 || bcTypes.isSymmetry())))
170 flux[dissipationEqIdx]
172 * (insideVolVars.dissipation() - outsideVolVars.dissipation())
173 * Extrusion::area(fvGeometry, scvf);
182 const Element& element,
183 const SubControlVolumeFace& scvf,
184 const FVElementGeometry& fvGeometry,
185 const ElementVolumeVariables& elemVolVars,
186 const ElementFaceVariables& elemFaceVars,
187 const GridFluxVariablesCache& gridFluxVarsCache)
189 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
191 return ParentType::computeFrontalMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
192 + ParentType::computeLateralMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, gridFluxVarsCache)
193 + 2.0 / ModelTraits::dim() * insideVolVars.density() * insideVolVars.turbulentKineticEnergy()
194 * Extrusion::area(fvGeometry, scvf) * scvf.directionSign() * insideVolVars.extrusionFactor();
CellCenterPrimaryVariables computeMassFlux(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const FluxVariablesCache &fluxVarsCache)
Computes the flux for the cell center residual.
Definition: freeflow/rans/twoeq/sst/staggered/fluxvariables.hh:76
FacePrimaryVariables computeMomentumFlux(const Problem &problem, const Element &element, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const GridFluxVariablesCache &gridFluxVarsCache)
Returns the momentum flux over all staggered faces.
Definition: freeflow/rans/twoeq/sst/staggered/fluxvariables.hh:181
The flux variables class for the SST model using the staggered grid discretization.
Definition: freeflow/rans/twoeq/sst/fluxvariables.hh:22
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
Base class for the flux variables living on a sub control volume face.
constexpr Scalar arithmeticMean(Scalar x, Scalar y, Scalar wx=1.0, Scalar wy=1.0) noexcept
Calculate the (weighted) arithmetic mean of two scalar values.
Definition: math.hh:38
static ctype distance(const Dune::FieldVector< ctype, dimWorld > &a, const Dune::FieldVector< ctype, dimWorld > &b)
Compute the shortest distance between two points.
Definition: distance.hh:282
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The available discretization methods in Dumux.
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:166
The available free flow turbulence models in Dumux.