3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
multidomain/boundary/stokesdarcy/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_STOKES_DARCY_COUPLINGMANAGER_HH
26#define DUMUX_STOKES_DARCY_COUPLINGMANAGER_HH
27
28#include <utility>
29#include <memory>
30
31#include <dune/common/float_cmp.hh>
32#include <dune/common/exceptions.hh>
36
37#include "couplingdata.hh"
38#include "couplingmapper.hh"
39
40namespace Dumux {
41
46template<class MDTraits>
48: public StaggeredCouplingManager<MDTraits>
49{
50 using Scalar = typename MDTraits::Scalar;
52
53public:
54 static constexpr auto stokesFaceIdx = typename MDTraits::template SubDomain<0>::Index();
55 static constexpr auto stokesCellCenterIdx = typename MDTraits::template SubDomain<1>::Index();
56 static constexpr auto stokesIdx = stokesFaceIdx;
57 static constexpr auto darcyIdx = typename MDTraits::template SubDomain<2>::Index();
58
59private:
60
61 using SolutionVector = typename MDTraits::SolutionVector;
62
63 // obtain the type tags of the sub problems
64 using StokesTypeTag = typename MDTraits::template SubDomain<0>::TypeTag;
65 using DarcyTypeTag = typename MDTraits::template SubDomain<2>::TypeTag;
66
67 using CouplingStencils = std::unordered_map<std::size_t, std::vector<std::size_t> >;
68 using CouplingStencil = CouplingStencils::mapped_type;
69
70 // the sub domain type tags
71 template<std::size_t id>
72 using SubDomainTypeTag = typename MDTraits::template SubDomain<id>::TypeTag;
73
74 static constexpr bool isCompositional = GetPropType<SubDomainTypeTag<0>, Properties::ModelTraits>::numFluidComponents()> 1;
75
76 template<std::size_t id> using GridView = typename GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>::GridView;
77 template<std::size_t id> using Problem = GetPropType<SubDomainTypeTag<id>, Properties::Problem>;
78 template<std::size_t id> using NumEqVector = GetPropType<SubDomainTypeTag<id>, Properties::NumEqVector>;
79 template<std::size_t id> using ElementVolumeVariables = typename GetPropType<SubDomainTypeTag<id>, Properties::GridVolumeVariables>::LocalView;
80 template<std::size_t id> using GridVolumeVariables = GetPropType<SubDomainTypeTag<id>, Properties::GridVolumeVariables>;
81 template<std::size_t id> using VolumeVariables = typename GetPropType<SubDomainTypeTag<id>, Properties::GridVolumeVariables>::VolumeVariables;
82 template<std::size_t id> using GridGeometry = GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>;
83 template<std::size_t id> using FVElementGeometry = typename GridGeometry<id>::LocalView;
84 template<std::size_t id> using ElementBoundaryTypes = GetPropType<SubDomainTypeTag<id>, Properties::ElementBoundaryTypes>;
85 template<std::size_t id> using ElementFluxVariablesCache = typename GetPropType<SubDomainTypeTag<id>, Properties::GridFluxVariablesCache>::LocalView;
86 template<std::size_t id> using GridVariables = GetPropType<SubDomainTypeTag<id>, Properties::GridVariables>;
87 template<std::size_t id> using Element = typename GridView<id>::template Codim<0>::Entity;
88 template<std::size_t id> using PrimaryVariables = typename MDTraits::template SubDomain<id>::PrimaryVariables;
89 template<std::size_t id> using SubControlVolumeFace = typename FVElementGeometry<id>::SubControlVolumeFace;
90
92
93 using VelocityVector = typename Element<stokesIdx>::Geometry::GlobalCoordinate;
94
96
97 struct StationaryStokesCouplingContext
98 {
99 Element<darcyIdx> element;
100 FVElementGeometry<darcyIdx> fvGeometry;
101 std::size_t darcyScvfIdx;
102 std::size_t stokesScvfIdx;
103 VolumeVariables<darcyIdx> volVars;
104 };
105
106 struct StationaryDarcyCouplingContext
107 {
108 Element<stokesIdx> element;
109 FVElementGeometry<stokesIdx> fvGeometry;
110 std::size_t stokesScvfIdx;
111 std::size_t darcyScvfIdx;
112 VelocityVector velocity;
113 VolumeVariables<stokesIdx> volVars;
114 };
115public:
116
120
122 StokesDarcyCouplingManager(std::shared_ptr<const GridGeometry<stokesIdx>> stokesFvGridGeometry,
123 std::shared_ptr<const GridGeometry<darcyIdx>> darcyFvGridGeometry)
124 { }
125
129 // \{
130
132 void init(std::shared_ptr<const Problem<stokesIdx>> stokesProblem,
133 std::shared_ptr<const Problem<darcyIdx>> darcyProblem,
134 const SolutionVector& curSol)
135 {
136 if (Dune::FloatCmp::ne(stokesProblem->gravity(), darcyProblem->spatialParams().gravity({})))
137 DUNE_THROW(Dune::InvalidStateException, "Both models must use the same gravity vector");
138
139 this->setSubProblems(std::make_tuple(stokesProblem, stokesProblem, darcyProblem));
140 this->curSol() = curSol;
141 couplingData_ = std::make_shared<CouplingData>(*this);
143 }
144
146 void update()
147 { }
148
149 // \}
150
152 void updateSolution(const SolutionVector& curSol)
153 { this->curSol() = curSol; }
154
157 {
158 couplingMapper_.computeCouplingMapsAndStencils(*this,
159 darcyToStokesCellCenterCouplingStencils_,
160 darcyToStokesFaceCouplingStencils_,
161 stokesCellCenterCouplingStencils_,
162 stokesFaceCouplingStencils_);
163
164 for(auto&& stencil : darcyToStokesCellCenterCouplingStencils_)
165 removeDuplicates_(stencil.second);
166 for(auto&& stencil : darcyToStokesFaceCouplingStencils_)
167 removeDuplicates_(stencil.second);
168 for(auto&& stencil : stokesCellCenterCouplingStencils_)
169 removeDuplicates_(stencil.second);
170 for(auto&& stencil : stokesFaceCouplingStencils_)
171 removeDuplicates_(stencil.second);
172 }
173
177 // \{
178
180
184 template<std::size_t i, class Assembler, std::enable_if_t<(i == stokesCellCenterIdx || i == stokesFaceIdx), int> = 0>
185 void bindCouplingContext(Dune::index_constant<i> domainI, const Element<stokesCellCenterIdx>& element, const Assembler& assembler) const
186 { bindCouplingContext(domainI, element); }
187
191 template<std::size_t i, std::enable_if_t<(i == stokesCellCenterIdx || i == stokesFaceIdx), int> = 0>
192 void bindCouplingContext(Dune::index_constant<i> domainI, const Element<stokesCellCenterIdx>& element) const
193 {
194 stokesCouplingContext_.clear();
195
196 const auto stokesElementIdx = this->problem(stokesIdx).gridGeometry().elementMapper().index(element);
197 boundStokesElemIdx_ = stokesElementIdx;
198
199 // do nothing if the element is not coupled to the other domain
200 if(!couplingMapper_.stokesElementToDarcyElementMap().count(stokesElementIdx))
201 return;
202
203 // prepare the coupling context
204 const auto& darcyIndices = couplingMapper_.stokesElementToDarcyElementMap().at(stokesElementIdx);
205 auto darcyFvGeometry = localView(this->problem(darcyIdx).gridGeometry());
206
207 for(auto&& indices : darcyIndices)
208 {
209 const auto& darcyElement = this->problem(darcyIdx).gridGeometry().boundingBoxTree().entitySet().entity(indices.eIdx);
210 darcyFvGeometry.bindElement(darcyElement);
211 const auto& scv = (*scvs(darcyFvGeometry).begin());
212
213 const auto darcyElemSol = elementSolution(darcyElement, this->curSol()[darcyIdx], this->problem(darcyIdx).gridGeometry());
214 VolumeVariables<darcyIdx> darcyVolVars;
215 darcyVolVars.update(darcyElemSol, this->problem(darcyIdx), darcyElement, scv);
216
217 // add the context
218 stokesCouplingContext_.push_back({darcyElement, darcyFvGeometry, indices.scvfIdx, indices.flipScvfIdx, darcyVolVars});
219 }
220 }
221
225 template<class Assembler>
226 void bindCouplingContext(Dune::index_constant<darcyIdx> domainI, const Element<darcyIdx>& element, const Assembler& assembler) const
227 { bindCouplingContext(domainI, element); }
228
232 void bindCouplingContext(Dune::index_constant<darcyIdx> domainI, const Element<darcyIdx>& element) const
233 {
234 darcyCouplingContext_.clear();
235
236 const auto darcyElementIdx = this->problem(darcyIdx).gridGeometry().elementMapper().index(element);
237 boundDarcyElemIdx_ = darcyElementIdx;
238
239 // do nothing if the element is not coupled to the other domain
240 if(!couplingMapper_.darcyElementToStokesElementMap().count(darcyElementIdx))
241 return;
242
243 // prepare the coupling context
244 const auto& stokesElementIndices = couplingMapper_.darcyElementToStokesElementMap().at(darcyElementIdx);
245 auto stokesFvGeometry = localView(this->problem(stokesIdx).gridGeometry());
246
247 for(auto&& indices : stokesElementIndices)
248 {
249 const auto& stokesElement = this->problem(stokesIdx).gridGeometry().boundingBoxTree().entitySet().entity(indices.eIdx);
250 stokesFvGeometry.bindElement(stokesElement);
251
252 VelocityVector faceVelocity(0.0);
253
254 for(const auto& scvf : scvfs(stokesFvGeometry))
255 {
256 if(scvf.index() == indices.scvfIdx)
257 faceVelocity[scvf.directionIndex()] = this->curSol()[stokesFaceIdx][scvf.dofIndex()];
258 }
259
260 using PriVarsType = typename VolumeVariables<stokesCellCenterIdx>::PrimaryVariables;
261 const auto& cellCenterPriVars = this->curSol()[stokesCellCenterIdx][indices.eIdx];
262 const auto elemSol = makeElementSolutionFromCellCenterPrivars<PriVarsType>(cellCenterPriVars);
263
264 VolumeVariables<stokesIdx> stokesVolVars;
265 for(const auto& scv : scvs(stokesFvGeometry))
266 stokesVolVars.update(elemSol, this->problem(stokesIdx), stokesElement, scv);
267
268 // add the context
269 darcyCouplingContext_.push_back({stokesElement, stokesFvGeometry, indices.scvfIdx, indices.flipScvfIdx, faceVelocity, stokesVolVars});
270 }
271 }
272
276 template<class LocalAssemblerI>
277 void updateCouplingContext(Dune::index_constant<darcyIdx> domainI,
278 const LocalAssemblerI& localAssemblerI,
279 Dune::index_constant<darcyIdx> domainJ,
280 std::size_t dofIdxGlobalJ,
281 const PrimaryVariables<darcyIdx>& priVarsJ,
282 int pvIdxJ)
283 {
284 this->curSol()[domainJ][dofIdxGlobalJ][pvIdxJ] = priVarsJ[pvIdxJ];
285 }
286
290 template<class LocalAssemblerI>
291 void updateCouplingContext(Dune::index_constant<darcyIdx> domainI,
292 const LocalAssemblerI& localAssemblerI,
293 Dune::index_constant<stokesCellCenterIdx> domainJ,
294 const std::size_t dofIdxGlobalJ,
295 const PrimaryVariables<stokesCellCenterIdx>& priVars,
296 int pvIdxJ)
297 {
298 this->curSol()[domainJ][dofIdxGlobalJ] = priVars;
299
300 for (auto& data : darcyCouplingContext_)
301 {
302 const auto stokesElemIdx = this->problem(stokesIdx).gridGeometry().elementMapper().index(data.element);
303
304 if(stokesElemIdx != dofIdxGlobalJ)
305 continue;
306
307 using PriVarsType = typename VolumeVariables<stokesCellCenterIdx>::PrimaryVariables;
308 const auto elemSol = makeElementSolutionFromCellCenterPrivars<PriVarsType>(priVars);
309
310 for(const auto& scv : scvs(data.fvGeometry))
311 data.volVars.update(elemSol, this->problem(stokesIdx), data.element, scv);
312 }
313 }
314
318 template<class LocalAssemblerI>
319 void updateCouplingContext(Dune::index_constant<darcyIdx> domainI,
320 const LocalAssemblerI& localAssemblerI,
321 Dune::index_constant<stokesFaceIdx> domainJ,
322 const std::size_t dofIdxGlobalJ,
323 const PrimaryVariables<stokesFaceIdx>& priVars,
324 int pvIdxJ)
325 {
326 this->curSol()[domainJ][dofIdxGlobalJ] = priVars;
327
328 for (auto& data : darcyCouplingContext_)
329 {
330 for(const auto& scvf : scvfs(data.fvGeometry))
331 {
332 if(scvf.dofIndex() == dofIdxGlobalJ)
333 data.velocity[scvf.directionIndex()] = priVars;
334 }
335 }
336 }
337
341 template<std::size_t i, class LocalAssemblerI, std::enable_if_t<(i == stokesCellCenterIdx || i == stokesFaceIdx), int> = 0>
342 void updateCouplingContext(Dune::index_constant<i> domainI,
343 const LocalAssemblerI& localAssemblerI,
344 Dune::index_constant<darcyIdx> domainJ,
345 const std::size_t dofIdxGlobalJ,
346 const PrimaryVariables<darcyIdx>& priVars,
347 int pvIdxJ)
348 {
349 this->curSol()[domainJ][dofIdxGlobalJ] = priVars;
350
351 for (auto& data : stokesCouplingContext_)
352 {
353 const auto darcyElemIdx = this->problem(darcyIdx).gridGeometry().elementMapper().index(data.element);
354
355 if(darcyElemIdx != dofIdxGlobalJ)
356 continue;
357
358 const auto darcyElemSol = elementSolution(data.element, this->curSol()[darcyIdx], this->problem(darcyIdx).gridGeometry());
359
360 for(const auto& scv : scvs(data.fvGeometry))
361 data.volVars.update(darcyElemSol, this->problem(darcyIdx), data.element, scv);
362 }
363 }
364
365 // \}
366
370 const auto& couplingData() const
371 {
372 return *couplingData_;
373 }
374
378 const auto& stokesCouplingContext(const Element<stokesIdx>& element, const SubControlVolumeFace<stokesIdx>& scvf) const
379 {
380 if (stokesCouplingContext_.empty() || boundStokesElemIdx_ != scvf.insideScvIdx())
382
383 for(const auto& context : stokesCouplingContext_)
384 {
385 if(scvf.index() == context.stokesScvfIdx)
386 return context;
387 }
388
389 DUNE_THROW(Dune::InvalidStateException, "No coupling context found at scvf " << scvf.center());
390 }
391
395 const auto& darcyCouplingContext(const Element<darcyIdx>& element, const SubControlVolumeFace<darcyIdx>& scvf) const
396 {
397 if (darcyCouplingContext_.empty() || boundDarcyElemIdx_ != scvf.insideScvIdx())
399
400 for(const auto& context : darcyCouplingContext_)
401 {
402 if(scvf.index() == context.darcyScvfIdx)
403 return context;
404 }
405
406 DUNE_THROW(Dune::InvalidStateException, "No coupling context found at scvf " << scvf.center());
407 }
408
412 // \{
413
417 const CouplingStencil& couplingStencil(Dune::index_constant<stokesCellCenterIdx> domainI,
418 const Element<stokesIdx>& element,
419 Dune::index_constant<darcyIdx> domainJ) const
420 {
421 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(element);
422 if(stokesCellCenterCouplingStencils_.count(eIdx))
423 return stokesCellCenterCouplingStencils_.at(eIdx);
424 else
425 return emptyStencil_;
426 }
427
432 template<std::size_t i, std::size_t j>
433 const CouplingStencil& couplingStencil(Dune::index_constant<i> domainI,
434 const Element<i>& element,
435 Dune::index_constant<j> domainJ) const
436 { return emptyStencil_; }
437
442 const CouplingStencil& couplingStencil(Dune::index_constant<darcyIdx> domainI,
443 const Element<darcyIdx>& element,
444 Dune::index_constant<stokesCellCenterIdx> domainJ) const
445 {
446 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(element);
447 if(darcyToStokesCellCenterCouplingStencils_.count(eIdx))
448 return darcyToStokesCellCenterCouplingStencils_.at(eIdx);
449 else
450 return emptyStencil_;
451 }
452
457 const CouplingStencil& couplingStencil(Dune::index_constant<darcyIdx> domainI,
458 const Element<darcyIdx>& element,
459 Dune::index_constant<stokesFaceIdx> domainJ) const
460 {
461 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(element);
462 if (darcyToStokesFaceCouplingStencils_.count(eIdx))
463 return darcyToStokesFaceCouplingStencils_.at(eIdx);
464 else
465 return emptyStencil_;
466 }
467
472 template<std::size_t i, std::size_t j>
473 const CouplingStencil& couplingStencil(Dune::index_constant<i> domainI,
474 const SubControlVolumeFace<stokesIdx>& scvf,
475 Dune::index_constant<j> domainJ) const
476 { return emptyStencil_; }
477
481 const CouplingStencil& couplingStencil(Dune::index_constant<stokesFaceIdx> domainI,
482 const SubControlVolumeFace<stokesIdx>& scvf,
483 Dune::index_constant<darcyIdx> domainJ) const
484 {
485 const auto faceDofIdx = scvf.dofIndex();
486 if(stokesFaceCouplingStencils_.count(faceDofIdx))
487 return stokesFaceCouplingStencils_.at(faceDofIdx);
488 else
489 return emptyStencil_;
490 }
491
492 // \}
493
497 template<class IdType>
498 const std::vector<std::size_t>& getAdditionalDofDependencies(IdType id, std::size_t stokesElementIdx) const
499 { return emptyStencil_; }
500
504 template<class IdType>
505 const std::vector<std::size_t>& getAdditionalDofDependenciesInverse(IdType id, std::size_t darcyElementIdx) const
506 { return emptyStencil_; }
507
511 bool isCoupledEntity(Dune::index_constant<stokesIdx>, const SubControlVolumeFace<stokesFaceIdx>& scvf) const
512 {
513 return stokesFaceCouplingStencils_.count(scvf.dofIndex());
514 }
515
519 bool isCoupledEntity(Dune::index_constant<darcyIdx>, const SubControlVolumeFace<darcyIdx>& scvf) const
520 {
521 return couplingMapper_.isCoupledDarcyScvf(scvf.index());
522 }
523
524protected:
525
527 std::vector<std::size_t>& emptyStencil()
528 { return emptyStencil_; }
529
530 void removeDuplicates_(std::vector<std::size_t>& stencil)
531 {
532 std::sort(stencil.begin(), stencil.end());
533 stencil.erase(std::unique(stencil.begin(), stencil.end()), stencil.end());
534 }
535
536private:
537
538 std::vector<bool> isCoupledDarcyDof_;
539 std::shared_ptr<CouplingData> couplingData_;
540
541 std::unordered_map<std::size_t, std::vector<std::size_t> > stokesCellCenterCouplingStencils_;
542 std::unordered_map<std::size_t, std::vector<std::size_t> > stokesFaceCouplingStencils_;
543 std::unordered_map<std::size_t, std::vector<std::size_t> > darcyToStokesCellCenterCouplingStencils_;
544 std::unordered_map<std::size_t, std::vector<std::size_t> > darcyToStokesFaceCouplingStencils_;
545 std::vector<std::size_t> emptyStencil_;
546
550 mutable std::vector<StationaryStokesCouplingContext> stokesCouplingContext_;
551 mutable std::vector<StationaryDarcyCouplingContext> darcyCouplingContext_;
552
553 mutable std::size_t boundStokesElemIdx_;
554 mutable std::size_t boundDarcyElemIdx_;
555
556 CouplingMapper couplingMapper_;
557};
558
559} // end namespace Dumux
560
561#endif
The interface of the coupling manager for multi domain problems.
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==DiscretizationMethod::box, BoxElementSolution< typename GridGeometry::LocalView, std::decay_t< decltype(std::declval< SolutionVector >()[0])> > >
Make an element solution for box schemes.
Definition: box/elementsolution.hh:115
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 vector of size number equations that can be used for Neumann fluxes, sources, residuals,...
Definition: common/properties.hh:51
Traits class encapsulating model specifications.
Definition: common/properties.hh:53
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:57
Stores the boundary types on an element.
Definition: common/properties.hh:98
Definition: common/properties.hh:101
The type for a global container for the volume variables.
Definition: common/properties.hh:108
The global vector of flux variable containers.
Definition: common/properties.hh:118
The grid variables object managing variable data on the grid (volvars/fluxvars cache)
Definition: common/properties.hh:122
Definition: couplingdata.hh:206
Coupling manager for Stokes and Darcy domains with equal dimension.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:49
void computeStencils()
Prepare the coupling stencils.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:156
void bindCouplingContext(Dune::index_constant< i > domainI, const Element< stokesCellCenterIdx > &element) const
prepares all data and variables that are necessary to evaluate the residual of an Darcy element (i....
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:192
const CouplingStencil & couplingStencil(Dune::index_constant< i > domainI, const Element< i > &element, Dune::index_constant< j > domainJ) const
The coupling stencil of domain I, i.e. which domain J DOFs the given domain I element's residual depe...
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:433
void updateCouplingContext(Dune::index_constant< darcyIdx > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< darcyIdx > domainJ, std::size_t dofIdxGlobalJ, const PrimaryVariables< darcyIdx > &priVarsJ, int pvIdxJ)
Update the coupling context for the Darcy residual w.r.t. Darcy DOFs.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:277
const CouplingStencil & couplingStencil(Dune::index_constant< stokesFaceIdx > domainI, const SubControlVolumeFace< stokesIdx > &scvf, Dune::index_constant< darcyIdx > domainJ) const
The coupling stencil of a Stokes face w.r.t. Darcy DOFs.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:481
std::vector< std::size_t > & emptyStencil()
Return a reference to an empty stencil.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:527
const CouplingStencil & couplingStencil(Dune::index_constant< stokesCellCenterIdx > domainI, const Element< stokesIdx > &element, Dune::index_constant< darcyIdx > domainJ) const
The coupling stencils.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:417
void updateCouplingContext(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< darcyIdx > domainJ, const std::size_t dofIdxGlobalJ, const PrimaryVariables< darcyIdx > &priVars, int pvIdxJ)
Update the coupling context for the Stokes cc residual w.r.t. the Darcy DOFs (FaceToDarcy)
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:342
const std::vector< std::size_t > & getAdditionalDofDependencies(IdType id, std::size_t stokesElementIdx) const
There are no additional dof dependencies.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:498
static constexpr auto darcyIdx
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:57
const auto & stokesCouplingContext(const Element< stokesIdx > &element, const SubControlVolumeFace< stokesIdx > &scvf) const
Access the coupling context needed for the Stokes domain.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:378
StokesDarcyCouplingManager(std::shared_ptr< const GridGeometry< stokesIdx > > stokesFvGridGeometry, std::shared_ptr< const GridGeometry< darcyIdx > > darcyFvGridGeometry)
Constructor.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:122
void update()
Update after the grid has changed.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:146
void bindCouplingContext(Dune::index_constant< darcyIdx > domainI, const Element< darcyIdx > &element) const
prepares all data and variables that are necessary to evaluate the residual of an Darcy element (i....
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:232
const auto & darcyCouplingContext(const Element< darcyIdx > &element, const SubControlVolumeFace< darcyIdx > &scvf) const
Access the coupling context needed for the Darcy domain.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:395
void removeDuplicates_(std::vector< std::size_t > &stencil)
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:530
void init(std::shared_ptr< const Problem< stokesIdx > > stokesProblem, std::shared_ptr< const Problem< darcyIdx > > darcyProblem, const SolutionVector &curSol)
Methods to be accessed by main.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:132
const CouplingStencil & couplingStencil(Dune::index_constant< darcyIdx > domainI, const Element< darcyIdx > &element, Dune::index_constant< stokesCellCenterIdx > domainJ) const
The coupling stencil of domain I, i.e. which domain J dofs the given domain I element's residual depe...
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:442
const auto & couplingData() const
Access the coupling data.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:370
void updateCouplingContext(Dune::index_constant< darcyIdx > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< stokesCellCenterIdx > domainJ, const std::size_t dofIdxGlobalJ, const PrimaryVariables< stokesCellCenterIdx > &priVars, int pvIdxJ)
Update the coupling context for the Darcy residual w.r.t. the Stokes cell-center DOFs (DarcyToCC)
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:291
const CouplingStencil & couplingStencil(Dune::index_constant< darcyIdx > domainI, const Element< darcyIdx > &element, Dune::index_constant< stokesFaceIdx > domainJ) const
The coupling stencil of domain I, i.e. which domain J dofs the given domain I element's residual depe...
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:457
bool isCoupledEntity(Dune::index_constant< stokesIdx >, const SubControlVolumeFace< stokesFaceIdx > &scvf) const
Returns whether a given free flow scvf is coupled to the other domain.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:511
void bindCouplingContext(Dune::index_constant< darcyIdx > domainI, const Element< darcyIdx > &element, const Assembler &assembler) const
prepares all data and variables that are necessary to evaluate the residual of an Darcy element (i....
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:226
static constexpr auto stokesCellCenterIdx
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:55
const std::vector< std::size_t > & getAdditionalDofDependenciesInverse(IdType id, std::size_t darcyElementIdx) const
There are no additional dof dependencies.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:505
void bindCouplingContext(Dune::index_constant< i > domainI, const Element< stokesCellCenterIdx > &element, const Assembler &assembler) const
prepares all data and variables that are necessary to evaluate the residual of an Darcy element (i....
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:185
void updateCouplingContext(Dune::index_constant< darcyIdx > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< stokesFaceIdx > domainJ, const std::size_t dofIdxGlobalJ, const PrimaryVariables< stokesFaceIdx > &priVars, int pvIdxJ)
Update the coupling context for the Darcy residual w.r.t. the Stokes face DOFs (DarcyToFace)
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:319
const CouplingStencil & couplingStencil(Dune::index_constant< i > domainI, const SubControlVolumeFace< stokesIdx > &scvf, Dune::index_constant< j > domainJ) const
The coupling stencil of domain I, i.e. which domain J DOFs the given domain I element's residual depe...
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:473
bool isCoupledEntity(Dune::index_constant< darcyIdx >, const SubControlVolumeFace< darcyIdx > &scvf) const
Returns whether a given free flow scvf is coupled to the other domain.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:519
void updateSolution(const SolutionVector &curSol)
Update the solution vector before assembly.
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:152
static constexpr auto stokesFaceIdx
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:54
static constexpr auto stokesIdx
Definition: multidomain/boundary/stokesdarcy/couplingmanager.hh:56
Coupling mapper for Stokes and Darcy domains with equal dimension.
Definition: boundary/stokesdarcy/couplingmapper.hh:42
bool isCoupledDarcyScvf(std::size_t darcyScvfIdx) const
Returns whether a Darcy scvf is coupled to the other domain.
Definition: boundary/stokesdarcy/couplingmapper.hh:133
void computeCouplingMapsAndStencils(const CouplingManager &couplingManager, Stencils &darcyToStokesCellCenterStencils, Stencils &darcyToStokesFaceStencils, Stencils &stokesCellCenterToDarcyStencils, Stencils &stokesFaceToDarcyStencils)
Main update routine.
Definition: boundary/stokesdarcy/couplingmapper.hh:56
const auto & darcyElementToStokesElementMap() const
A map that returns all Stokes elements coupled to a Darcy element.
Definition: boundary/stokesdarcy/couplingmapper.hh:141
const auto & stokesElementToDarcyElementMap() const
A map that returns all Darcy elements coupled to a Stokes element.
Definition: boundary/stokesdarcy/couplingmapper.hh:149
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
SolutionVector & curSol()
the solution vector of the coupled problem
Definition: multidomain/couplingmanager.hh:278
Base coupling manager for the staggered discretization.
Definition: staggeredcouplingmanager.hh:42
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
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
Declares all properties used in Dumux.
The local element solution class for staggered methods.