3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
facet/cellcentered/tpfa/couplingmapper.hh
Go to the documentation of this file.
1/*****************************************************************************
2 * See the file COPYING for full copying permissions. *
3 * *
4 * This program is free software: you can redistribute it and/or modify *
5 * it under the terms of the GNU General Public License as published by *
6 * the Free Software Foundation, either version 3 of the License, or *
7 * (at your option) any later version. *
8 * *
9 * This program is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
12 * GNU General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU General Public License *
15 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
16 *****************************************************************************/
22#ifndef DUMUX_CCTPFA_FACETCOUPLING_MAPPER_HH
23#define DUMUX_CCTPFA_FACETCOUPLING_MAPPER_HH
24
25#include <dune/common/indices.hh>
26#include <dune/common/float_cmp.hh>
27
32
33namespace Dumux {
34
47template<class BulkFVG, class LowDimFVG, std::size_t bulkId, std::size_t lowDimId>
48class FacetCouplingMapper<BulkFVG, LowDimFVG, bulkId, lowDimId, DiscretizationMethod::cctpfa>
49: public virtual FacetCouplingMapperBase<BulkFVG, LowDimFVG, bulkId, lowDimId>
50{
52 using LowDimElement = typename LowDimFVG::GridView::template Codim<0>::Entity;
53 using GlobalPosition = typename LowDimElement::Geometry::GlobalCoordinate;
54
55public:
57 using ParentType::bulkGridId;
58 using ParentType::facetGridId;
59
68 template< class Embeddings >
69 void update(const BulkFVG& bulkFvGridGeometry,
70 const LowDimFVG& lowDimFvGridGeometry,
71 std::shared_ptr<const Embeddings> embeddings)
72 {
73 // define the policy how to add map entries for given lowdim element and adjoined entity indices
74 auto addCouplingEntryPolicy = [&] (auto&& adjoinedEntityIndices,
75 const LowDimElement& lowDimElement,
76 const LowDimFVG& lowDimFvGridGeometry,
77 const BulkFVG& bulkFvGridGeometry)
78 {
79 using LowDimIndexType = typename IndexTraits<typename LowDimFVG::GridView>::GridIndex;
80 using BulkIndexType = typename IndexTraits<typename BulkFVG::GridView>::GridIndex;
81
82 const auto lowDimGeometry = lowDimElement.geometry();
83 const auto lowDimElemIdx = lowDimFvGridGeometry.elementMapper().index(lowDimElement);
84 auto& lowDimData = this->couplingMap_(facetGridId, bulkGridId)[lowDimElemIdx];
85
86 // find the scvfs in the adjoined entities coinciding with the low dim element
87 // since the bulk domain uses tpfa, there is always only going to be one scvf
88 for (auto bulkElemIdx : adjoinedEntityIndices)
89 {
90 const auto bulkElement = bulkFvGridGeometry.element(bulkElemIdx);
91
92 auto fvGeometry = localView(bulkFvGridGeometry);
93 fvGeometry.bindElement(bulkElement);
94
95 bool found = false;
96 BulkIndexType embeddedScvfIdx;
97 for (const auto& scvf : scvfs(fvGeometry))
98 {
99 // for non-boundary faces, it suffices to check if one
100 // of the outside scv indices is in the set of embedments
101 if (!scvf.boundary())
102 {
103 if ( std::find(adjoinedEntityIndices.begin(),
104 adjoinedEntityIndices.end(),
105 scvf.outsideScvIdx()) != adjoinedEntityIndices.end() )
106 {
107 embeddedScvfIdx = scvf.index();
108 found = true; break;
109 }
110 }
111
112 // otherwise, do float comparison of element and scvf center
113 else
114 {
115 const auto eps = lowDimGeometry.volume()*1e-8;
116 const auto diffVec = lowDimGeometry.center()-scvf.center();
117
118 if ( Dune::FloatCmp::eq<GlobalPosition, Dune::FloatCmp::CmpStyle::absolute>(diffVec, GlobalPosition(0.0), eps) )
119 {
120 embeddedScvfIdx = scvf.index();
121 found = true; break;
122 }
123 }
124 }
125
126 // Error tracking. The boundary scvf detection might has to be improved for very fine grids!?
127 if (!found)
128 DUNE_THROW(Dune::InvalidStateException, "Could not find coupling scvf in embedment");
129
130 // add each dof in the low dim element to coupling stencil of the bulk element
131 auto& bulkData = this->couplingMap_(bulkGridId, facetGridId)[bulkElemIdx];
132 const auto lowDimElementDofs = LowDimFVG::discMethod == DiscretizationMethod::box
133 ? this->extractNodalDofs_(lowDimElement, lowDimFvGridGeometry)
134 : std::vector<LowDimIndexType>( {lowDimElemIdx} );
135
136 for (auto dofIdx : lowDimElementDofs)
137 {
138 bulkData.couplingStencil.push_back( dofIdx );
139 bulkData.dofToCouplingScvfMap[dofIdx].push_back( embeddedScvfIdx );
140 }
141
142 // add info on which scvfs coincide with which low dim element
143 bulkData.couplingElementStencil.push_back(lowDimElemIdx);
144 bulkData.elementToScvfMap[lowDimElemIdx].push_back( embeddedScvfIdx );
145
146 // add embedment (coupling stencil will be done below)
147 lowDimData.embedments.emplace_back( bulkElemIdx, std::vector<BulkIndexType>({embeddedScvfIdx}) );
148 }
149
150 // adjoint entity indices = coupling stencil for tpfa
151 lowDimData.couplingStencil = std::move(adjoinedEntityIndices);
152 };
153
154 // let the parent do the update subject to the execution policy defined above
155 ParentType::update_(bulkFvGridGeometry, lowDimFvGridGeometry, embeddings, addCouplingEntryPolicy);
156
157 // coupling stencils might not be unique if box is used in lowdim domain
158 if (LowDimFVG::discMethod == DiscretizationMethod::box)
159 {
160 auto makeStencilUnique = [] (auto& data)
161 {
162 auto& cs = data.second.couplingStencil;
163 std::sort(cs.begin(), cs.end());
164 cs.erase( std::unique(cs.begin(), cs.end()), cs.end() );
165 };
166
167 auto& bulkCouplingData = this->couplingMap_(bulkGridId, facetGridId);
168 std::for_each(bulkCouplingData.begin(), bulkCouplingData.end(), makeStencilUnique);
169 }
170 }
171};
172
173} // end namespace Dumux
174
175#endif
Defines the index types used for grid and local indices.
The available discretization methods in Dumux.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
DiscretizationMethod
The available discretization methods in Dumux.
Definition: method.hh:37
Definition: adapt.hh:29
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:39
void update(const BulkFVG &bulkFvGridGeometry, const LowDimFVG &lowDimFvGridGeometry, std::shared_ptr< const Embeddings > embeddings)
Update coupling maps. This is the standard interface required by any mapper implementation.
Definition: facet/cellcentered/tpfa/couplingmapper.hh:69
Implementation for the coupling mapper that sets up and stores the coupling maps between two domains ...
Definition: facet/couplingmapper.hh:52
Base class for the coupling mapper that sets up and stores the coupling maps between two domains of d...
Definition: couplingmapperbase.hh:51