GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/multidomain/facet/couplingmanager.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 21 41 51.2%
Functions: 13 34 38.2%
Branches: 15 40 37.5%

Line Branch Exec Source
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 //
7 /*!
8 * \file
9 * \ingroup FacetCoupling
10 * \copydoc Dumux::FacetCouplingManager
11 */
12 #ifndef DUMUX_FACETCOUPLING_MANAGER_HH
13 #define DUMUX_FACETCOUPLING_MANAGER_HH
14
15 #include <dumux/common/properties.hh>
16 #include <dumux/common/indextraits.hh>
17
18 #include <dumux/discretization/method.hh>
19 #include <dumux/discretization/elementsolution.hh>
20 #include <dumux/discretization/evalsolution.hh>
21
22 namespace Dumux {
23 namespace FacetCoupling{
24
25 /*!
26 * \ingroup FacetCoupling
27 * \brief Free function that allows the creation of a volume variables object
28 * interpolated to a given position within an element. This is the standard
29 * implementation which simply interpolates the solution to the given position
30 * and then performs a volume variables update with the interpolated solution.
31 *
32 * \note This assumes element-wise constant parameters for the computation of secondary
33 * variables. For heteregeneous parameter distributions a default implementation
34 * cannot be defined and an adequate overload of this function has to be provided!
35 * \note For cell-centered schemes this is an unnecessary overhead because all variables
36 * are constant within the cells and a volume variables update can usually be realized
37 * more efficiently. This function is mainly to be used for the box scheme!
38 */
39 template<class VolumeVariables, class Problem, class SolutionVector, class FVGeometry>
40 922988 void makeInterpolatedVolVars(VolumeVariables& volVars,
41 const Problem& problem,
42 const SolutionVector& sol,
43 const FVGeometry& fvGeometry,
44 const typename FVGeometry::GridGeometry::GridView::template Codim<0>::Entity& element,
45 const typename FVGeometry::GridGeometry::GridView::template Codim<0>::Entity::Geometry& elemGeom,
46 const typename FVGeometry::GridGeometry::GridView::template Codim<0>::Entity::Geometry::GlobalCoordinate& pos)
47 {
48 // interpolate solution and set it for each entry in element solution
49 1801240 auto elemSol = elementSolution(element, sol, fvGeometry.gridGeometry());
50 1801240 const auto centerSol = evalSolution(element, elemGeom, fvGeometry.gridGeometry(), elemSol, pos);
51
4/4
✓ Branch 0 taken 1188728 times.
✓ Branch 1 taken 583180 times.
✓ Branch 2 taken 1188728 times.
✓ Branch 3 taken 583180 times.
4826716 for (unsigned int i = 0; i < fvGeometry.numScv(); ++i)
52 3781424 elemSol[i] = centerSol;
53
54 // Update volume variables with the interpolated solution. Note that this standard
55 // implementation only works for element-wise constant parameters as we simply use
56 // the first element scv for the vol var update. For heterogeneities within the element
57 // or more complex models (e.g. 2p with interface solver) a corresponding overload
58 // of this function has to be provided!
59 1845976 volVars.update(elemSol, problem, element, *scvs(fvGeometry).begin());
60 922988 }
61 } // end namespace FacetCoupling
62
63 /*!
64 * \ingroup FacetCoupling
65 * \brief Implementation for the coupling manager between two domains of dimension d
66 * and (d-1) for models considering coupling across the bulk domain element facets.
67 * The implementations are specificto the discretization method used in the bulk
68 * domain, which is extracted automatically from the grid geometry corresponding
69 * to the provided bulk domain id. Implementations for the different methods have
70 * to be provided and included at the end of this file.
71 *
72 * \tparam MDTraits The multidomain traits containing the types on all sub-domains
73 * \tparam CouplingMapper Class containing maps on the coupling between dofs of different grids
74 * \tparam bulkDomainId The domain id of the bulk problem
75 * \tparam lowDimDomainId The domain id of the lower-dimensional problem
76 * \tparam bulkDM Discretization method used in the bulk domain
77 */
78 template< class MDTraits,
79 class CouplingMapper,
80 std::size_t bulkDomainId = 0,
81 std::size_t lowDimDomainId = 1,
82 class DiscretizationMethod = typename GetPropType<typename MDTraits::template SubDomain<bulkDomainId>::TypeTag, Properties::GridGeometry>::DiscretizationMethod >
83 class FacetCouplingManager;
84
85 /*!
86 * \ingroup FacetCoupling
87 * \brief Class that handles the coupling between three sub-domains in models where
88 * the coupling between the two occurs across the facets of the d- and (d-1)-
89 * dimensional domains.
90 *
91 * \tparam MDTraits The multidomain traits containing the types on all sub-domains
92 * \tparam CouplingMapper Class containing maps on the coupling between dofs of different grids
93 * \tparam bulkDomainId The domain id of the d-dimensional bulk problem
94 * \tparam facetDomainId The domain id of the (d-1)-dimensional problem living on the bulk grid facets
95 * \tparam facetDomainId The domain id of the (d-2)-dimensional problem living on the bulk grid edges
96 */
97 template< class MDTraits,
98 class CouplingMapper,
99 std::size_t bulkDomainId = 0,
100 std::size_t facetDomainId = 1,
101 std::size_t edgeDomainId = 2 >
102 class FacetCouplingThreeDomainManager
103 : public FacetCouplingManager<MDTraits, CouplingMapper, bulkDomainId, facetDomainId>
104 , public FacetCouplingManager<MDTraits, CouplingMapper, facetDomainId, edgeDomainId>
105 {
106 using BulkFacetManager = FacetCouplingManager<MDTraits, CouplingMapper, bulkDomainId, facetDomainId>;
107 using FacetEdgeManager = FacetCouplingManager<MDTraits, CouplingMapper, facetDomainId, edgeDomainId>;
108
109 // convenience aliases and instances of the domain ids
110 using BulkIdType = typename MDTraits::template SubDomain<bulkDomainId>::Index;
111 using FacetIdType = typename MDTraits::template SubDomain<facetDomainId>::Index;
112 using EdgeIdType = typename MDTraits::template SubDomain<edgeDomainId>::Index;
113 static constexpr auto bulkId = BulkIdType();
114 static constexpr auto facetId = FacetIdType();
115 static constexpr auto edgeId = EdgeIdType();
116
117 // the sub-domain type tags
118 template<std::size_t id> using SubDomainTypeTag = typename MDTraits::template SubDomain<id>::TypeTag;
119
120 // further types specific to the sub-problems
121 template<std::size_t id> using PrimaryVariables = GetPropType<SubDomainTypeTag<id>, Properties::PrimaryVariables>;
122 template<std::size_t id> using Problem = GetPropType<SubDomainTypeTag<id>, Properties::Problem>;
123
124 template<std::size_t id> using GridGeometry = GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>;
125 template<std::size_t id> using FVElementGeometry = typename GridGeometry<id>::LocalView;
126 template<std::size_t id> using SubControlVolumeFace = typename GridGeometry<id>::SubControlVolumeFace;
127 template<std::size_t id> using GridView = typename GridGeometry<id>::GridView;
128 template<std::size_t id> using GridIndexType = typename IndexTraits<GridView<id>>::GridIndex;
129 template<std::size_t id> using Element = typename GridView<id>::template Codim<0>::Entity;
130
131 template<std::size_t id> using GridVariables = GetPropType<SubDomainTypeTag<id>, Properties::GridVariables>;
132 template<std::size_t id> using ElementVolumeVariables = typename GridVariables<id>::GridVolumeVariables::LocalView;
133 template<std::size_t id> using ElementFluxVariablesCache = typename GridVariables<id>::GridFluxVariablesCache::LocalView;
134
135 // helper function to check if a domain uses mpfa
136 template<std::size_t id>
137 static constexpr bool usesMpfa(Dune::index_constant<id> domainId)
138 { return GridGeometry<domainId>::discMethod == DiscretizationMethods::ccmpfa; }
139
140 public:
141 //! types used for coupling stencils
142 template<std::size_t i, std::size_t j>
143 using CouplingStencilType = typename std::conditional< (j == edgeDomainId),
144 typename FacetEdgeManager::template CouplingStencilType<i, j>,
145 typename BulkFacetManager::template CouplingStencilType<i, j> >::type;
146
147 //! the type of the solution vector
148 using SolutionVector = typename MDTraits::SolutionVector;
149
150 /*!
151 * \brief Initialize the coupling manager.
152 *
153 * \param bulkProblem The problem to be solved on the (3d) bulk domain
154 * \param facetProblem The problem to be solved on the (2d) facet domain
155 * \param edgeProblem The problem to be solved on the (1d) edge domain
156 * \param couplingMapper The mapper object containing the connectivity between the domains
157 * \param curSol The current solution
158 */
159 2 void init(std::shared_ptr< Problem<bulkId> > bulkProblem,
160 std::shared_ptr< Problem<facetId> > facetProblem,
161 std::shared_ptr< Problem<edgeId> > edgeProblem,
162 std::shared_ptr< CouplingMapper > couplingMapper,
163 const SolutionVector& curSol)
164 {
165
4/14
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
6 BulkFacetManager::init(bulkProblem, facetProblem, couplingMapper, curSol);
166
4/14
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
6 FacetEdgeManager::init(facetProblem, edgeProblem, couplingMapper, curSol);
167 2 }
168
169 //! Pull up functionalities from the parent classes
170 using BulkFacetManager::couplingStencil;
171 using FacetEdgeManager::couplingStencil;
172
173 using BulkFacetManager::isCoupled;
174 using FacetEdgeManager::isCoupled;
175
176 using BulkFacetManager::isOnInteriorBoundary;
177 using FacetEdgeManager::isOnInteriorBoundary;
178
179 using BulkFacetManager::getLowDimVolVars;
180 using FacetEdgeManager::getLowDimVolVars;
181
182 using BulkFacetManager::getLowDimElement;
183 using FacetEdgeManager::getLowDimElement;
184
185 using BulkFacetManager::getLowDimElementIndex;
186 using FacetEdgeManager::getLowDimElementIndex;
187
188 using BulkFacetManager::evalSourcesFromBulk;
189 using FacetEdgeManager::evalSourcesFromBulk;
190
191 using BulkFacetManager::evalCouplingResidual;
192 using FacetEdgeManager::evalCouplingResidual;
193
194 using BulkFacetManager::bindCouplingContext;
195 using FacetEdgeManager::bindCouplingContext;
196
197 using BulkFacetManager::updateCouplingContext;
198 using FacetEdgeManager::updateCouplingContext;
199
200 using BulkFacetManager::updateCoupledVariables;
201 using FacetEdgeManager::updateCoupledVariables;
202
203 using BulkFacetManager::extendJacobianPattern;
204 using FacetEdgeManager::extendJacobianPattern;
205
206 using BulkFacetManager::evalAdditionalDomainDerivatives;
207 using FacetEdgeManager::evalAdditionalDomainDerivatives;
208
209 // extension of the jacobian pattern for the facet domain only occurs
210 // within the bulk-facet coupling & for mpfa being used in the bulk domain.
211 template<class JacobianPattern>
212 void extendJacobianPattern(FacetIdType, JacobianPattern& pattern) const
213
1/4
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
2 { BulkFacetManager::extendJacobianPattern(facetId, pattern); }
214
215 template<class FacetLocalAssembler, class JacobianMatrixDiagBlock, class GridVariables>
216 void evalAdditionalDomainDerivatives(FacetIdType,
217 const FacetLocalAssembler& facetLocalAssembler,
218 const typename FacetLocalAssembler::LocalResidual::ElementResidualVector& origResiduals,
219 JacobianMatrixDiagBlock& A,
220 GridVariables& gridVariables)
221 616 { BulkFacetManager::evalAdditionalDomainDerivatives(facetId, facetLocalAssembler, origResiduals, A, gridVariables); }
222
223 /*!
224 * \brief The coupling stencil of the bulk with the edge domain (empty stencil).
225 */
226 const CouplingStencilType<bulkId, edgeId>& couplingStencil(BulkIdType domainI,
227 const Element<bulkId>& element,
228 EdgeIdType domainJ) const
229
2/4
✓ Branch 3 taken 1927 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1927 times.
✗ Branch 7 not taken.
34686 { return FacetEdgeManager::getEmptyStencil(edgeId); }
230
231 /*!
232 * \brief The coupling stencil of the edge with the bulk domain (empty stencil).
233 */
234 const CouplingStencilType<edgeId, bulkId>& couplingStencil(EdgeIdType domainI,
235 const Element<edgeId>& element,
236 BulkIdType domainJ) const
237 72 { return BulkFacetManager::getEmptyStencil(bulkId); }
238
239 /*!
240 * \brief updates the current solution. We have to overload this here
241 * to avoid ambiguity and update the solution in both managers
242 */
243 21 void updateSolution(const SolutionVector& sol)
244 {
245 21 BulkFacetManager::updateSolution(sol);
246 21 FacetEdgeManager::updateSolution(sol);
247 21 }
248
249 /*!
250 * \brief Interface for evaluating the coupling residual between the bulk and the edge domain.
251 * This is always zero as coupling only occurs between grids of codimension one. These
252 * overloads are provided by the two parent classes. However, we need this overload in
253 * order for the overload resolution not to fail.
254 */
255 template<std::size_t i,
256 std::size_t j,
257 class LocalAssembler,
258 std::enable_if_t<((i==bulkId && j==edgeId) || ((i==edgeId && j==bulkId))), int> = 0>
259 typename LocalAssembler::LocalResidual::ElementResidualVector
260 evalCouplingResidual(Dune::index_constant<i> domainI,
261 const LocalAssembler& localAssembler,
262 Dune::index_constant<j> domainJ,
263 GridIndexType<j> dofIdxGlobalJ)
264 {
265 typename LocalAssembler::LocalResidual::ElementResidualVector res(1);
266 res = 0.0;
267 return res;
268 }
269
270 /*!
271 * \brief Interface for binding the coupling context for the facet domain. In this case
272 * we have to bind both the facet -> bulk and the facet -> edge coupling context.
273 */
274 template< class Assembler >
275 void bindCouplingContext(FacetIdType, const Element<facetId>& element, const Assembler& assembler)
276 {
277 BulkFacetManager::bindCouplingContext(facetId, element, assembler);
278 FacetEdgeManager::bindCouplingContext(facetId, element, assembler);
279 }
280
281 /*!
282 * \brief Interface for updating the coupling context of the facet domain. In this case
283 * we have to update both the facet -> bulk and the facet -> edge coupling context.
284 */
285 template< class FacetLocalAssembler >
286 void updateCouplingContext(FacetIdType domainI,
287 const FacetLocalAssembler& facetLocalAssembler,
288 FacetIdType domainJ,
289 GridIndexType<facetId> dofIdxGlobalJ,
290 const PrimaryVariables<facetId>& priVarsJ,
291 unsigned int pvIdxJ)
292 {
293 BulkFacetManager::updateCouplingContext(domainI, facetLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ);
294 FacetEdgeManager::updateCouplingContext(domainI, facetLocalAssembler, domainJ, dofIdxGlobalJ, priVarsJ, pvIdxJ);
295 }
296
297 /*!
298 * \brief Interface for updating the coupling context between the bulk and the edge domain.
299 * We do nothing here because coupling only occurs between grids of codimension one.
300 * We have to provide this overload as the overload resolution fails otherwise.
301 */
302 template<std::size_t i,
303 std::size_t j,
304 class LocalAssembler,
305 std::enable_if_t<((i==bulkId && j==edgeId) || (i==edgeId && j==bulkId)), int> = 0>
306 void updateCouplingContext(Dune::index_constant<i> domainI,
307 const LocalAssembler& localAssembler,
308 Dune::index_constant<j> domainJ,
309 GridIndexType<j> dofIdxGlobalJ,
310 const PrimaryVariables<j>& priVarsJ,
311 unsigned int pvIdxJ)
312 { /*do nothing here*/ }
313
314 /*!
315 * \brief Interface for updating the local views of the facet domain after updateCouplingContext
316 * the coupling context. In this case we have to forward the both managers as the facet
317 * domain is a part in both.
318 */
319 template< class FacetLocalAssembler, class UpdatableElementVolVars, class UpdatableFluxVarCache>
320 void updateCoupledVariables(FacetIdType domainI,
321 const FacetLocalAssembler& facetLocalAssembler,
322 UpdatableElementVolVars& elemVolVars,
323 UpdatableFluxVarCache& elemFluxVarsCache)
324 {
325 5152 BulkFacetManager::updateCoupledVariables(domainI, facetLocalAssembler, elemVolVars, elemFluxVarsCache);
326 5152 FacetEdgeManager::updateCoupledVariables(domainI, facetLocalAssembler, elemVolVars, elemFluxVarsCache);
327 }
328
329 //! Return a const reference to bulk or facet problem
330 template<std::size_t id, std::enable_if_t<(id == bulkId || id == facetId), int> = 0>
331 const Problem<id>& problem() const { return BulkFacetManager::template problem<id>(); }
332
333 //! Return a reference to bulk or facet problem
334 template<std::size_t id, std::enable_if_t<(id == bulkId || id == facetId), int> = 0>
335 Problem<id>& problem() { return BulkFacetManager::template problem<id>(); }
336
337 //! Return a const reference to edge problem
338 template<std::size_t id, std::enable_if_t<(id == edgeId), int> = 0>
339 const Problem<id>& problem() const { return FacetEdgeManager::template problem<id>(); }
340
341 //! Return a reference to edge problem
342 template<std::size_t id, std::enable_if_t<(id == edgeId), int> = 0>
343 Problem<id>& problem() { return FacetEdgeManager::template problem<id>(); }
344 };
345
346 } // end namespace Dumux
347
348 // Here, we have to include all available implementations
349 #include <dumux/multidomain/facet/box/couplingmanager.hh>
350 #include <dumux/multidomain/facet/cellcentered/tpfa/couplingmanager.hh>
351 #include <dumux/multidomain/facet/cellcentered/mpfa/couplingmanager.hh>
352
353 #endif
354