24#ifndef DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
25#define DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
30#include <dune/common/hybridutilities.hh>
36template<
class TypeTag, DiscretizationMethod discMethod>
43template<
class TypeTag>
52 using GridVolumeVariables =
typename GridVariables::GridVolumeVariables;
53 using ElementVolumeVariables =
typename GridVolumeVariables::LocalView;
54 using VolumeVariables =
typename GridVolumeVariables::VolumeVariables;
56 using GridFluxVariablesCache =
typename GridVariables::GridFluxVariablesCache;
57 using ElementFluxVariablesCache =
typename GridFluxVariablesCache::LocalView;
59 using GridFaceVariables =
typename GridVariables::GridFaceVariables;
60 using ElementFaceVariables =
typename GridFaceVariables::LocalView;
66 using FVElementGeometry =
typename GridGeometry::LocalView;
67 using GridView =
typename GridGeometry::GridView;
68 using Element =
typename GridView::template Codim<0>::Entity;
69 using SubControlVolume =
typename FVElementGeometry::SubControlVolume;
70 using SubControlVolumeFace =
typename FVElementGeometry::SubControlVolumeFace;
77 using CellCenterResidual = CellCenterPrimaryVariables;
78 using FaceResidual = FacePrimaryVariables;
86 static constexpr auto cellCenterOffset = ModelTraits::numEq() - CellCenterPrimaryVariables::dimension;
87 static_assert(cellCenterOffset == ModelTraits::dim(),
"cellCenterOffset must equal dim for staggered NavierStokes");
90 using ParentType::ParentType;
94 const Element &element,
95 const FVElementGeometry& fvGeometry,
96 const ElementVolumeVariables& elemVolVars,
97 const ElementFaceVariables& elemFaceVars,
98 const SubControlVolumeFace &scvf,
99 const ElementFluxVariablesCache& elemFluxVarsCache)
const
101 FluxVariables fluxVars;
102 CellCenterPrimaryVariables flux = fluxVars.computeMassFlux(problem, element, fvGeometry, elemVolVars,
103 elemFaceVars, scvf, elemFluxVarsCache[scvf]);
105 EnergyLocalResidual::heatFlux(flux, problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf);
112 const Element &element,
113 const FVElementGeometry& fvGeometry,
114 const ElementVolumeVariables& elemVolVars,
115 const ElementFaceVariables& elemFaceVars,
116 const SubControlVolume &scv)
const
118 CellCenterPrimaryVariables result(0.0);
121 const auto sourceValues = problem.source(element, fvGeometry, elemVolVars, elemFaceVars, scv);
124 for (
int i = 0; i < result.size(); ++i)
125 result[i] = sourceValues[i + cellCenterOffset];
133 const SubControlVolume& scv,
134 const VolumeVariables& volVars)
const
136 CellCenterPrimaryVariables storage;
137 storage[Indices::conti0EqIdx - ModelTraits::dim()] = volVars.density();
139 EnergyLocalResidual::fluidPhaseStorage(storage, volVars);
146 const SubControlVolumeFace& scvf,
147 const VolumeVariables& volVars,
148 const ElementFaceVariables& elemFaceVars)
const
150 FacePrimaryVariables storage(0.0);
151 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
152 storage[0] = volVars.density() * velocity;
158 const Element& element,
159 const FVElementGeometry& fvGeometry,
160 const SubControlVolumeFace& scvf,
161 const ElementVolumeVariables& elemVolVars,
162 const ElementFaceVariables& elemFaceVars)
const
164 FacePrimaryVariables source(0.0);
165 const auto& insideVolVars = elemVolVars[scvf.insideScvIdx()];
166 source += problem.gravity()[scvf.directionIndex()] * insideVolVars.density();
167 source += problem.source(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[Indices::velocity(scvf.directionIndex())];
174 const Element& element,
175 const SubControlVolumeFace& scvf,
176 const FVElementGeometry& fvGeometry,
177 const ElementVolumeVariables& elemVolVars,
178 const ElementFaceVariables& elemFaceVars,
179 const ElementFluxVariablesCache& elemFluxVarsCache)
const
181 FluxVariables fluxVars;
182 return fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
189 const Element& element,
190 const FVElementGeometry& fvGeometry,
191 const SubControlVolumeFace& scvf,
192 const ElementVolumeVariables& elemVolVars,
193 const ElementFaceVariables& elemFaceVars,
194 const ElementBoundaryTypes& elemBcTypes,
195 const ElementFluxVariablesCache& elemFluxVarsCache)
const
197 CellCenterResidual result(0.0);
201 const auto bcTypes = problem.boundaryTypes(element, scvf);
204 if (bcTypes.isSymmetry())
208 result = computeFluxForCellCenter(problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf, elemFluxVarsCache);
211 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
212 if (bcTypes.hasNeumann())
214 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
215 const auto neumannFluxes = problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
217 for (
int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
219 if (bcTypes.isNeumann(eqIdx + cellCenterOffset))
220 result[eqIdx] = neumannFluxes[eqIdx + cellCenterOffset] * extrusionFactor * scvf.area();
225 incorporateWallFunction_(result, problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars);
234 const Problem& problem,
235 const Element& element,
236 const FVElementGeometry& fvGeometry,
237 const SubControlVolumeFace& scvf,
238 const ElementVolumeVariables& elemVolVars,
239 const ElementFaceVariables& elemFaceVars,
240 const ElementBoundaryTypes& elemBcTypes,
241 const ElementFluxVariablesCache& elemFluxVarsCache)
const
246 const auto bcTypes = problem.boundaryTypes(element, scvf);
248 if(bcTypes.isDirichlet(Indices::velocity(scvf.directionIndex())))
251 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
252 const Scalar dirichletValue = problem.dirichlet(element, scvf)[Indices::velocity(scvf.directionIndex())];
253 residual = velocity - dirichletValue;
255 else if(bcTypes.isSymmetry())
259 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
260 const Scalar fixedValue = 0.0;
261 residual = velocity - fixedValue;
270 const Element& element,
271 const FVElementGeometry& fvGeometry,
272 const SubControlVolumeFace& scvf,
273 const ElementVolumeVariables& elemVolVars,
274 const ElementFaceVariables& elemFaceVars,
275 const ElementBoundaryTypes& elemBcTypes,
276 const ElementFluxVariablesCache& elemFluxVarsCache)
const
278 FaceResidual result(0.0);
282 FluxVariables fluxVars;
285 const auto bcTypes = problem.boundaryTypes(element, scvf);
286 if (bcTypes.isNeumann(Indices::velocity(scvf.directionIndex())))
290 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
291 result = problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[Indices::velocity(scvf.directionIndex())]
292 * extrusionFactor * scvf.area();
295 result += fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
297 else if(bcTypes.isDirichlet(Indices::pressureIdx))
301 result = fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
304 result += fluxVars.inflowOutflowBoundaryFlux(problem, element, scvf, elemVolVars, elemFaceVars);
313 template<
class ...Args,
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<!turbulenceModel, int> = 0>
314 void incorporateWallFunction_(Args&&... args)
const
318 template<
bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<turbulenceModel,
int> = 0>
319 void incorporateWallFunction_(CellCenterResidual& boundaryFlux,
320 const Problem& problem,
321 const Element& element,
322 const FVElementGeometry& fvGeometry,
323 const SubControlVolumeFace& scvf,
324 const ElementVolumeVariables& elemVolVars,
325 const ElementFaceVariables& elemFaceVars)
const
327 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
328 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
331 for(
int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
334 if(problem.useWallFunction(element, scvf, eqIdx + cellCenterOffset))
336 boundaryFlux[eqIdx] = problem.wallFunction(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[eqIdx]
337 * extrusionFactor * scvf.area();
343 Implementation &asImp_()
344 {
return *
static_cast<Implementation *
>(
this); }
347 const Implementation &asImp_()
const
348 {
return *
static_cast<const Implementation *
>(
this); }
Calculates the element-wise residual for the staggered FV scheme.
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
Calculates the element-wise residual for the staggered FV scheme.
Definition: staggeredlocalresidual.hh:40
Definition: freeflow/navierstokes/localresidual.hh:35
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:111
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:233
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:132
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 boundary fluxes for a face dof.
Definition: freeflow/navierstokes/staggered/localresidual.hh:269
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:93
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:173
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:188
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:157
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:145
Definition: freeflow/nonisothermal/localresidual.hh:34
Declares all properties used in Dumux.