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 bool fixedSize(int, int) const
82 { return true; }
83
84 template <class EntityType>
85 size_t size(const EntityType&) const
86 { return 2; }
87
88 template <class MessageBufferImp, class EntityType>
89 void gather(MessageBufferImp& buff, const EntityType& e) const
90 {
91 unsigned myIdx = static_cast<unsigned>(map_.index(e));
92 buff.write(static_cast<unsigned>(gridView_.comm().rank()));
93 buff.write(myIdx);
94 }
95
96 template <class MessageBufferImp>
97 void scatter(MessageBufferImp& buff,
98 const Element& e,
99 size_t)
100 {
101 // discard the index if it is not on the process boundary
102 bool isInteriorNeighbor = false;
103 auto isIt = gridView_.ibegin(e);
104 const auto& isEndIt = gridView_.iend(e);
105 for (; isIt != isEndIt; ++isIt) {
106 if (!isIt->neighbor())
107 continue;
108 else if (isIt->outside().partitionType() == Dune::InteriorEntity) {
109 isInteriorNeighbor = true;
110 break;
111 }
112 }
114 return;
115
117
118 bIdx.localIdx = static_cast<Index>(map_.index(e));
119 {
120 ProcessRank tmp;
121 buff.read(tmp);
122 bIdx.peerRank = tmp;
123 peerSet_.insert(tmp);
124 }
125 {
126 unsigned tmp;
127 buff.read(tmp);
128 bIdx.peerIdx = static_cast<Index>(tmp);
129 }
130 bIdx.borderDistance = 1;
131
132 borderList_.push_back(bIdx);
133 }
134
135 // this template method is needed because the above one only works for codim-0
136 // entities (i.e., elements) but the dune grid uses some code which causes the
137 // compiler to invoke the scatter method for every codim...
138 template <class MessageBufferImp, class EntityType>
139 void scatter(MessageBufferImp&,
140 const EntityType&,
141 size_t)
142 { }
143
144 const std::set<ProcessRank>& peerSet() const
145 { return peerSet_; }
146
147 private:
148 GridView gridView_;
149 const ElementMapper& map_;
150 std::set<ProcessRank> peerSet_;
151 BlackList& blackList_;
152 BorderList& borderList_;
153 };
154
155 class PeerBlackListHandle_
156 : public Dune::CommDataHandleIF<PeerBlackListHandle_, int>
157 {
158 public:
159 PeerBlackListHandle_(const GridView& gridView,
160 const ElementMapper& map,
161 PeerBlackLists& peerBlackLists)
162 : gridView_(gridView)
163 , map_(map)
164 , peerBlackLists_(peerBlackLists)
165 {}
166
167 // data handle methods
168 bool contains(int, int codim) const
169 { return codim == 0; }
170
171 bool fixedSize(int, int) const
172 { return true; }
173
174 template <class EntityType>
175 size_t size(const EntityType&) const
176 { return 2; }
177
178 template <class MessageBufferImp, class EntityType>
179 void gather(MessageBufferImp& buff, const EntityType& e) const
180 {
181 buff.write(static_cast<int>(gridView_.comm().rank()));
182 buff.write(static_cast<int>(map_.index(e)));
183 }
184
185 template <class MessageBufferImp, class EntityType>
186 void scatter(MessageBufferImp& buff, const EntityType& e, size_t)
187 {
188 int peerRank;
189 int peerIdx;
190 Index localIdx;
191
192 buff.read(peerRank);
193 buff.read(peerIdx);
194 localIdx = static_cast<Index>(map_.index(e));
195
197 pIdx.nativeIndexOfPeer = static_cast<Index>(peerIdx);
198 pIdx.myOwnNativeIndex = static_cast<Index>(localIdx);
199
200 peerBlackLists_[static_cast<ProcessRank>(peerRank)].push_back(pIdx);
201 }
202
203 const PeerBlackList& peerBlackList(ProcessRank peerRank) const
204 { return peerBlackLists_.at(peerRank); }
205
206 private:
207 GridView gridView_;
208 const ElementMapper& map_;
209 PeerBlackLists peerBlackLists_;
210 };
211
212public:
213 ElementBorderListFromGrid(const GridView& gridView, const ElementMapper& map)
214 : gridView_(gridView)
215 , map_(map)
216 {
217 BorderListHandle_ blh(gridView, map, blackList_, borderList_);
218 gridView.communicate(blh,
219 Dune::InteriorBorder_All_Interface,
220 Dune::BackwardCommunication);
221
222 PeerBlackListHandle_ pblh(gridView, map, peerBlackLists_);
223 gridView.communicate(pblh,
224 Dune::InteriorBorder_All_Interface,
225 Dune::BackwardCommunication);
226
227 auto peerIt = blh.peerSet().begin();
228 const auto& peerEndIt = blh.peerSet().end();
229 for (; peerIt != peerEndIt; ++peerIt) {
230 blackList_.setPeerList(*peerIt, pblh.peerBlackList(*peerIt));
231 }
232 }
233
234 // Access to the border list.
235 const BorderList& borderList() const
236 { return borderList_; }
237
238 // Access to the black-list indices.
239 const BlackList& blackList() const
240 { return blackList_; }
241
242private:
243 const GridView gridView_;
244 const ElementMapper& map_;
245
246 BorderList borderList_;
247
248 BlackList blackList_;
249 PeerBlackLists peerBlackLists_;
250};
251
252} // namespace Linear
253} // namespace Opm
254
255#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
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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