version 3.10.0
Loading...
Searching...
No Matches
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// SPDX-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
13
14#ifndef DUMUX_STAGGERED_COUPLING_MANAGER_HH
15#define DUMUX_STAGGERED_COUPLING_MANAGER_HH
16
22
23namespace Dumux {
24
30template<class MDTraits>
31class StaggeredCouplingManager: virtual public CouplingManager<MDTraits>
32{
33 using ParentType = CouplingManager<MDTraits>;
34
35 template<std::size_t id> using GridGeometry = typename MDTraits::template SubDomain<id>::GridGeometry;
36 template<std::size_t id> using GridView = typename GridGeometry<id>::GridView;
37
38 using FVElementGeometry = typename GridGeometry<0>::LocalView;
39 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
40 using Element = typename GridView<0>::template Codim<0>::Entity;
41
42 using GridIndexType = typename IndexTraits< GridView<0> >::GridIndex;
43 using CouplingStencil = std::vector<GridIndexType>;
44
45public:
46
48
49 using Traits = MDTraits;
50
51 static constexpr auto cellCenterIdx = GridGeometry<0>::cellCenterIdx();
52 static constexpr auto faceIdx = GridGeometry<0>::faceIdx();
53
57 template<std::size_t i, std::size_t j, class LocalAssemblerI, class PriVarsJ>
58 void updateCouplingContext(Dune::index_constant<i> domainI,
59 const LocalAssemblerI& localAssemblerI,
60 Dune::index_constant<j> domainJ,
61 const std::size_t dofIdxGlobalJ,
62 const PriVarsJ& priVarsJ,
63 int pvIdxJ)
64 {
65 auto& curSol = this->curSol(domainJ);
66
67 // only proceed if the solution vector has actually been set in the base class
68 // (which is not the case if the staggered model is not coupled to another domain)
69 if (curSol.size() > 0)
70 curSol[dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
71 }
72
83 const CouplingStencil& couplingStencil(Dune::index_constant<cellCenterIdx> domainI,
84 const Element& elementI,
85 Dune::index_constant<faceIdx> domainJ) const
86 {
87 const auto& connectivityMap = this->problem(domainI).gridGeometry().connectivityMap();
88 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(elementI);
89 return connectivityMap(domainI, domainJ, eIdx);
90 }
91
102 template<std::size_t i, std::size_t j>
103 const CouplingStencil couplingStencil(Dune::index_constant<i> domainI,
104 const SubControlVolumeFace& scvfI,
105 Dune::index_constant<j> domainJ) const
106 {
107 static_assert(i != j, "Domain i cannot be coupled to itself!");
108 static_assert(AlwaysFalse<Dune::index_constant<i>>::value,
109 "The coupling manager does not implement the couplingStencil() function" );
110
111 return CouplingStencil(); // suppress compiler warning of function having no return statement
112 }
113
124 const CouplingStencil& couplingStencil(Dune::index_constant<faceIdx> domainI,
125 const SubControlVolumeFace& scvfI,
126 Dune::index_constant<cellCenterIdx> domainJ) const
127 {
128 const auto& connectivityMap = this->problem(domainI).gridGeometry().connectivityMap();
129 return connectivityMap(domainI, domainJ, scvfI.index());
130 }
131
138 template<class LocalAssemblerI, std::size_t j>
139 decltype(auto) evalCouplingResidual(Dune::index_constant<cellCenterIdx> domainI,
140 const LocalAssemblerI& localAssemblerI,
141 Dune::index_constant<j> domainJ,
142 std::size_t dofIdxGlobalJ) const
143 {
144 static_assert(domainI != domainJ, "Domain i cannot be coupled to itself!");
145 return localAssemblerI.evalLocalResidualForCellCenter();
146 }
147
163 template<class LocalAssemblerI, std::size_t j>
164 decltype(auto) evalCouplingResidual(Dune::index_constant<faceIdx> domainI,
165 const SubControlVolumeFace& scvfI,
166 const LocalAssemblerI& localAssemblerI,
167 Dune::index_constant<j> domainJ,
168 std::size_t dofIdxGlobalJ) const
169 {
170 static_assert(domainI != domainJ, "Domain i cannot be coupled to itself!");
171 return localAssemblerI.evalLocalResidualForFace(scvfI);
172 }
173
178 template<std::size_t i, typename std::enable_if_t<(GridGeometry<i>::discMethod != DiscretizationMethods::staggered), int> = 0>
179 decltype(auto) numericEpsilon(Dune::index_constant<i> id,
180 const std::string& paramGroup) const
181 {
182 return ParentType::numericEpsilon(id, paramGroup);
183 }
184
189 template<std::size_t i, typename std::enable_if_t<(GridGeometry<i>::discMethod == DiscretizationMethods::staggered), int> = 0>
190 decltype(auto) numericEpsilon(Dune::index_constant<i>,
191 const std::string& paramGroup) const
192 {
193 constexpr std::size_t numEqCellCenter = Traits::template SubDomain<cellCenterIdx>::PrimaryVariables::dimension;
194 constexpr std::size_t numEqFace = Traits::template SubDomain<faceIdx>::PrimaryVariables::dimension;
195 constexpr bool isCellCenter = GridGeometry<i>::isCellCenter();
196 constexpr std::size_t numEq = isCellCenter ? numEqCellCenter : numEqFace;
197 constexpr auto prefix = isCellCenter ? "CellCenter" : "Face";
198
199 try {
200 if(paramGroup.empty())
202 else
203 return NumericEpsilon<typename Traits::Scalar, numEq>(paramGroup + "." + prefix);
204 }
205 catch (Dune::RangeError& e)
206 {
207 DUNE_THROW(Dumux::ParameterException, "For the staggered model, you can to specify \n\n"
208 " CellCenter.Assembly.NumericDifference.PriVarMagnitude = mCC\n"
209 " Face.Assembly.NumericDifference.PriVarMagnitude = mFace\n"
210 " CellCenter.Assembly.NumericDifference.BaseEpsilon = eCC_0 ... eCC_numEqCellCenter-1\n"
211 " Face.Assembly.NumericDifference.BaseEpsilon = eFace_0 ... eFace_numEqFace-1\n\n"
212 "Wrong number of values set for " << prefix << " (has " << numEq << " primary variable(s))\n\n" << e);
213 }
214
215 }
216
217};
218
219} //end namespace Dumux
220
221#endif
const Problem< i > & problem(Dune::index_constant< i > domainIdx) const
Definition multidomain/couplingmanager.hh:297
SubSolutionVector< i > & curSol(Dune::index_constant< i > domainIdx)
Definition multidomain/couplingmanager.hh:326
decltype(auto) numericEpsilon(Dune::index_constant< i >, const std::string &paramGroup) const
Definition multidomain/couplingmanager.hh:263
CouplingManager()
Definition multidomain/couplingmanager.hh:70
A helper class for local assemblers using numeric differentiation to determine the epsilon.
Definition numericepsilon.hh:29
Exception thrown if a run-time parameter is not specified correctly.
Definition exceptions.hh:48
Base coupling manager for the staggered discretization.
Definition staggeredcouplingmanager.hh:32
const CouplingStencil & couplingStencil(Dune::index_constant< faceIdx > domainI, const SubControlVolumeFace &scvfI, Dune::index_constant< cellCenterIdx > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition staggeredcouplingmanager.hh:124
const CouplingStencil couplingStencil(Dune::index_constant< i > domainI, const SubControlVolumeFace &scvfI, Dune::index_constant< j > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition staggeredcouplingmanager.hh:103
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:139
const CouplingStencil & couplingStencil(Dune::index_constant< cellCenterIdx > domainI, const Element &elementI, Dune::index_constant< faceIdx > domainJ) const
returns an iterable container of all indices of degrees of freedom of domain j that couple with / inf...
Definition staggeredcouplingmanager.hh:83
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:179
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:58
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:164
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:190
static constexpr auto cellCenterIdx
Definition staggeredcouplingmanager.hh:51
MDTraits Traits
Definition staggeredcouplingmanager.hh:49
static constexpr auto faceIdx
Definition staggeredcouplingmanager.hh:52
Type traits.
decltype(auto) evalCouplingResidual(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ) const
Definition multidomain/couplingmanager.hh:237
Defines the index types used for grid and local indices.
The available discretization methods in Dumux.
The interface of the coupling manager for multi domain problems.
Definition adapt.hh:17
An adapter class for local assemblers using numeric differentiation.
Template which always yields a false value.
Definition common/typetraits/typetraits.hh:24
Structure to define the index types used for grid and local indices.
Definition indextraits.hh:26