3.1-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
freeflow/navierstokes/staggered/localresidual.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*****************************************************************************
4 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
25#define DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
26
30#include <dune/common/hybridutilities.hh>
32
33namespace Dumux {
34
35// forward declaration
36template<class TypeTag, DiscretizationMethod discMethod>
38
43template<class TypeTag>
45: public StaggeredLocalResidual<TypeTag>
46{
48 friend class StaggeredLocalResidual<TypeTag>;
49
51
52 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
53 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
54 using VolumeVariables = typename GridVolumeVariables::VolumeVariables;
55
56 using GridFluxVariablesCache = typename GridVariables::GridFluxVariablesCache;
57 using ElementFluxVariablesCache = typename GridFluxVariablesCache::LocalView;
58
59 using GridFaceVariables = typename GridVariables::GridFaceVariables;
60 using ElementFaceVariables = typename GridFaceVariables::LocalView;
61
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;
72 using CellCenterPrimaryVariables = GetPropType<TypeTag, Properties::CellCenterPrimaryVariables>;
76
77 using CellCenterResidual = CellCenterPrimaryVariables;
78 using FaceResidual = FacePrimaryVariables;
79
81
82public:
83 using EnergyLocalResidual = FreeFlowEnergyLocalResidual<GridGeometry, FluxVariables, ModelTraits::enableEnergyBalance(), (ModelTraits::numFluidComponents() > 1)>;
84
85 // account for the offset of the cell center privars within the PrimaryVariables container
86 static constexpr auto cellCenterOffset = ModelTraits::numEq() - CellCenterPrimaryVariables::dimension;
87 static_assert(cellCenterOffset == ModelTraits::dim(), "cellCenterOffset must equal dim for staggered NavierStokes");
88
90 using ParentType::ParentType;
91
93 CellCenterPrimaryVariables computeFluxForCellCenter(const Problem& problem,
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
100 {
101 FluxVariables fluxVars;
102 CellCenterPrimaryVariables flux = fluxVars.computeMassFlux(problem, element, fvGeometry, elemVolVars,
103 elemFaceVars, scvf, elemFluxVarsCache[scvf]);
104
105 EnergyLocalResidual::heatFlux(flux, problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf);
106
107 return flux;
108 }
109
111 CellCenterPrimaryVariables computeSourceForCellCenter(const Problem& problem,
112 const Element &element,
113 const FVElementGeometry& fvGeometry,
114 const ElementVolumeVariables& elemVolVars,
115 const ElementFaceVariables& elemFaceVars,
116 const SubControlVolume &scv) const
117 {
118 CellCenterPrimaryVariables result(0.0);
119
120 // get the values from the problem
121 const auto sourceValues = problem.source(element, fvGeometry, elemVolVars, elemFaceVars, scv);
122
123 // copy the respective cell center related values to the result
124 for (int i = 0; i < result.size(); ++i)
125 result[i] = sourceValues[i + cellCenterOffset];
126
127 return result;
128 }
129
130
132 CellCenterPrimaryVariables computeStorageForCellCenter(const Problem& problem,
133 const SubControlVolume& scv,
134 const VolumeVariables& volVars) const
135 {
136 CellCenterPrimaryVariables storage;
137 storage[Indices::conti0EqIdx - ModelTraits::dim()] = volVars.density();
138
139 EnergyLocalResidual::fluidPhaseStorage(storage, volVars);
140
141 return storage;
142 }
143
145 FacePrimaryVariables computeStorageForFace(const Problem& problem,
146 const SubControlVolumeFace& scvf,
147 const VolumeVariables& volVars,
148 const ElementFaceVariables& elemFaceVars) const
149 {
150 FacePrimaryVariables storage(0.0);
151 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
152 storage[0] = volVars.density() * velocity;
153 return storage;
154 }
155
157 FacePrimaryVariables computeSourceForFace(const Problem& problem,
158 const Element& element,
159 const FVElementGeometry& fvGeometry,
160 const SubControlVolumeFace& scvf,
161 const ElementVolumeVariables& elemVolVars,
162 const ElementFaceVariables& elemFaceVars) const
163 {
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())];
168
169 return source;
170 }
171
173 FacePrimaryVariables computeFluxForFace(const Problem& problem,
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
180 {
181 FluxVariables fluxVars;
182 return fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
183 }
184
188 CellCenterResidual computeBoundaryFluxForCellCenter(const Problem& problem,
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
196 {
197 CellCenterResidual result(0.0);
198
199 if (scvf.boundary())
200 {
201 const auto bcTypes = problem.boundaryTypes(element, scvf);
202
203 // no fluxes occur over symmetry boundaries
204 if (bcTypes.isSymmetry())
205 return result;
206
207 // treat Dirichlet and outflow BCs
208 result = computeFluxForCellCenter(problem, element, fvGeometry, elemVolVars, elemFaceVars, scvf, elemFluxVarsCache);
209
210 // treat Neumann BCs, i.e. overwrite certain fluxes by user-specified values
211 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
212 if (bcTypes.hasNeumann())
213 {
214 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
215 const auto neumannFluxes = problem.neumann(element, fvGeometry, elemVolVars, elemFaceVars, scvf);
216
217 for (int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
218 {
219 if (bcTypes.isNeumann(eqIdx + cellCenterOffset))
220 result[eqIdx] = neumannFluxes[eqIdx + cellCenterOffset] * extrusionFactor * scvf.area();
221 }
222 }
223
224 // account for wall functions, if used
225 incorporateWallFunction_(result, problem, element, fvGeometry, scvf, elemVolVars, elemFaceVars);
226 }
227 return result;
228 }
229
233 void evalDirichletBoundariesForFace(FaceResidual& residual,
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
242 {
243 if (scvf.boundary())
244 {
245 // handle the actual boundary conditions:
246 const auto bcTypes = problem.boundaryTypes(element, scvf);
247
248 if(bcTypes.isDirichlet(Indices::velocity(scvf.directionIndex())))
249 {
250 // set a fixed value for the velocity for Dirichlet boundary conditions
251 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
252 const Scalar dirichletValue = problem.dirichlet(element, scvf)[Indices::velocity(scvf.directionIndex())];
253 residual = velocity - dirichletValue;
254 }
255 else if(bcTypes.isSymmetry())
256 {
257 // for symmetry boundary conditions, there is no flow accross the boundary and
258 // we therefore treat it like a Dirichlet boundary conditions with zero velocity
259 const Scalar velocity = elemFaceVars[scvf].velocitySelf();
260 const Scalar fixedValue = 0.0;
261 residual = velocity - fixedValue;
262 }
263 }
264 }
265
269 FaceResidual computeBoundaryFluxForFace(const Problem& problem,
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
277 {
278 FaceResidual result(0.0);
279
280 if (scvf.boundary())
281 {
282 FluxVariables fluxVars;
283
284 // handle the actual boundary conditions:
285 const auto bcTypes = problem.boundaryTypes(element, scvf);
286 if (bcTypes.isNeumann(Indices::velocity(scvf.directionIndex())))
287 {
288 // the source term has already been accounted for, here we
289 // add a given Neumann flux for the face on the boundary itself ...
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();
293
294 // ... and treat the fluxes of the remaining (frontal and lateral) faces of the staggered control volume
295 result += fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
296 }
297 else if(bcTypes.isDirichlet(Indices::pressureIdx))
298 {
299 // we are at an "fixed pressure" boundary for which the resdiual of the momentum balance needs to be assembled
300 // as if it where inside the domain and not on the boundary (source term has already been acounted for)
301 result = fluxVars.computeMomentumFlux(problem, element, scvf, fvGeometry, elemVolVars, elemFaceVars, elemFluxVarsCache.gridFluxVarsCache());
302
303 // incorporate the inflow or outflow contribution
304 result += fluxVars.inflowOutflowBoundaryFlux(problem, element, scvf, elemVolVars, elemFaceVars);
305 }
306 }
307 return result;
308 }
309
310private:
311
313 template<class ...Args, bool turbulenceModel = ModelTraits::usesTurbulenceModel(), std::enable_if_t<!turbulenceModel, int> = 0>
314 void incorporateWallFunction_(Args&&... args) const
315 {}
316
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
326 {
327 static constexpr auto numEqCellCenter = CellCenterResidual::dimension;
328 const auto extrusionFactor = elemVolVars[scvf.insideScvIdx()].extrusionFactor();
329
330 // account for wall functions, if used
331 for(int eqIdx = 0; eqIdx < numEqCellCenter; ++eqIdx)
332 {
333 // use a wall function
334 if(problem.useWallFunction(element, scvf, eqIdx + cellCenterOffset))
335 {
336 boundaryFlux[eqIdx] = problem.wallFunction(element, fvGeometry, elemVolVars, elemFaceVars, scvf)[eqIdx]
337 * extrusionFactor * scvf.area();
338 }
339 }
340 }
341
343 Implementation &asImp_()
344 { return *static_cast<Implementation *>(this); }
345
347 const Implementation &asImp_() const
348 { return *static_cast<const Implementation *>(this); }
349};
350
351} // end namespace Dumux
352
353#endif // DUMUX_STAGGERED_NAVIERSTOKES_LOCAL_RESIDUAL_HH
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.