3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
couplingmanagerbase.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 *****************************************************************************/
26#ifndef DUMUX_MULTIDOMAIN_EMBEDDED_COUPLINGMANAGERBASE_HH
27#define DUMUX_MULTIDOMAIN_EMBEDDED_COUPLINGMANAGERBASE_HH
28
29#include <iostream>
30#include <fstream>
31#include <string>
32#include <utility>
33#include <unordered_map>
34
35#include <dune/common/timer.hh>
36#include <dune/geometry/quadraturerules.hh>
37
46
47namespace Dumux {
48
50template<class MDTraits>
52{
53private:
54 template<std::size_t i> using SubDomainTypeTag = typename MDTraits::template SubDomain<i>::TypeTag;
55 template<std::size_t i> using GridGeometry = GetPropType<SubDomainTypeTag<i>, Properties::GridGeometry>;
56 template<std::size_t i> using NumEqVector = GetPropType<SubDomainTypeTag<i>, Properties::NumEqVector>;
57public:
59 template<std::size_t i>
61
63 template<std::size_t i>
65
68};
69
75template<class MDTraits, class Implementation, class PSTraits = DefaultPointSourceTraits<MDTraits>>
77: public CouplingManager<MDTraits>
78{
80 using Scalar = typename MDTraits::Scalar;
81 static constexpr auto bulkIdx = typename MDTraits::template SubDomain<0>::Index();
82 static constexpr auto lowDimIdx = typename MDTraits::template SubDomain<1>::Index();
83 using SolutionVector = typename MDTraits::SolutionVector;
84 using PointSourceData = typename PSTraits::PointSourceData;
85
86 // the sub domain type tags
87 template<std::size_t id> using PointSource = typename PSTraits::template PointSource<id>;
88 template<std::size_t id> using SubDomainTypeTag = typename MDTraits::template SubDomain<id>::TypeTag;
89 template<std::size_t id> using Problem = GetPropType<SubDomainTypeTag<id>, Properties::Problem>;
90 template<std::size_t id> using PrimaryVariables = GetPropType<SubDomainTypeTag<id>, Properties::PrimaryVariables>;
91 template<std::size_t id> using GridGeometry = GetPropType<SubDomainTypeTag<id>, Properties::GridGeometry>;
92 template<std::size_t id> using GridView = typename GridGeometry<id>::GridView;
93 template<std::size_t id> using ElementMapper = typename GridGeometry<id>::ElementMapper;
94 template<std::size_t id> using Element = typename GridView<id>::template Codim<0>::Entity;
95 template<std::size_t id> using GridIndex = typename IndexTraits<GridView<id>>::GridIndex;
96 template<std::size_t id> using CouplingStencil = std::vector<GridIndex<id>>;
97
98 static constexpr int bulkDim = GridView<bulkIdx>::dimension;
99 static constexpr int lowDimDim = GridView<lowDimIdx>::dimension;
100 static constexpr int dimWorld = GridView<bulkIdx>::dimensionworld;
101
102 template<std::size_t id>
103 static constexpr bool isBox()
104 { return GridGeometry<id>::discMethod == DiscretizationMethod::box; }
105
106 using GlobalPosition = typename Element<bulkIdx>::Geometry::GlobalCoordinate;
107 using GlueType = MultiDomainGlue<GridView<bulkIdx>, GridView<lowDimIdx>, ElementMapper<bulkIdx>, ElementMapper<lowDimIdx>>;
108
109public:
111 using MultiDomainTraits = MDTraits;
113 using PointSourceTraits = PSTraits;
115 template<std::size_t id> using CouplingStencils = std::unordered_map<GridIndex<id>, CouplingStencil<id>>;
116
120 void updateAfterGridAdaption(std::shared_ptr<const GridGeometry<bulkIdx>> bulkGridGeometry,
121 std::shared_ptr<const GridGeometry<lowDimIdx>> lowDimGridGeometry)
122 {
123 glue_ = std::make_shared<GlueType>();
124 }
125
129 EmbeddedCouplingManagerBase(std::shared_ptr<const GridGeometry<bulkIdx>> bulkGridGeometry,
130 std::shared_ptr<const GridGeometry<lowDimIdx>> lowDimGridGeometry)
131 {
132 updateAfterGridAdaption(bulkGridGeometry, lowDimGridGeometry);
133 }
134
138 // \{
139
140 void init(std::shared_ptr<Problem<bulkIdx>> bulkProblem,
141 std::shared_ptr<Problem<lowDimIdx>> lowDimProblem,
142 const SolutionVector& curSol)
143 {
144 this->updateSolution(curSol);
145 this->setSubProblems(std::make_tuple(bulkProblem, lowDimProblem));
146
147 integrationOrder_ = getParam<int>("MixedDimension.IntegrationOrder", 1);
148 asImp_().computePointSourceData(integrationOrder_);
149 }
150
151 // \}
152
156 // \{
157
172 template<std::size_t i, std::size_t j>
173 const CouplingStencil<j>& couplingStencil(Dune::index_constant<i> domainI,
174 const Element<i>& element,
175 Dune::index_constant<j> domainJ) const
176 {
177 static_assert(i != j, "A domain cannot be coupled to itself!");
178
179 const auto eIdx = this->problem(domainI).gridGeometry().elementMapper().index(element);
180 if (couplingStencils(domainI).count(eIdx))
181 return couplingStencils(domainI).at(eIdx);
182 else
183 return emptyStencil(domainI);
184 }
185
198 template<std::size_t i, std::size_t j, class LocalAssemblerI>
199 decltype(auto) evalCouplingResidual(Dune::index_constant<i> domainI,
200 const LocalAssemblerI& localAssemblerI,
201 Dune::index_constant<j> domainJ,
202 std::size_t dofIdxGlobalJ)
203 {
204 static_assert(i != j, "A domain cannot be coupled to itself!");
205
206 typename LocalAssemblerI::LocalResidual::ElementResidualVector residual;
207
208 const auto& element = localAssemblerI.element();
209 const auto& fvGeometry = localAssemblerI.fvGeometry();
210 const auto& curElemVolVars = localAssemblerI.curElemVolVars();
211
212 residual.resize(fvGeometry.numScv());
213 for (const auto& scv : scvs(fvGeometry))
214 {
215 auto couplingSource = this->problem(domainI).scvPointSources(element, fvGeometry, curElemVolVars, scv);
216 couplingSource += this->problem(domainI).source(element, fvGeometry, curElemVolVars, scv);
217 couplingSource *= -GridGeometry<i>::Extrusion::volume(scv)*curElemVolVars[scv].extrusionFactor();
218 residual[scv.indexInElement()] = couplingSource;
219 }
220 return residual;
221 }
222
223 // \}
224
225 /* \brief Compute integration point point sources and associated data
226 *
227 * This method uses grid glue to intersect the given grids. Over each intersection
228 * we later need to integrate a source term. This method places point sources
229 * at each quadrature point and provides the point source with the necessary
230 * information to compute integrals (quadrature weight and integration element)
231 * \param order The order of the quadrature rule for integration of sources over an intersection
232 * \param verbose If the point source computation is verbose
233 */
234 void computePointSourceData(std::size_t order = 1, bool verbose = false)
235 {
236 // initilize the maps
237 // do some logging and profiling
238 Dune::Timer watch;
239 std::cout << "Initializing the point sources..." << std::endl;
240
241 // clear all internal members like pointsource vectors and stencils
242 // initializes the point source id counter
243 clear();
244
245 // precompute the vertex indices for efficiency for the box method
246 this->precomputeVertexIndices(bulkIdx);
247 this->precomputeVertexIndices(lowDimIdx);
248
249 const auto& bulkGridGeometry = this->problem(bulkIdx).gridGeometry();
250 const auto& lowDimGridGeometry = this->problem(lowDimIdx).gridGeometry();
251
252 // intersect the bounding box trees
253 glueGrids();
254
255 pointSourceData_.reserve(this->glue().size());
256 averageDistanceToBulkCell_.reserve(this->glue().size());
257 for (const auto& is : intersections(this->glue()))
258 {
259 // all inside elements are identical...
260 const auto& inside = is.targetEntity(0);
261 // get the intersection geometry for integrating over it
262 const auto intersectionGeometry = is.geometry();
263
264 // get the Gaussian quadrature rule for the local intersection
265 const auto& quad = Dune::QuadratureRules<Scalar, lowDimDim>::rule(intersectionGeometry.type(), order);
266 const std::size_t lowDimElementIdx = lowDimGridGeometry.elementMapper().index(inside);
267
268 // iterate over all quadrature points
269 for (auto&& qp : quad)
270 {
271 // compute the coupling stencils
272 for (std::size_t outsideIdx = 0; outsideIdx < is.numDomainNeighbors(); ++outsideIdx)
273 {
274 const auto& outside = is.domainEntity(outsideIdx);
275 const std::size_t bulkElementIdx = bulkGridGeometry.elementMapper().index(outside);
276
277 // each quadrature point will be a point source for the sub problem
278 const auto globalPos = intersectionGeometry.global(qp.position());
279 const auto id = idCounter_++;
280 const auto qpweight = qp.weight();
281 const auto ie = intersectionGeometry.integrationElement(qp.position());
282 pointSources(bulkIdx).emplace_back(globalPos, id, qpweight, ie, std::vector<std::size_t>({bulkElementIdx}));
283 pointSources(bulkIdx).back().setEmbeddings(is.numDomainNeighbors());
284 pointSources(lowDimIdx).emplace_back(globalPos, id, qpweight, ie, std::vector<std::size_t>({lowDimElementIdx}));
285 pointSources(lowDimIdx).back().setEmbeddings(is.numDomainNeighbors());
286
287 // pre compute additional data used for the evaluation of
288 // the actual solution dependent source term
289 PointSourceData psData;
290
291 if constexpr (isBox<lowDimIdx>())
292 {
293 using ShapeValues = std::vector<Dune::FieldVector<Scalar, 1> >;
294 const auto lowDimGeometry = this->problem(lowDimIdx).gridGeometry().element(lowDimElementIdx).geometry();
295 ShapeValues shapeValues;
296 this->getShapeValues(lowDimIdx, this->problem(lowDimIdx).gridGeometry(), lowDimGeometry, globalPos, shapeValues);
297 psData.addLowDimInterpolation(shapeValues, this->vertexIndices(lowDimIdx, lowDimElementIdx), lowDimElementIdx);
298 }
299 else
300 {
301 psData.addLowDimInterpolation(lowDimElementIdx);
302 }
303
304 // add data needed to compute integral over the circle
305 if constexpr (isBox<bulkIdx>())
306 {
307 using ShapeValues = std::vector<Dune::FieldVector<Scalar, 1> >;
308 const auto bulkGeometry = this->problem(bulkIdx).gridGeometry().element(bulkElementIdx).geometry();
309 ShapeValues shapeValues;
310 this->getShapeValues(bulkIdx, this->problem(bulkIdx).gridGeometry(), bulkGeometry, globalPos, shapeValues);
311 psData.addBulkInterpolation(shapeValues, this->vertexIndices(bulkIdx, bulkElementIdx), bulkElementIdx);
312 }
313 else
314 {
315 psData.addBulkInterpolation(bulkElementIdx);
316 }
317
318 // publish point source data in the global vector
319 this->pointSourceData().emplace_back(std::move(psData));
320
321 // compute average distance to bulk cell
322 averageDistanceToBulkCell_.push_back(averageDistancePointGeometry(globalPos, outside.geometry()));
323
324 // export the lowdim coupling stencil
325 // we insert all vertices / elements and make it unique later
326 if (isBox<bulkIdx>())
327 {
328 const auto& vertices = this->vertexIndices(bulkIdx, bulkElementIdx);
329 this->couplingStencils(lowDimIdx)[lowDimElementIdx].insert(this->couplingStencils(lowDimIdx)[lowDimElementIdx].end(),
330 vertices.begin(), vertices.end());
331 }
332 else
333 {
334 this->couplingStencils(lowDimIdx)[lowDimElementIdx].push_back(bulkElementIdx);
335 }
336
337 // export the bulk coupling stencil
338 // we insert all vertices / elements and make it unique later
339 if (isBox<lowDimIdx>())
340 {
341 const auto& vertices = this->vertexIndices(lowDimIdx, lowDimElementIdx);
342 this->couplingStencils(bulkIdx)[bulkElementIdx].insert(this->couplingStencils(bulkIdx)[bulkElementIdx].end(),
343 vertices.begin(), vertices.end());
344
345 }
346 else
347 {
348 this->couplingStencils(bulkIdx)[bulkElementIdx].push_back(lowDimElementIdx);
349 }
350 }
351 }
352 }
353
354 // make stencils unique
355 using namespace Dune::Hybrid;
356 forEach(integralRange(std::integral_constant<std::size_t, MDTraits::numSubDomains>{}), [&](const auto domainIdx)
357 {
358 for (auto&& stencil : this->couplingStencils(domainIdx))
359 {
360 std::sort(stencil.second.begin(), stencil.second.end());
361 stencil.second.erase(std::unique(stencil.second.begin(), stencil.second.end()), stencil.second.end());
362 }
363 });
364
365 std::cout << "took " << watch.elapsed() << " seconds." << std::endl;
366 }
367
371 // \{
372
374 const PointSourceData& pointSourceData(std::size_t id) const
375 { return pointSourceData_[id]; }
376
378 template<std::size_t id>
379 const GridView<id>& gridView(Dune::index_constant<id> domainIdx) const
380 { return this->problem(domainIdx).gridGeometry().gridView(); }
381
383 PrimaryVariables<bulkIdx> bulkPriVars(std::size_t id) const
384 { return pointSourceData_[id].interpolateBulk(this->curSol()[bulkIdx]); }
385
387 PrimaryVariables<lowDimIdx> lowDimPriVars(std::size_t id) const
388 { return pointSourceData_[id].interpolateLowDim(this->curSol()[lowDimIdx]); }
389
391 Scalar averageDistance(std::size_t id) const
392 { return averageDistanceToBulkCell_[id]; }
393
395 const std::vector<PointSource<bulkIdx>>& bulkPointSources() const
396 { return std::get<bulkIdx>(pointSources_); }
397
399 const std::vector<PointSource<lowDimIdx>>& lowDimPointSources() const
400 { return std::get<lowDimIdx>(pointSources_); }
401
403 template<std::size_t i>
404 const std::vector<PointSource<i>>& pointSources(Dune::index_constant<i> dom) const
405 { return std::get<i>(pointSources_); }
406
408 template<std::size_t i>
409 const CouplingStencils<i>& couplingStencils(Dune::index_constant<i> dom) const
410 { return std::get<i>(couplingStencils_); }
411
413 const std::vector<PointSourceData>& pointSourceData() const
414 { return pointSourceData_; }
415
417 template<std::size_t i>
418 const CouplingStencil<i>& emptyStencil(Dune::index_constant<i> dom) const
419 { return std::get<i>(emptyStencil_); }
420
421protected:
422
424 template<std::size_t id>
425 void precomputeVertexIndices(Dune::index_constant<id> domainIdx)
426 {
427 // fill helper structure for box discretization
428 if constexpr (isBox<domainIdx>())
429 {
430 this->vertexIndices(domainIdx).resize(gridView(domainIdx).size(0));
431 for (const auto& element : elements(gridView(domainIdx)))
432 {
433 constexpr int dim = GridView<domainIdx>::dimension;
434 const auto eIdx = this->problem(domainIdx).gridGeometry().elementMapper().index(element);
435 this->vertexIndices(domainIdx, eIdx).resize(element.subEntities(dim));
436 for (int i = 0; i < element.subEntities(dim); ++i)
437 this->vertexIndices(domainIdx, eIdx)[i] = this->problem(domainIdx).gridGeometry().vertexMapper().subIndex(element, i, dim);
438 }
439 }
440 }
441
443 template<std::size_t i, class FVGG, class Geometry, class ShapeValues>
444 void getShapeValues(Dune::index_constant<i> domainI, const FVGG& gridGeometry, const Geometry& geo, const GlobalPosition& globalPos, ShapeValues& shapeValues)
445 {
446 if constexpr (FVGG::discMethod == DiscretizationMethod::box)
447 {
448 const auto ipLocal = geo.local(globalPos);
449 const auto& localBasis = this->problem(domainI).gridGeometry().feCache().get(geo.type()).localBasis();
450 localBasis.evaluateFunction(ipLocal, shapeValues);
451 }
452 else
453 DUNE_THROW(Dune::InvalidStateException, "Shape values requested for other discretization than box!");
454 }
455
457 void clear()
458 {
459 pointSources(bulkIdx).clear();
460 pointSources(lowDimIdx).clear();
461 couplingStencils(bulkIdx).clear();
462 couplingStencils(lowDimIdx).clear();
463 vertexIndices(bulkIdx).clear();
464 vertexIndices(lowDimIdx).clear();
465 pointSourceData_.clear();
466 averageDistanceToBulkCell_.clear();
467
468 idCounter_ = 0;
469 }
470
473 {
474 const auto& bulkGridGeometry = this->problem(bulkIdx).gridGeometry();
475 const auto& lowDimGridGeometry = this->problem(lowDimIdx).gridGeometry();
476
477 // intersect the bounding box trees
478 glue_->build(bulkGridGeometry.boundingBoxTree(), lowDimGridGeometry.boundingBoxTree());
479 }
480
482 std::vector<PointSourceData>& pointSourceData()
483 { return pointSourceData_; }
484
486 std::vector<Scalar>& averageDistanceToBulkCell()
487 { return averageDistanceToBulkCell_; }
488
490 template<std::size_t i>
491 std::vector<PointSource<i>>& pointSources(Dune::index_constant<i> dom)
492 { return std::get<i>(pointSources_); }
493
495 template<std::size_t i>
496 CouplingStencils<i>& couplingStencils(Dune::index_constant<i> dom)
497 { return std::get<i>(couplingStencils_); }
498
500 template<std::size_t i>
501 std::vector<GridIndex<i>>& vertexIndices(Dune::index_constant<i> dom, GridIndex<i> eIdx)
502 { return std::get<i>(vertexIndices_)[eIdx]; }
503
505 template<std::size_t i>
506 std::vector<std::vector<GridIndex<i>>>& vertexIndices(Dune::index_constant<i> dom)
507 { return std::get<i>(vertexIndices_); }
508
509 const GlueType& glue() const
510 { return *glue_; }
511
513 Implementation &asImp_()
514 { return *static_cast<Implementation *>(this); }
515
517 const Implementation &asImp_() const
518 { return *static_cast<const Implementation *>(this); }
519
521 std::size_t idCounter_ = 0;
522
523private:
524
526 std::tuple<std::vector<PointSource<bulkIdx>>, std::vector<PointSource<lowDimIdx>>> pointSources_;
527 std::vector<PointSourceData> pointSourceData_;
528 std::vector<Scalar> averageDistanceToBulkCell_;
529
531 std::tuple<std::vector<std::vector<GridIndex<bulkIdx>>>,
532 std::vector<std::vector<GridIndex<lowDimIdx>>>> vertexIndices_;
533 std::tuple<CouplingStencils<bulkIdx>, CouplingStencils<lowDimIdx>> couplingStencils_;
534 std::tuple<CouplingStencil<bulkIdx>, CouplingStencil<lowDimIdx>> emptyStencil_;
535
537 std::shared_ptr<GlueType> glue_;
538
540 int integrationOrder_ = 1;
541};
542
543} // end namespace Dumux
544
545#endif
The available discretization methods in Dumux.
An integration point source class, i.e. sources located at a single point in space associated with a ...
Data associated with a point source.
A class glueing two grids of potentially different dimension geometrically. Intersections are compute...
Geometry::ctype averageDistancePointGeometry(const typename Geometry::GlobalCoordinate &p, const Geometry &geometry, std::size_t integrationOrder=2)
Compute the average distance from a point to a geometry by integration.
Definition: geometry/distance.hh:36
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
Struture to define the index types used for grid and local indices.
Definition: indextraits.hh:38
A vector of primary variables.
Definition: common/properties.hh:49
A vector of size number equations that can be used for Neumann fluxes, sources, residuals,...
Definition: common/properties.hh:51
Property to specify the type of a problem which has to be solved.
Definition: common/properties.hh:57
Definition: common/properties.hh:101
A class representing the intersection entites two geometric entity sets.
Definition: geometry/intersectionentityset.hh:55
Definition: multidomain/couplingmanager.hh:46
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
void updateSolution(const SolutionVector &curSol)
Updates the entire solution vector, e.g. before assembly or after grid adaption.
Definition: multidomain/couplingmanager.hh:186
the default point source traits
Definition: couplingmanagerbase.hh:52
Manages the coupling between bulk elements and lower dimensional elements Point sources on each integ...
Definition: couplingmanagerbase.hh:78
std::vector< GridIndex< i > > & vertexIndices(Dune::index_constant< i > dom, GridIndex< i > eIdx)
Return a reference to the vertex indices.
Definition: couplingmanagerbase.hh:501
Scalar averageDistance(std::size_t id) const
return the average distance to the coupled bulk cell center
Definition: couplingmanagerbase.hh:391
Implementation & asImp_()
Returns the implementation of the problem (i.e. static polymorphism)
Definition: couplingmanagerbase.hh:513
const std::vector< PointSource< i > > & pointSources(Dune::index_constant< i > dom) const
Return the point source if domain i.
Definition: couplingmanagerbase.hh:404
PSTraits PointSourceTraits
export the point source traits
Definition: couplingmanagerbase.hh:113
EmbeddedCouplingManagerBase(std::shared_ptr< const GridGeometry< bulkIdx > > bulkGridGeometry, std::shared_ptr< const GridGeometry< lowDimIdx > > lowDimGridGeometry)
Constructor.
Definition: couplingmanagerbase.hh:129
const CouplingStencils< i > & couplingStencils(Dune::index_constant< i > dom) const
Return reference to bulk coupling stencil member of domain i.
Definition: couplingmanagerbase.hh:409
void updateAfterGridAdaption(std::shared_ptr< const GridGeometry< bulkIdx > > bulkGridGeometry, std::shared_ptr< const GridGeometry< lowDimIdx > > lowDimGridGeometry)
call this after grid adaption
Definition: couplingmanagerbase.hh:120
const GridView< id > & gridView(Dune::index_constant< id > domainIdx) const
Return a reference to the bulk problem.
Definition: couplingmanagerbase.hh:379
const Implementation & asImp_() const
Returns the implementation of the problem (i.e. static polymorphism)
Definition: couplingmanagerbase.hh:517
std::vector< Scalar > & averageDistanceToBulkCell()
Return reference to average distances to bulk cell.
Definition: couplingmanagerbase.hh:486
const CouplingStencil< i > & emptyStencil(Dune::index_constant< i > dom) const
Return a reference to an empty stencil.
Definition: couplingmanagerbase.hh:418
const std::vector< PointSource< bulkIdx > > & bulkPointSources() const
Return reference to bulk point sources.
Definition: couplingmanagerbase.hh:395
const PointSourceData & pointSourceData(std::size_t id) const
Methods to be accessed by the subproblems.
Definition: couplingmanagerbase.hh:374
std::vector< std::vector< GridIndex< i > > > & vertexIndices(Dune::index_constant< i > dom)
Return a reference to the vertex indices container.
Definition: couplingmanagerbase.hh:506
void glueGrids()
compute the intersections between the two grids
Definition: couplingmanagerbase.hh:472
MDTraits MultiDomainTraits
export traits
Definition: couplingmanagerbase.hh:111
void computePointSourceData(std::size_t order=1, bool verbose=false)
Definition: couplingmanagerbase.hh:234
std::vector< PointSourceData > & pointSourceData()
Return reference to point source data vector member.
Definition: couplingmanagerbase.hh:482
const CouplingStencil< j > & couplingStencil(Dune::index_constant< i > domainI, const Element< i > &element, Dune::index_constant< j > domainJ) const
Methods to be accessed by the assembly.
Definition: couplingmanagerbase.hh:173
PrimaryVariables< lowDimIdx > lowDimPriVars(std::size_t id) const
Return data for a low dim point source with the identifier id.
Definition: couplingmanagerbase.hh:387
void getShapeValues(Dune::index_constant< i > domainI, const FVGG &gridGeometry, const Geometry &geo, const GlobalPosition &globalPos, ShapeValues &shapeValues)
compute the shape function for a given point and geometry
Definition: couplingmanagerbase.hh:444
const std::vector< PointSource< lowDimIdx > > & lowDimPointSources() const
Return reference to low dim point sources.
Definition: couplingmanagerbase.hh:399
void clear()
Clear all internal data members.
Definition: couplingmanagerbase.hh:457
const GlueType & glue() const
Definition: couplingmanagerbase.hh:509
decltype(auto) evalCouplingResidual(Dune::index_constant< i > domainI, const LocalAssemblerI &localAssemblerI, Dune::index_constant< j > domainJ, std::size_t dofIdxGlobalJ)
evaluates the element residual of a coupled element of domain i which depends on the variables at the...
Definition: couplingmanagerbase.hh:199
void init(std::shared_ptr< Problem< bulkIdx > > bulkProblem, std::shared_ptr< Problem< lowDimIdx > > lowDimProblem, const SolutionVector &curSol)
Methods to be accessed by main.
Definition: couplingmanagerbase.hh:140
CouplingStencils< i > & couplingStencils(Dune::index_constant< i > dom)
Return reference to bulk coupling stencil member of domain i.
Definition: couplingmanagerbase.hh:496
std::vector< PointSource< i > > & pointSources(Dune::index_constant< i > dom)
Return the point source if domain i.
Definition: couplingmanagerbase.hh:491
std::unordered_map< GridIndex< id >, CouplingStencil< id > > CouplingStencils
export stencil types
Definition: couplingmanagerbase.hh:115
PrimaryVariables< bulkIdx > bulkPriVars(std::size_t id) const
Return data for a bulk point source with the identifier id.
Definition: couplingmanagerbase.hh:383
void precomputeVertexIndices(Dune::index_constant< id > domainIdx)
computes the vertex indices per element for the box method
Definition: couplingmanagerbase.hh:425
std::size_t idCounter_
id generator for point sources
Definition: couplingmanagerbase.hh:521
const std::vector< PointSourceData > & pointSourceData() const
Return reference to point source data vector member.
Definition: couplingmanagerbase.hh:413
An integration point source class with an identifier to attach data and a quadrature weight and integ...
Definition: integrationpointsource.hh:45
A helper class calculating a DOF-index to point source map.
Definition: integrationpointsource.hh:135
A point source data class used for integration in multidimension models.
Definition: pointsourcedata.hh:43
Helper functions for distance queries.
Algorithms that finds which geometric entites intersect.
Declares all properties used in Dumux.
The interface of the coupling manager for multi domain problems.