My Project
Loading...
Searching...
No Matches
elementborderlistfromgrid.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 This file is part of the Open Porous Media project (OPM).
5
6 OPM 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 2 of the License, or
9 (at your option) any later version.
10
11 OPM 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 OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef EWOMS_ELEMENT_BORDER_LIST_FROM_GRID_HH
28#define EWOMS_ELEMENT_BORDER_LIST_FROM_GRID_HH
29
30#include "overlaptypes.hh"
31#include "blacklist.hh"
32
33#include <dune/grid/common/datahandleif.hh>
34#include <dune/grid/common/gridenums.hh>
35#include <dune/grid/common/partitionset.hh>
36#include <dune/common/version.hh>
37
38#include <set>
39
40namespace Opm {
41namespace Linear {
47template <class GridView, class ElementMapper>
49{
51 using PeerBlackList = BlackList::PeerBlackList;
52 using PeerBlackLists = BlackList::PeerBlackLists;
53
54 using Element = typename GridView::template Codim<0>::Entity;
55
56 class BorderListHandle_
57 : public Dune::CommDataHandleIF<BorderListHandle_, ProcessRank>
58 {
59 public:
60 BorderListHandle_(const GridView& gridView,
61 const ElementMapper& map,
62 BlackList& blackList,
63 BorderList& borderList)
64 : gridView_(gridView)
65 , map_(map)
66 , blackList_(blackList)
67 , borderList_(borderList)
68 {
69 for (const auto& elem : elements(gridView_)) {
70 if (elem.partitionType() != Dune::InteriorEntity) {
71 Index elemIdx = static_cast<Index>(map_.index(elem));
72 blackList_.addIndex(elemIdx);
73 }
74 }
75 }
76
77 // data handle methods
78 bool contains(int, int codim) const
79 { return codim == 0; }
80
81#if DUNE_VERSION_LT(DUNE_GRID, 2, 8)
82 bool fixedsize(int, int) const
83#else
84 bool fixedSize(int, int) const
85#endif
86 { return true; }
87
88 template <class EntityType>
89 size_t size(const EntityType&) const
90 { return 2; }
91
92 template <class MessageBufferImp, class EntityType>
93 void gather(MessageBufferImp& buff, const EntityType& e) const
94 {
95 unsigned myIdx = static_cast<unsigned>(map_.index(e));
96 buff.write(static_cast<unsigned>(gridView_.comm().rank()));
97 buff.write(myIdx);
98 }
99
100 template <class MessageBufferImp>
101 void scatter(MessageBufferImp& buff,
102 const Element& e,
103 size_t)
104 {
105 // discard the index if it is not on the process boundary
106 bool isInteriorNeighbor = false;
107 auto isIt = gridView_.ibegin(e);
108 const auto& isEndIt = gridView_.iend(e);
109 for (; isIt != isEndIt; ++isIt) {
110 if (!isIt->neighbor())
111 continue;
112 else if (isIt->outside().partitionType() == Dune::InteriorEntity) {
113 isInteriorNeighbor = true;
114 break;
115 }
116 }
117 if (!isInteriorNeighbor)
118 return;
119
120 BorderIndex bIdx;
121
122 bIdx.localIdx = static_cast<Index>(map_.index(e));
123 {
124 ProcessRank tmp;
125 buff.read(tmp);
126 bIdx.peerRank = tmp;
127 peerSet_.insert(tmp);
128 }
129 {
130 unsigned tmp;
131 buff.read(tmp);
132 bIdx.peerIdx = static_cast<Index>(tmp);
133 }
134 bIdx.borderDistance = 1;
135
136 borderList_.push_back(bIdx);
137 }
138
139 // this template method is needed because the above one only works for codim-0
140 // entities (i.e., elements) but the dune grid uses some code which causes the
141 // compiler to invoke the scatter method for every codim...
142 template <class MessageBufferImp, class EntityType>
143 void scatter(MessageBufferImp&,
144 const EntityType&,
145 size_t)
146 { }
147
148 const std::set<ProcessRank>& peerSet() const
149 { return peerSet_; }
150
151 private:
152 GridView gridView_;
153 const ElementMapper& map_;
154 std::set<ProcessRank> peerSet_;
155 BlackList& blackList_;
156 BorderList& borderList_;
157 };
158
159 class PeerBlackListHandle_
160 : public Dune::CommDataHandleIF<PeerBlackListHandle_, int>
161 {
162 public:
163 PeerBlackListHandle_(const GridView& gridView,
164 const ElementMapper& map,
165 PeerBlackLists& peerBlackLists)
166 : gridView_(gridView)
167 , map_(map)
168 , peerBlackLists_(peerBlackLists)
169 {}
170
171 // data handle methods
172 bool contains(int, int codim) const
173 { return codim == 0; }
174
175#if DUNE_VERSION_LT(DUNE_GRID, 2, 8)
176 bool fixedsize(int, int) const
177#else
178 bool fixedSize(int, int) const
179#endif
180 { return true; }
181
182 template <class EntityType>
183 size_t size(const EntityType&) const
184 { return 2; }
185
186 template <class MessageBufferImp, class EntityType>
187 void gather(MessageBufferImp& buff, const EntityType& e) const
188 {
189 buff.write(static_cast<int>(gridView_.comm().rank()));
190 buff.write(static_cast<int>(map_.index(e)));
191 }
192
193 template <class MessageBufferImp, class EntityType>
194 void scatter(MessageBufferImp& buff, const EntityType& e, size_t)
195 {
196 int peerRank;
197 int peerIdx;
198 Index localIdx;
199
200 buff.read(peerRank);
201 buff.read(peerIdx);
202 localIdx = static_cast<Index>(map_.index(e));
203
205 pIdx.nativeIndexOfPeer = static_cast<Index>(peerIdx);
206 pIdx.myOwnNativeIndex = static_cast<Index>(localIdx);
207
208 peerBlackLists_[static_cast<ProcessRank>(peerRank)].push_back(pIdx);
209 }
210
211 const PeerBlackList& peerBlackList(ProcessRank peerRank) const
212 { return peerBlackLists_.at(peerRank); }
213
214 private:
215 GridView gridView_;
216 const ElementMapper& map_;
217 PeerBlackLists peerBlackLists_;
218 };
219
220public:
221 ElementBorderListFromGrid(const GridView& gridView, const ElementMapper& map)
222 : gridView_(gridView)
223 , map_(map)
224 {
225 BorderListHandle_ blh(gridView, map, blackList_, borderList_);
226 gridView.communicate(blh,
227 Dune::InteriorBorder_All_Interface,
228 Dune::BackwardCommunication);
229
230 PeerBlackListHandle_ pblh(gridView, map, peerBlackLists_);
231 gridView.communicate(pblh,
232 Dune::InteriorBorder_All_Interface,
233 Dune::BackwardCommunication);
234
235 auto peerIt = blh.peerSet().begin();
236 const auto& peerEndIt = blh.peerSet().end();
237 for (; peerIt != peerEndIt; ++peerIt) {
238 blackList_.setPeerList(*peerIt, pblh.peerBlackList(*peerIt));
239 }
240 }
241
242 // Access to the border list.
243 const BorderList& borderList() const
244 { return borderList_; }
245
246 // Access to the black-list indices.
247 const BlackList& blackList() const
248 { return blackList_; }
249
250private:
251 const GridView gridView_;
252 const ElementMapper& map_;
253
254 BorderList borderList_;
255
256 BlackList blackList_;
257 PeerBlackLists peerBlackLists_;
258};
259
260} // namespace Linear
261} // namespace Opm
262
263#endif
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Definition blacklist.hh:49
Uses communication on the grid to find the initial seed list of indices for methods which use element...
Definition elementborderlistfromgrid.hh:49
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
This files provides several data structures for storing tuples of indices of remote and/or local proc...
unsigned ProcessRank
The type of the rank of a process.
Definition overlaptypes.hh:49
int Index
The type of an index of a degree of freedom.
Definition overlaptypes.hh:44
std::list< BorderIndex > BorderList
This class managages a list of indices which are on the border of a process' partition of the grid.
Definition overlaptypes.hh:120
A single index intersecting with the process boundary.
Definition overlaptypes.hh:102
Index localIdx
Index of the entity for the local process.
Definition overlaptypes.hh:104
BorderDistance borderDistance
Distance to the process border for the peer (in hops)
Definition overlaptypes.hh:113
ProcessRank peerRank
Rank of the peer process.
Definition overlaptypes.hh:110
Index peerIdx
Index of the entity for the peer process.
Definition overlaptypes.hh:107