3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
facecentered/staggered/elementsolution.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_DISCRETIZATION_FACECENTERED_STAGGERED_ELEMENT_SOLUTION_HH
25#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_ELEMENT_SOLUTION_HH
26
27#include <type_traits>
28#include <array>
29#include <utility>
30
33
34namespace Dumux {
35
40template<class FVElementGeometry, class PV>
42{
43 using GridGeometry = typename FVElementGeometry::GridGeometry;
44 using GridView = typename GridGeometry::GridView;
45 using Element = typename GridView::template Codim<0>::Entity;
46 using SmallLocalIndexType = typename IndexTraits<GridView>::SmallLocalIndex;
47
48 static constexpr auto numScvsPerElement = GridView::Grid::dimension * 2;
49
50public:
52 using PrimaryVariables = PV;
53
55
57 template<class SolutionVector>
59 const SolutionVector& sol,
60 const GridGeometry& gridGeometry)
61 {
62 const auto fvGeometry = localView(gridGeometry).bindElement(element);
63
64 for (const auto& scv : scvs(fvGeometry))
65 priVars_[scv.indexInElement()] = sol[scv.dofIndex()];
66 }
67
69 template<class ElementVolumeVariables>
71 const ElementVolumeVariables& elemVolVars,
72 const FVElementGeometry& fvGeometry)
73 {
74 for (const auto& scv : scvs(fvGeometry))
75 priVars_ = elemVolVars[scv].priVars();
76 }
77
80 {
81 priVars_[0] = std::move(priVars);
82 for (int i = 1; i < priVars_.size(); ++i)
83 priVars_[i] = priVars_[0];
84 }
85
88 {
89 priVars_[0] = priVars;
90 for (int i = 1; i < priVars_.size(); ++i)
91 priVars_[i] = priVars_[0];
92 }
93
95 template<class SolutionVector>
96 void update(const Element& element, const SolutionVector& sol,
97 const GridGeometry& gridGeometry)
98 {
99 const auto fvGeometry = localView(gridGeometry).bindElement(element);
100
101 for (const auto& scv : scvs(fvGeometry))
102 priVars_[scv.indexInElement()] = sol[scv.dofIndex()];
103 }
104
106 const PrimaryVariables& operator [](SmallLocalIndexType localScvIdx) const
107 {
108 return priVars_[localScvIdx];
109 }
110
112 PrimaryVariables& operator [](SmallLocalIndexType localScvIdx)
113 {
114 return priVars_[localScvIdx];
115 }
116
118 static constexpr std::size_t size()
119 { return numScvsPerElement; }
120
121private:
122 std::array<PrimaryVariables, numScvsPerElement> priVars_;
123};
124
129template<class Element, class SolutionVector, class GridGeometry>
130auto elementSolution(const Element& element, const SolutionVector& sol, const GridGeometry& gg)
131-> std::enable_if_t<
132 GridGeometry::discMethod == DiscretizationMethods::fcstaggered,
133 FaceCenteredStaggeredElementSolution<
134 typename GridGeometry::LocalView,
135 std::decay_t<decltype(std::declval<SolutionVector>()[0])>
136 >
137>
138{ return { element, sol, gg }; }
139
144template<class Element, class ElementVolumeVariables, class FVElementGeometry>
145auto elementSolution(const Element& element, const ElementVolumeVariables& elemVolVars, const FVElementGeometry& gg)
146-> std::enable_if_t<
147 FVElementGeometry::GridGeometry::discMethod == DiscretizationMethods::fcstaggered,
148 FaceCenteredStaggeredElementSolution<
149 FVElementGeometry,
150 typename ElementVolumeVariables::VolumeVariables::PrimaryVariables
151 >
152>
153{ return { element, elemVolVars, gg }; }
154
160template<class FVElementGeometry, class PrimaryVariables>
161auto elementSolution(PrimaryVariables&& priVars)
162-> std::enable_if_t<
163 FVElementGeometry::GridGeometry::discMethod == DiscretizationMethods::fcstaggered,
164 FaceCenteredStaggeredElementSolution<
165 FVElementGeometry,
166 std::decay_t<PrimaryVariables>
167 >
168>
169{ return { std::forward<PrimaryVariables>(priVars) }; }
170
171} // end namespace Dumux
172
173#endif
Defines the index types used for grid and local indices.
The available discretization methods in Dumux.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
auto elementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gg) -> std::enable_if_t< GridGeometry::discMethod==DiscretizationMethods::box, BoxElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for box schemes.
Definition: box/elementsolution.hh:118
Definition: adapt.hh:29
constexpr FCStaggered fcstaggered
Definition: method.hh:142
std::uint_least8_t SmallLocalIndex
Definition: indextraits.hh:41
The global face variables class for staggered models.
Definition: facecentered/staggered/elementsolution.hh:42
FaceCenteredStaggeredElementSolution(PrimaryVariables &&priVars)
Constructor with a primary variable object.
Definition: facecentered/staggered/elementsolution.hh:79
FaceCenteredStaggeredElementSolution(const Element &element, const ElementVolumeVariables &elemVolVars, const FVElementGeometry &fvGeometry)
Constructor with element, element volume variables and fv element geometry.
Definition: facecentered/staggered/elementsolution.hh:70
FaceCenteredStaggeredElementSolution(const Element &element, const SolutionVector &sol, const GridGeometry &gridGeometry)
Constructor with element, solution vector and grid geometry.
Definition: facecentered/staggered/elementsolution.hh:58
PV PrimaryVariables
export the primary variables type
Definition: facecentered/staggered/elementsolution.hh:52
const PrimaryVariables & operator[](SmallLocalIndexType localScvIdx) const
bracket operator const access
Definition: facecentered/staggered/elementsolution.hh:106
static constexpr std::size_t size()
return the size of the element solution
Definition: facecentered/staggered/elementsolution.hh:118
FaceCenteredStaggeredElementSolution(const PrimaryVariables &priVars)
Constructor with a primary variable object.
Definition: facecentered/staggered/elementsolution.hh:87
void update(const Element &element, const SolutionVector &sol, const GridGeometry &gridGeometry)
extract the element solution from the solution vector using a mapper
Definition: facecentered/staggered/elementsolution.hh:96