version 3.10-dev
staggeredfvproblem.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// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
13#ifndef DUMUX_STAGGERD_FV_PROBLEM_HH
14#define DUMUX_STAGGERD_FV_PROBLEM_HH
15
16#include <dune/common/rangeutilities.hh>
17
21
22namespace Dumux {
23
34template<class TypeTag>
36{
38 using Implementation = GetPropType<TypeTag, Properties::Problem>;
40 using Element = typename GridView::template Codim<0>::Entity;
41
43 using GridVolumeVariables = typename GridVariables::GridVolumeVariables;
44 using ElementVolumeVariables = typename GridVolumeVariables::LocalView;
45 using GridFaceVariables = typename GridVariables::GridFaceVariables;
46 using ElementFaceVariables = typename GridFaceVariables::LocalView;
47
50 using FVElementGeometry = typename GridGeometry::LocalView;
51 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
52 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
54
55 using CoordScalar = typename GridView::ctype;
56 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
57
58 static constexpr auto cellCenterIdx = GridGeometry::cellCenterIdx();
59 static constexpr auto faceIdx = GridGeometry::faceIdx();
60
61 static constexpr auto numEqCellCenter = getPropValue<TypeTag, Properties::NumEqCellCenter>();
62 static constexpr auto numEqFace = getPropValue<TypeTag, Properties::NumEqFace>();
63
64public:
70 StaggeredFVProblem(std::shared_ptr<const GridGeometry> gridGeometry,
71 const std::string& paramGroup = "")
73 { }
74
83 bool isDirichletCell(const Element& element,
84 const FVElementGeometry& fvGeometry,
85 const SubControlVolume& scv,
86 int pvIdx) const
87 { return false; }
88
107 template<class ElementVolumeVariables, class ElementFaceVariables, class Entity>
108 NumEqVector source(const Element &element,
109 const FVElementGeometry& fvGeometry,
110 const ElementVolumeVariables& elemVolVars,
111 const ElementFaceVariables& elementFaceVars,
112 const Entity &e) const
113 {
114 // forward to solution independent, fully-implicit specific interface
115 return this->asImp_().sourceAtPos(e.center());
116 }
117
133 NumEqVector neumann(const Element& element,
134 const FVElementGeometry& fvGeometry,
135 const ElementVolumeVariables& elemVolVars,
136 const ElementFaceVariables& elemFaceVars,
137 const SubControlVolumeFace& scvf) const
138 {
139 // forward it to the interface with only the global position
140 return this->asImp_().neumannAtPos(scvf.ipGlobal());
141 }
142
149 template<class Entity>
150 PrimaryVariables initial(const Entity& entity) const
151 {
152 return this->asImp_().initialAtPos(entity.center());
153 }
154
159 template<class SolutionVector>
160 void applyInitialSolution(SolutionVector& sol) const
161 {
162 sol[cellCenterIdx].resize(this->gridGeometry().numCellCenterDofs());
163 sol[faceIdx].resize(this->gridGeometry().numFaceDofs());
164
165 auto fvGeometry = localView(this->gridGeometry());
166 for (const auto& element : elements(this->gridGeometry().gridView()))
167 {
168 fvGeometry.bindElement(element);
169
170 // loop over sub control volumes
171 for (auto&& scv : scvs(fvGeometry))
172 {
173 // let the problem do the dirty work of nailing down
174 // the initial solution.
175 auto initPriVars = this->asImp_().initial(scv);
176 this->asImp_().applyInitialCellCenterSolution(sol, scv, initPriVars);
177 }
178
179 // loop over faces
180 for(auto&& scvf : scvfs(fvGeometry))
181 {
182 auto initPriVars = this->asImp_().initial(scvf);
183 this->asImp_().applyInitialFaceSolution(sol, scvf, initPriVars);
184 }
185 }
186 }
187
188
190 template<class SolutionVector>
191 void applyInitialCellCenterSolution(SolutionVector& sol,
192 const SubControlVolume& scv,
193 const PrimaryVariables& initSol) const
194 {
195 // while the container within the actual solution vector holds numEqCellCenter
196 // elements, we need to specify an offset to get the correct entry of the initial solution
197 static constexpr auto offset = PrimaryVariables::dimension - numEqCellCenter;
198
199 for(int pvIdx = 0; pvIdx < numEqCellCenter; ++pvIdx)
200 sol[cellCenterIdx][scv.dofIndex()][pvIdx] = initSol[pvIdx + offset];
201 }
202
204 template<class SolutionVector>
205 void applyInitialFaceSolution(SolutionVector& sol,
206 const SubControlVolumeFace& scvf,
207 const PrimaryVariables& initSol) const
208 {
209 for(int pvIdx = 0; pvIdx < numEqFace; ++pvIdx)
210 sol[faceIdx][scvf.dofIndex()][pvIdx] = initSol[pvIdx];
211 }
212
213};
214
215} // end namespace Dumux
216
217#endif
Base class for all finite-volume problems.
Definition: common/fvproblem.hh:43
const std::string & paramGroup() const
The parameter group in which to retrieve runtime parameters.
Definition: common/fvproblem.hh:524
const GridGeometry & gridGeometry() const
The finite volume grid geometry.
Definition: common/fvproblem.hh:520
Implementation & asImp_()
Returns the implementation of the problem (i.e. static polymorphism)
Definition: common/fvproblem.hh:529
Base class for all finite-volume problems using spatial parameters.
Definition: fvproblemwithspatialparams.hh:29
Base class for all staggered finite-volume problems.
Definition: staggeredfvproblem.hh:36
NumEqVector neumann(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elemFaceVars, const SubControlVolumeFace &scvf) const
Evaluate the boundary conditions for a neumann boundary segment.
Definition: staggeredfvproblem.hh:133
bool isDirichletCell(const Element &element, const FVElementGeometry &fvGeometry, const SubControlVolume &scv, int pvIdx) const
Returns whether a fixed Dirichlet value shall be used at a given cell.
Definition: staggeredfvproblem.hh:83
StaggeredFVProblem(std::shared_ptr< const GridGeometry > gridGeometry, const std::string &paramGroup="")
Constructor.
Definition: staggeredfvproblem.hh:70
PrimaryVariables initial(const Entity &entity) const
Evaluate the initial value for an element (for cell-centered primary variables) or face (for velociti...
Definition: staggeredfvproblem.hh:150
void applyInitialFaceSolution(SolutionVector &sol, const SubControlVolumeFace &scvf, const PrimaryVariables &initSol) const
Applies the initial face solution.
Definition: staggeredfvproblem.hh:205
void applyInitialCellCenterSolution(SolutionVector &sol, const SubControlVolume &scv, const PrimaryVariables &initSol) const
Applies the initial cell center solution.
Definition: staggeredfvproblem.hh:191
void applyInitialSolution(SolutionVector &sol) const
Applies the initial solution for all degrees of freedom of the grid.
Definition: staggeredfvproblem.hh:160
NumEqVector source(const Element &element, const FVElementGeometry &fvGeometry, const ElementVolumeVariables &elemVolVars, const ElementFaceVariables &elementFaceVars, const Entity &e) const
Evaluate the source term for all phases within a given sub-control-volume (-face).
Definition: staggeredfvproblem.hh:108
Defines all properties used in Dumux.
Base class for all finite volume problems that are parameterized.
typename NumEqVectorTraits< PrimaryVariables >::type NumEqVector
A vector with the same size as numbers of equations This is the default implementation and has to be ...
Definition: numeqvector.hh:34
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:26
typename GetProp< TypeTag, Property >::type GetPropType
get the type alias defined in the property
Definition: propertysystem.hh:296
Definition: adapt.hh:17
A helper to deduce a vector with the same size as numbers of equations.