27 #ifndef EWOMS_ELEMENT_BORDER_LIST_FROM_GRID_HH
28 #define EWOMS_ELEMENT_BORDER_LIST_FROM_GRID_HH
33 #include <opm/common/Unused.hpp>
35 #include <dune/grid/common/datahandleif.hh>
36 #include <dune/grid/common/gridenums.hh>
37 #include <dune/istl/bcrsmatrix.hh>
38 #include <dune/istl/scalarproducts.hh>
39 #include <dune/istl/operators.hh>
40 #include <dune/common/version.hh>
51 template <
class Gr
idView,
class ElementMapper>
55 typedef BlackList::PeerBlackList PeerBlackList;
56 typedef BlackList::PeerBlackLists PeerBlackLists;
58 typedef typename GridView::template Codim<0>::Entity Element;
60 class BorderListHandle_
61 :
public Dune::CommDataHandleIF<BorderListHandle_, ProcessRank>
64 BorderListHandle_(
const GridView& gridView,
65 const ElementMapper& map,
67 BorderList& borderList)
70 , blackList_(blackList)
71 , borderList_(borderList)
73 auto elemIt = gridView_.template begin<0>();
74 const auto& elemEndIt = gridView_.template end<0>();
75 for (; elemIt != elemEndIt; ++elemIt) {
76 if (elemIt->partitionType() != Dune::InteriorEntity) {
77 Index elemIdx =
static_cast<Index
>(map_.index(*elemIt));
78 blackList_.addIndex(elemIdx);
84 bool contains(
int dim OPM_UNUSED,
int codim)
const
85 {
return codim == 0; }
87 bool fixedsize(
int dim OPM_UNUSED,
int codim OPM_UNUSED)
const
90 template <
class EntityType>
91 size_t size(
const EntityType& e OPM_UNUSED)
const
94 template <
class MessageBufferImp,
class EntityType>
95 void gather(MessageBufferImp& buff,
const EntityType& e)
const
97 unsigned myIdx =
static_cast<unsigned>(map_.index(e));
98 buff.write(static_cast<unsigned>(gridView_.comm().rank()));
102 template <
class MessageBufferImp>
103 void scatter(MessageBufferImp& buff,
108 bool isInteriorNeighbor =
false;
109 auto isIt = gridView_.ibegin(e);
110 const auto& isEndIt = gridView_.iend(e);
111 for (; isIt != isEndIt; ++isIt) {
112 if (!isIt->neighbor())
114 else if (isIt->outside().partitionType() == Dune::InteriorEntity) {
115 isInteriorNeighbor =
true;
119 if (!isInteriorNeighbor)
124 bIdx.
localIdx =
static_cast<Index
>(map_.index(e));
129 peerSet_.insert(tmp);
134 bIdx.
peerIdx =
static_cast<Index
>(tmp);
138 borderList_.push_back(bIdx);
144 template <
class MessageBufferImp,
class EntityType>
145 void scatter(MessageBufferImp& buff OPM_UNUSED,
146 const EntityType& e OPM_UNUSED,
150 const std::set<ProcessRank>& peerSet()
const
155 const ElementMapper& map_;
156 std::set<ProcessRank> peerSet_;
158 BorderList& borderList_;
161 class PeerBlackListHandle_
162 :
public Dune::CommDataHandleIF<PeerBlackListHandle_, int>
165 PeerBlackListHandle_(
const GridView& gridView,
166 const ElementMapper& map,
167 PeerBlackLists& peerBlackLists)
168 : gridView_(gridView)
170 , peerBlackLists_(peerBlackLists)
174 bool contains(
int dim OPM_UNUSED,
int codim)
const
175 {
return codim == 0; }
177 bool fixedsize(
int dim OPM_UNUSED,
int codim OPM_UNUSED)
const
180 template <
class EntityType>
181 size_t size(
const EntityType& e OPM_UNUSED)
const
184 template <
class MessageBufferImp,
class EntityType>
185 void gather(MessageBufferImp& buff,
const EntityType& e)
const
187 buff.write(static_cast<int>(gridView_.comm().rank()));
188 buff.write(static_cast<int>(map_.index(e)));
191 template <
class MessageBufferImp,
class EntityType>
192 void scatter(MessageBufferImp& buff,
const EntityType& e,
size_t n OPM_UNUSED)
200 localIdx =
static_cast<Index
>(map_.index(e));
203 pIdx.nativeIndexOfPeer =
static_cast<Index
>(peerIdx);
204 pIdx.myOwnNativeIndex =
static_cast<Index
>(localIdx);
206 peerBlackLists_[
static_cast<ProcessRank
>(peerRank)].push_back(pIdx);
209 const PeerBlackList& peerBlackList(ProcessRank peerRank)
const
210 {
return peerBlackLists_.at(peerRank); }
214 const ElementMapper& map_;
215 PeerBlackLists peerBlackLists_;
220 : gridView_(gridView)
223 BorderListHandle_ blh(gridView, map, blackList_, borderList_);
224 gridView.communicate(blh,
225 Dune::InteriorBorder_All_Interface,
226 Dune::BackwardCommunication);
228 PeerBlackListHandle_ pblh(gridView, map, peerBlackLists_);
229 gridView.communicate(pblh,
230 Dune::InteriorBorder_All_Interface,
231 Dune::BackwardCommunication);
233 auto peerIt = blh.peerSet().begin();
234 const auto& peerEndIt = blh.peerSet().end();
235 for (; peerIt != peerEndIt; ++peerIt) {
236 blackList_.setPeerList(*peerIt, pblh.peerBlackList(*peerIt));
241 const BorderList& borderList()
const
242 {
return borderList_; }
246 {
return blackList_; }
249 const GridView gridView_;
250 const ElementMapper& map_;
252 BorderList borderList_;
255 PeerBlackLists peerBlackLists_;
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Index peerIdx
Index of the entity for the peer process.
Definition: overlaptypes.hh:107
Index localIdx
Index of the entity for the local process.
Definition: overlaptypes.hh:104
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Definition: blacklist.hh:47
ProcessRank peerRank
Rank of the peer process.
Definition: overlaptypes.hh:110
This files provides several data structures for storing tuples of indices of remote and/or local proc...
BorderDistance borderDistance
Distance to the process border for the peer (in hops)
Definition: overlaptypes.hh:113
Definition: blacklist.hh:50
Uses communication on the grid to find the initial seed list of indices for methods which use element...
Definition: elementborderlistfromgrid.hh:52
A single index intersecting with the process boundary.
Definition: overlaptypes.hh:101