3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
multidomain/couplingmanager.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_MULTIDOMAIN_COUPLING_MANAGER_HH
26#define DUMUX_MULTIDOMAIN_COUPLING_MANAGER_HH
27
28#include <memory>
29#include <tuple>
30#include <vector>
31#include <dune/common/exceptions.hh>
32#include <dune/common/indices.hh>
36
37namespace Dumux {
38
44template<class Traits>
46{
47 template<std::size_t id> using SubDomainTypeTag = typename Traits::template SubDomain<id>::TypeTag;
48 template<std::size_t id> using PrimaryVariables = GetPropType<SubDomainTypeTag<id>, Properties::PrimaryVariables>;
49 template<std::size_t id> using GridView = typename GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>::GridView;
50 template<std::size_t id> using Element = typename GridView<id>::template Codim<0>::Entity;
51 template<std::size_t id> using Problem = GetPropType<SubDomainTypeTag<id>, Properties::Problem>;
52 template<std::size_t id> using ProblemWeakPtr = std::weak_ptr<const Problem<id>>;
53 using Problems = typename Traits::template Tuple<ProblemWeakPtr>;
54
55public:
57 template<std::size_t i, std::size_t j>
58 using CouplingStencilType = std::vector<std::size_t>;
59
61 using SolutionVector = typename Traits::SolutionVector;
62
66 // \{
67
82 template<std::size_t i, std::size_t j>
83 const CouplingStencilType<i, j>& couplingStencil(Dune::index_constant<i> domainI,
84 const Element<i>& elementI,
85 Dune::index_constant<j> domainJ) const
86 {
87 static_assert(i != j, "Domain i cannot be coupled to itself!");
88 static_assert(AlwaysFalse<Dune::index_constant<i>>::value,
89 "The coupling manager does not implement the couplingStencil() function");
90 }
91
100 template<std::size_t id, class JacobianPattern>
101 void extendJacobianPattern(Dune::index_constant<id> domainI, JacobianPattern& pattern) const
102 {}
103
104 // \}
105
109 // \{
110
125 template<std::size_t i, class Assembler>
126 void bindCouplingContext(Dune::index_constant<i> domainI,
127 const Element<i>& elementI,
128 const Assembler& assembler)
129 {}
130
131
151 template<std::size_t i, std::size_t j, class LocalAssemblerI>
152 void updateCouplingContext(Dune::index_constant<i> domainI,
153 const LocalAssemblerI& localAssemblerI,
154 Dune::index_constant<j> domainJ,
155 std::size_t dofIdxGlobalJ,
156 const PrimaryVariables<j>& priVarsJ,
157 int pvIdxJ)
158 {
159 curSol_[domainJ][dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
160 }
161
176 template<std::size_t i, class LocalAssemblerI, class UpdatableElementVolVars, class UpdatableFluxVarCache>
177 void updateCoupledVariables(Dune::index_constant<i> domainI,
178 const LocalAssemblerI& localAssemblerI,
179 UpdatableElementVolVars& elemVolVars,
180 UpdatableFluxVarCache& elemFluxVarsCache)
181 {}
182
187 { curSol_ = curSol; }
188
189 // \}
190
208 template<std::size_t i, std::size_t j, class LocalAssemblerI>
209 decltype(auto) evalCouplingResidual(Dune::index_constant<i> domainI,
210 const LocalAssemblerI& localAssemblerI,
211 Dune::index_constant<j> domainJ,
212 std::size_t dofIdxGlobalJ) const
213 {
214 return localAssemblerI.evalLocalResidual();
215 }
216
223 template<std::size_t i, class LocalAssemblerI, class JacobianMatrixDiagBlock, class GridVariables>
224 void evalAdditionalDomainDerivatives(Dune::index_constant<i> domainI,
225 const LocalAssemblerI& localAssemblerI,
226 const typename LocalAssemblerI::LocalResidual::ElementResidualVector& origResiduals,
227 JacobianMatrixDiagBlock& A,
228 GridVariables& gridVariables)
229 {}
230
234 template<std::size_t i>
235 decltype(auto) numericEpsilon(Dune::index_constant<i>,
236 const std::string& paramGroup) const
237 {
238 constexpr auto numEq = PrimaryVariables<i>::dimension;
240 }
241
246 template<typename... SubProblems>
247 void setSubProblems(const std::tuple<std::shared_ptr<SubProblems>...>& problems)
248 { problems_ = problems; }
249
255 template<class SubProblem, std::size_t i>
256 void setSubProblem(std::shared_ptr<SubProblem> problem, Dune::index_constant<i> domainIdx)
257 { std::get<i>(problems_) = problem; }
258
263 template<std::size_t i>
264 const Problem<i>& problem(Dune::index_constant<i> domainIdx) const
265 {
266 if (!std::get<i>(problems_).expired())
267 return *std::get<i>(problems_).lock();
268 else
269 DUNE_THROW(Dune::InvalidStateException, "The problem pointer was not set or has already expired. Use setSubProblems() before calling this function");
270 }
271
272protected:
273
279 { return curSol_; }
280
285 const SolutionVector& curSol() const
286 { return curSol_; }
287
288private:
293 SolutionVector curSol_;
294
299 Problems problems_;
300
301};
302
303} //end namespace Dumux
304
305#endif
An adapter class for local assemblers using numeric differentiation.
void updateCouplingContext(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ, const PrimaryVariables< j > &priVarsJ, int pvIdxJ)
updates all data and variables that are necessary to evaluate the residual of the element of domain i...
Definition: multidomain/couplingmanager.hh:152
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
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
A helper class for local assemblers using numeric differentiation to determine the epsilon.
Definition: numericepsilon.hh:41
A vector of primary variables.
Definition: common/properties.hh:49
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:57
Definition: common/properties.hh:101
Template which always yields a false value.
Definition: typetraits.hh:36
Definition: multidomain/couplingmanager.hh:46
void updateCoupledVariables(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, UpdatableElementVolVars &elemVolVars, UpdatableFluxVarCache &elemFluxVarsCache)
update variables of domain i that depend on variables in domain j after the coupling context has been...
Definition: multidomain/couplingmanager.hh:177
void extendJacobianPattern(Dune::index_constant< id > domainI, JacobianPattern &pattern) const
extend the jacobian pattern of the diagonal block of domain i by those entries that are not already i...
Definition: multidomain/couplingmanager.hh:101
void setSubProblem(std::shared_ptr< SubProblem > problem, Dune::index_constant< i > domainIdx)
set a pointer to one of the sub problems
Definition: multidomain/couplingmanager.hh:256
void setSubProblems(const std::tuple< std::shared_ptr< SubProblems >... > &problems)
set the pointers to the sub problems
Definition: multidomain/couplingmanager.hh:247
const Problem< i > & problem(Dune::index_constant< i > domainIdx) const
Return a reference to the sub problem.
Definition: multidomain/couplingmanager.hh:264
void bindCouplingContext(Dune::index_constant< i > domainI, const Element< i > &elementI, const Assembler &assembler)
prepares all data and variables that are necessary to evaluate the residual of the element of domain ...
Definition: multidomain/couplingmanager.hh:126
std::vector< std::size_t > CouplingStencilType
default type used for coupling element stencils
Definition: multidomain/couplingmanager.hh:58
SolutionVector & curSol()
the solution vector of the coupled problem
Definition: multidomain/couplingmanager.hh:278
void updateSolution(const SolutionVector &curSol)
Updates the entire solution vector, e.g. before assembly or after grid adaption.
Definition: multidomain/couplingmanager.hh:186
void evalAdditionalDomainDerivatives(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, const typename LocalAssemblerI::LocalResidual::ElementResidualVector &origResiduals, JacobianMatrixDiagBlock &A, GridVariables &gridVariables)
evaluate additional derivatives of the element residual of a domain with respect to dofs in the same ...
Definition: multidomain/couplingmanager.hh:224
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
const CouplingStencilType< i, j > & couplingStencil(Dune::index_constant< i > domainI, const Element< i > &elementI, 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: multidomain/couplingmanager.hh:83
const SolutionVector & curSol() const
the solution vector of the coupled problem
Definition: multidomain/couplingmanager.hh:285
typename Traits::SolutionVector SolutionVector
the type of the solution vector
Definition: multidomain/couplingmanager.hh:61
Declares all properties used in Dumux.