3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
staggeredcouplingmanager.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 *****************************************************************************/
25#ifndef DUMUX_STAGGERED_COUPLING_MANAGER_HH
26#define DUMUX_STAGGERED_COUPLING_MANAGER_HH
27
32
33namespace Dumux {
34
40template<class MDTraits>
41class StaggeredCouplingManager: virtual public CouplingManager<MDTraits>
42{
44
45 template<std::size_t id> using GridGeometry = typename MDTraits::template SubDomain<id>::GridGeometry;
46 template<std::size_t id> using GridView = typename GridGeometry<id>::GridView;
47
48 using FVElementGeometry = typename GridGeometry<0>::LocalView;
49 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
50 using Element = typename GridView<0>::template Codim<0>::Entity;
51
52 using GridIndexType = typename IndexTraits< GridView<0> >::GridIndex;
53 using CouplingStencil = std::vector<GridIndexType>;
54
55public:
56
58
59 using Traits = MDTraits;
60
61 static constexpr auto cellCenterIdx = GridGeometry<0>::cellCenterIdx();
62 static constexpr auto faceIdx = GridGeometry<0>::faceIdx();
63
67 template<std::size_t i, std::size_t j, class LocalAssemblerI, class PriVarsJ>
68 void updateCouplingContext(Dune::index_constant<i> domainI,
69 const LocalAssemblerI& localAssemblerI,
70 Dune::index_constant<j> domainJ,
71 const std::size_t dofIdxGlobalJ,
72 const PriVarsJ& priVarsJ,
73 int pvIdxJ)
74 {
75 auto& curSol = this->curSol()[domainJ];
76
77 // only proceed if the solution vector has actually been set in the base class
78 // (which is not the case if the staggered model is not coupled to another domain)
79 if (curSol.size() > 0)
80 curSol[dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
81 }
82
93 const CouplingStencil& couplingStencil(Dune::index_constant<cellCenterIdx> domainI,
94 const Element& elementI,
95 Dune::index_constant<faceIdx> domainJ) const
96 {
97 const auto& connectivityMap = this->problem(domainI).gridGeometry().connectivityMap();
98 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(elementI);
99 return connectivityMap(domainI, domainJ, eIdx);
100 }
101
112 template<std::size_t i, std::size_t j>
113 const CouplingStencil couplingStencil(Dune::index_constant<i> domainI,
114 const SubControlVolumeFace& scvfI,
115 Dune::index_constant<j> domainJ) const
116 {
117 static_assert(i != j, "Domain i cannot be coupled to itself!");
118 static_assert(AlwaysFalse<Dune::index_constant<i>>::value,
119 "The coupling manager does not implement the couplingStencil() function" );
120
121 return CouplingStencil(); // supress compiler warning of function having no return statement
122 }
123
134 const CouplingStencil& couplingStencil(Dune::index_constant<faceIdx> domainI,
135 const SubControlVolumeFace& scvfI,
136 Dune::index_constant<cellCenterIdx> domainJ) const
137 {
138 const auto& connectivityMap = this->problem(domainI).gridGeometry().connectivityMap();
139 return connectivityMap(domainI, domainJ, scvfI.index());
140 }
141
148 template<class LocalAssemblerI, std::size_t j>
149 decltype(auto) evalCouplingResidual(Dune::index_constant<cellCenterIdx> domainI,
150 const LocalAssemblerI& localAssemblerI,
151 Dune::index_constant<j> domainJ,
152 std::size_t dofIdxGlobalJ) const
153 {
154 static_assert(domainI != domainJ, "Domain i cannot be coupled to itself!");
155 return localAssemblerI.evalLocalResidualForCellCenter();
156 }
157
174 template<class LocalAssemblerI, std::size_t j>
175 decltype(auto) evalCouplingResidual(Dune::index_constant<faceIdx> domainI,
176 const SubControlVolumeFace& scvfI,
177 const LocalAssemblerI& localAssemblerI,
178 Dune::index_constant<j> domainJ,
179 std::size_t dofIdxGlobalJ) const
180 {
181 static_assert(domainI != domainJ, "Domain i cannot be coupled to itself!");
182 return localAssemblerI.evalLocalResidualForFace(scvfI);
183 }
184
189 template<std::size_t i, typename std::enable_if_t<(GridGeometry<i>::discMethod != DiscretizationMethod::staggered), int> = 0>
190 decltype(auto) numericEpsilon(Dune::index_constant<i> id,
191 const std::string& paramGroup) const
192 {
193 return ParentType::numericEpsilon(id, paramGroup);
194 }
195
200 template<std::size_t i, typename std::enable_if_t<(GridGeometry<i>::discMethod == DiscretizationMethod::staggered), int> = 0>
201 decltype(auto) numericEpsilon(Dune::index_constant<i>,
202 const std::string& paramGroup) const
203 {
204 constexpr std::size_t numEqCellCenter = Traits::template SubDomain<cellCenterIdx>::PrimaryVariables::dimension;
205 constexpr std::size_t numEqFace = Traits::template SubDomain<faceIdx>::PrimaryVariables::dimension;
206 constexpr bool isCellCenter = GridGeometry<i>::isCellCenter();
207 constexpr std::size_t numEq = isCellCenter ? numEqCellCenter : numEqFace;
208 constexpr auto prefix = isCellCenter ? "CellCenter" : "Face";
209
210 try {
211 if(paramGroup.empty())
213 else
214 return NumericEpsilon<typename Traits::Scalar, numEq>(paramGroup + "." + prefix);
215 }
216 catch (Dune::RangeError& e)
217 {
218 DUNE_THROW(Dumux::ParameterException, "For the staggered model, you can to specify \n\n"
219 " CellCenter.Assembly.NumericDifference.PriVarMagnitude = mCC\n"
220 " Face.Assembly.NumericDifference.PriVarMagnitude = mFace\n"
221 " CellCenter.Assembly.NumericDifference.BaseEpsilon = eCC_0 ... eCC_numEqCellCenter-1\n"
222 " Face.Assembly.NumericDifference.BaseEpsilon = eFace_0 ... eFace_numEqFace-1\n\n"
223 "Wrong numer of values set for " << prefix << " (has " << numEq << " primary variable(s))\n\n" << e);
224 }
225
226 }
227
228};
229
230} //end namespace Dumux
231
232#endif
An adapter class for local assemblers using numeric differentiation.
Defines the index types used for grid and local indices.
The available discretization methods in Dumux.
decltype(auto) evalCouplingResidual(Dune::index_constant< faceIdx > domainI, const SubControlVolumeFace &scvfI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
evaluates the face residual of a coupled face of domain i which depends on the variables at the degre...
Definition: staggeredcouplingmanager.hh:175
decltype(auto) evalCouplingResidual(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
evaluates the element residual of a coupled element of domain i which depends on the variables at the...
Definition: multidomain/couplingmanager.hh:209
Definition: adapt.hh:29
A helper class for local assemblers using numeric differentiation to determine the epsilon.
Definition: numericepsilon.hh:41
Exception thrown if a run-time parameter is not specified correctly.
Definition: exceptions.hh:60
Struture to define the index types used for grid and local indices.
Definition: indextraits.hh:38
Template which always yields a false value.
Definition: typetraits.hh:36
Definition: multidomain/couplingmanager.hh:46
const Problem< i > & problem(Dune::index_constant< i > domainIdx) const
Return a reference to the sub problem.
Definition: multidomain/couplingmanager.hh:264
SolutionVector & curSol()
the solution vector of the coupled problem
Definition: multidomain/couplingmanager.hh:278
decltype(auto) numericEpsilon(Dune::index_constant< i >, const std::string &paramGroup) const
return the numeric epsilon used for deflecting primary variables of coupled domain i
Definition: multidomain/couplingmanager.hh:235
Base coupling manager for the staggered discretization.
Definition: staggeredcouplingmanager.hh:42
const CouplingStencil & couplingStencil(Dune::index_constant< faceIdx > domainI, const SubControlVolumeFace &scvfI, Dune::index_constant< cellCenterIdx > domainJ) const
returns an iteratable container of all indices of degrees of freedom of domain j that couple with / i...
Definition: staggeredcouplingmanager.hh:134
const CouplingStencil couplingStencil(Dune::index_constant< i > domainI, const SubControlVolumeFace &scvfI, Dune::index_constant< j > domainJ) const
returns an iteratable container of all indices of degrees of freedom of domain j that couple with / i...
Definition: staggeredcouplingmanager.hh:113
decltype(auto) evalCouplingResidual(Dune::index_constant< cellCenterIdx > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
evaluates the element residual of a coupled element of domain i which depends on the variables at the...
Definition: staggeredcouplingmanager.hh:149
const CouplingStencil & couplingStencil(Dune::index_constant< cellCenterIdx > domainI, const Element &elementI, Dune::index_constant< faceIdx > domainJ) const
returns an iteratable container of all indices of degrees of freedom of domain j that couple with / i...
Definition: staggeredcouplingmanager.hh:93
decltype(auto) numericEpsilon(Dune::index_constant< i > id, const std::string &paramGroup) const
return the numeric epsilon used for deflecting primary variables of coupled domain i.
Definition: staggeredcouplingmanager.hh:190
void updateCouplingContext(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, const std::size_t dofIdxGlobalJ, const PriVarsJ &priVarsJ, int pvIdxJ)
updates all data and variables that are necessary to evaluate the residual of the element of domain i...
Definition: staggeredcouplingmanager.hh:68
decltype(auto) numericEpsilon(Dune::index_constant< i >, const std::string &paramGroup) const
return the numeric epsilon used for deflecting primary variables of coupled domain i.
Definition: staggeredcouplingmanager.hh:201
static constexpr auto cellCenterIdx
Definition: staggeredcouplingmanager.hh:61
MDTraits Traits
Definition: staggeredcouplingmanager.hh:59
static constexpr auto faceIdx
Definition: staggeredcouplingmanager.hh:62
The interface of the coupling manager for multi domain problems.