12#ifndef DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
13#define DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
15#include <dune/common/hybridutilities.hh>
29template<
int radialAxis>
30static constexpr bool isRotationalExtrusion<RotationalExtrusion<radialAxis>> =
true;
39template<
class TypeTag,
class DiscretizationMethod>
42template<
class TypeTag>
51 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
52 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
53 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
55 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
56 using ElementFluxVariablesCache =
typename GridFluxVariablesCache::LocalView;
58 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
59 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
65 using FVElementGeometry =
typename GridGeometry::LocalView;
66 using GridView =
typename GridGeometry::GridView;
67 using Element =
typename GridView::template Codim<0>::Entity;
68 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
69 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
77 static constexpr bool normalizePressure = getPropValue<TypeTag, Properties::NormalizePressure>();
79 using CellCenterResidual = CellCenterPrimaryVariables;
80 using FaceResidual = FacePrimaryVariables;
88 static constexpr auto cellCenterOffset = ModelTraits::numEq() - CellCenterPrimaryVariables::dimension;
89 static_assert(cellCenterOffset == ModelTraits::dim(),
"cellCenterOffset must equal dim for staggered NavierStokes");
92 using ParentType::ParentType;
96 const Element &element,
97 const FVElementGeometry& fvGeometry,
98 const ElementVolumeVariables& elemVolVars,
99 const ElementFaceVariables& elemFaceVars,
100 const SubControlVolumeFace &scvf,
101 const ElementFluxVariablesCache& elemFluxVarsCache)
const
103 FluxVariables fluxVars;
104 CellCenterPrimaryVariables flux = fluxVars.computeMassFlux(problem, element, fvGeometry, elemVolVars,
105 elemFaceVars, scvf, elemFluxVarsCache[scvf]);
107 EnergyLocalResidual::heatFlux(flux, problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf);
114 const Element &element,
115 const FVElementGeometry& fvGeometry,
116 const ElementVolumeVariables& elemVolVars,
117 const ElementFaceVariables& elemFaceVars,
118 const SubControlVolume &scv)
const
120 CellCenterPrimaryVariables result(0.0);
123 const auto sourceValues = problem.source(element, fvGeometry, elemVolVars, elemFaceVars, scv);
126 for (
int i = 0; i < result.size(); ++i)
127 result[i] = sourceValues[i + cellCenterOffset];
135 const SubControlVolume& scv,
136 const VolumeVariables& volVars)
const
138 CellCenterPrimaryVariables storage;
139 storage[Indices::conti0EqIdx - ModelTraits::dim()] = volVars.density();
141 EnergyLocalResidual::fluidPhaseStorage(storage, volVars);
148 const SubControlVolumeFace& scvf,
149 const VolumeVariables& volVars,
150 const ElementFaceVariables& elemFaceVars)
const
152 FacePrimaryVariables storage(0.0);
153 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
154 storage[0] = volVars.density() * velocity;
160 const Element& element,
161 const FVElementGeometry& fvGeometry,
162 const SubControlVolumeFace& scvf,
163 const ElementVolumeVariables& elemVolVars,
164 const ElementFaceVariables& elemFaceVars)
const
166 FacePrimaryVariables source(0.0);
167 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
168 source += problem.gravity()[scvf.directionIndex()] * insideVolVars.density();
169 source += problem.source(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[Indices::velocity(scvf.directionIndex())];
174 if constexpr (ModelTraits::dim() == 2 && Impl::isRotationalExtrusion<Extrusion>)
176 if (scvf.directionIndex() == Extrusion::radialAxis)
179 const auto r = scvf.center()[scvf.directionIndex()];
180 source -= -2.0*insideVolVars.effectiveViscosity() * elemFaceVars[scvf].velocitySelf() / (r*r);
186 normalizePressure ? insideVolVars.pressure() - problem.initial(scvf)[Indices::pressureIdx]
187 : insideVolVars.pressure();
198 const Element& element,
199 const SubControlVolumeFace& scvf,
200 const FVElementGeometry& fvGeometry,
201 const ElementVolumeVariables& elemVolVars,
202 const ElementFaceVariables& elemFaceVars,
203 const ElementFluxVariablesCache& elemFluxVarsCache)
const
205 FluxVariables fluxVars;
206 return fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
213 const Element& element,
214 const FVElementGeometry& fvGeometry,
215 const SubControlVolumeFace& scvf,
216 const ElementVolumeVariables& elemVolVars,
217 const ElementFaceVariables& elemFaceVars,
218 const ElementBoundaryTypes& elemBcTypes,
219 const ElementFluxVariablesCache& elemFluxVarsCache)
const
221 CellCenterResidual result(0.0);
225 const auto bcTypes = problem.boundaryTypes(element, scvf);
228 if (bcTypes.isSymmetry())
232 result = computeFluxForCellCenter(problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf, elemFluxVarsCache);
235 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
236 if (bcTypes.hasNeumann())
238 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
239 const auto neumannFluxes = problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
241 for (
int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
243 if (bcTypes.isNeumann(eqIdx + cellCenterOffset))
244 result[eqIdx] = neumannFluxes[eqIdx + cellCenterOffset] * extrusionFactor * Extrusion::area(fvGeometry, scvf);
249 incorporateWallFunction_(result, problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars);
258 const Problem& problem,
259 const Element& element,
260 const FVElementGeometry& fvGeometry,
261 const SubControlVolumeFace& scvf,
262 const ElementVolumeVariables& elemVolVars,
263 const ElementFaceVariables& elemFaceVars,
264 const ElementBoundaryTypes& elemBcTypes,
265 const ElementFluxVariablesCache& elemFluxVarsCache)
const
270 const auto bcTypes = problem.boundaryTypes(element, scvf);
272 if(bcTypes.isDirichlet(Indices::velocity(scvf.directionIndex())))
275 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
276 const Scalar dirichletValue = problem.dirichlet(element, scvf)[Indices::velocity(scvf.directionIndex())];
277 residual = velocity - dirichletValue;
279 else if(bcTypes.isSymmetry())
283 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
284 const Scalar fixedValue = 0.0;
285 residual = velocity - fixedValue;
294 const Element& element,
295 const FVElementGeometry& fvGeometry,
296 const SubControlVolumeFace& scvf,
297 const ElementVolumeVariables& elemVolVars,
298 const ElementFaceVariables& elemFaceVars,
299 const ElementBoundaryTypes& elemBcTypes,
300 const ElementFluxVariablesCache& elemFluxVarsCache)
const
302 FaceResidual result(0.0);
306 FluxVariables fluxVars;
309 const auto bcTypes = problem.boundaryTypes(element, scvf);
310 if (bcTypes.isNeumann(Indices::velocity(scvf.directionIndex())))
314 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
315 result = problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[Indices::velocity(scvf.directionIndex())]
316 * extrusionFactor * Extrusion::area(fvGeometry, scvf);
319 result += fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
321 else if(bcTypes.isDirichlet(Indices::pressureIdx))
325 result = fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
328 result += fluxVars.inflowOutflowBoundaryFlux(problem, scvf, fvGeometry, elemVolVars, elemFaceVars);
337 template<
class ...Args,
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<!turbulenceModel, int> = 0>
338 void incorporateWallFunction_(Args&&... args)
const
342 template<
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<turbulenceModel,
int> = 0>
343 void incorporateWallFunction_(CellCenterResidual& boundaryFlux,
344 const Problem& problem,
345 const Element& element,
346 const FVElementGeometry& fvGeometry,
347 const SubControlVolumeFace& scvf,
348 const ElementVolumeVariables& elemVolVars,
349 const ElementFaceVariables& elemFaceVars)
const
351 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
352 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
355 for(
int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
358 if(problem.useWallFunction(element, scvf, eqIdx + cellCenterOffset))
360 boundaryFlux[eqIdx] = problem.wallFunction(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[eqIdx]
361 * extrusionFactor * Extrusion::area(fvGeometry, scvf);
367 Implementation &asImp_()
368 {
return *
static_cast<Implementation *
>(
this); }
371 const Implementation &asImp_()
const
372 {
return *
static_cast<const Implementation *
>(
this); }
Definition: freeflow/nonisothermal/localresidual.hh:24
FacePrimaryVariables computeFluxForFace(const Problem &problem, const Element &element, const SubControlVolumeFace &scvf, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate the momentum flux for the face control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:197
void evalDirichletBoundariesForFace(FaceResidual &residual, const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate Dirichlet (fixed value) boundary conditions for a face dof.
Definition: freeflow/navierstokes/staggered/localresidual.hh:257
FacePrimaryVariables computeStorageForFace(const Problem &problem, const SubControlVolumeFace &scvf, const VolumeVariables &volVars, const ElementFaceVariables &elemFaceVars) const
Evaluate the storage term for the face control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:147
FacePrimaryVariables computeSourceForFace(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars) const
Evaluate the source term for the face control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:159
FaceResidual computeBoundaryFluxForFace(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate boundary fluxes for a face dof.
Definition: freeflow/navierstokes/staggered/localresidual.hh:293
CellCenterResidual computeBoundaryFluxForCellCenter(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const ElementBoundaryTypes &elemBcTypes, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate boundary conditions for a cell center dof.
Definition: freeflow/navierstokes/staggered/localresidual.hh:212
CellCenterPrimaryVariables computeFluxForCellCenter(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf, const ElementFluxVariablesCache &elemFluxVarsCache) const
Evaluate fluxes entering or leaving the cell center control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:95
CellCenterPrimaryVariables computeStorageForCellCenter(const Problem &problem, const SubControlVolume &scv, const VolumeVariables &volVars) const
Evaluate the storage term for the cell center control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:134
CellCenterPrimaryVariables computeSourceForCellCenter(const Problem &problem, const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolume &scv) const
Evaluate the source term for the cell center control volume.
Definition: freeflow/navierstokes/staggered/localresidual.hh:113
Element-wise calculation of the Navier-Stokes residual for models using the staggered discretization.
Definition: freeflow/navierstokes/localresidual.hh:23
Calculates the element-wise residual for the staggered FV scheme.
Definition: staggeredlocalresidual.hh:29
Defines all properties used in Dumux.
Helper classes to compute the integration elements.
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
The available discretization methods in Dumux.
std::string pressure(int phaseIdx) noexcept
I/O name of pressure for multiphase systems.
Definition: name.hh:22
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:166
constexpr bool isRotationalExtrusion
Convenience trait to check whether the extrusion is rotational.
Definition: extrusion.hh:172
Calculates the element-wise residual for the staggered FV scheme.