51 using PeerBlackList = BlackList::PeerBlackList;
52 using PeerBlackLists = BlackList::PeerBlackLists;
54 using Element =
typename GridView::template Codim<0>::Entity;
56 class BorderListHandle_
57 :
public Dune::CommDataHandleIF<BorderListHandle_, ProcessRank>
60 BorderListHandle_(
const GridView& gridView,
61 const ElementMapper& map,
66 , blackList_(blackList)
67 , borderList_(borderList)
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);
78 bool contains(
int,
int codim)
const
79 {
return codim == 0; }
81#if DUNE_VERSION_LT(DUNE_GRID, 2, 8)
82 bool fixedsize(
int,
int)
const
84 bool fixedSize(
int,
int)
const
88 template <
class EntityType>
89 size_t size(
const EntityType&)
const
92 template <
class MessageBufferImp,
class EntityType>
93 void gather(MessageBufferImp& buff,
const EntityType& e)
const
95 unsigned myIdx =
static_cast<unsigned>(map_.index(e));
96 buff.write(
static_cast<unsigned>(gridView_.comm().rank()));
100 template <
class MessageBufferImp>
101 void scatter(MessageBufferImp& buff,
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())
112 else if (isIt->outside().partitionType() == Dune::InteriorEntity) {
113 isInteriorNeighbor =
true;
117 if (!isInteriorNeighbor)
127 peerSet_.insert(tmp);
136 borderList_.push_back(bIdx);
142 template <
class MessageBufferImp,
class EntityType>
143 void scatter(MessageBufferImp&,
148 const std::set<ProcessRank>& peerSet()
const
153 const ElementMapper& map_;
154 std::set<ProcessRank> peerSet_;
159 class PeerBlackListHandle_
160 :
public Dune::CommDataHandleIF<PeerBlackListHandle_, int>
163 PeerBlackListHandle_(
const GridView& gridView,
164 const ElementMapper& map,
165 PeerBlackLists& peerBlackLists)
166 : gridView_(gridView)
168 , peerBlackLists_(peerBlackLists)
172 bool contains(
int,
int codim)
const
173 {
return codim == 0; }
175#if DUNE_VERSION_LT(DUNE_GRID, 2, 8)
176 bool fixedsize(
int,
int)
const
178 bool fixedSize(
int,
int)
const
182 template <
class EntityType>
183 size_t size(
const EntityType&)
const
186 template <
class MessageBufferImp,
class EntityType>
187 void gather(MessageBufferImp& buff,
const EntityType& e)
const
189 buff.write(
static_cast<int>(gridView_.comm().rank()));
190 buff.write(
static_cast<int>(map_.index(e)));
193 template <
class MessageBufferImp,
class EntityType>
194 void scatter(MessageBufferImp& buff,
const EntityType& e,
size_t)
202 localIdx =
static_cast<Index>(map_.index(e));
205 pIdx.nativeIndexOfPeer =
static_cast<Index>(peerIdx);
206 pIdx.myOwnNativeIndex =
static_cast<Index>(localIdx);
208 peerBlackLists_[
static_cast<ProcessRank>(peerRank)].push_back(pIdx);
211 const PeerBlackList& peerBlackList(
ProcessRank peerRank)
const
212 {
return peerBlackLists_.at(peerRank); }
216 const ElementMapper& map_;
217 PeerBlackLists peerBlackLists_;
222 : gridView_(gridView)
225 BorderListHandle_ blh(gridView, map, blackList_, borderList_);
226 gridView.communicate(blh,
227 Dune::InteriorBorder_All_Interface,
228 Dune::BackwardCommunication);
230 PeerBlackListHandle_ pblh(gridView, map, peerBlackLists_);
231 gridView.communicate(pblh,
232 Dune::InteriorBorder_All_Interface,
233 Dune::BackwardCommunication);
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));
244 {
return borderList_; }
248 {
return blackList_; }
251 const GridView gridView_;
252 const ElementMapper& map_;
257 PeerBlackLists peerBlackLists_;