vertexborderlistfromgrid.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_VERTEX_BORDER_LIST_FROM_GRID_HH
28#define EWOMS_VERTEX_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/common/version.hh>
36
37namespace Opm {
38namespace Linear {
48template <class GridView, class VertexMapper>
50 : public Dune::CommDataHandleIF<VertexBorderListFromGrid<GridView, VertexMapper>,
51 int>
52{
53 static const int dimWorld = GridView::dimensionworld;
54
55public:
56 VertexBorderListFromGrid(const GridView& gridView, const VertexMapper& map)
57 : gridView_(gridView), map_(map)
58 {
59 gridView.communicate(*this,
60 Dune::InteriorBorder_InteriorBorder_Interface,
61 Dune::ForwardCommunication);
62
63 auto vIt = gridView.template begin<dimWorld>();
64 const auto& vEndIt = gridView.template end<dimWorld >();
65 for (; vIt != vEndIt; ++vIt) {
66 if (vIt->partitionType() != Dune::InteriorEntity
67 && vIt->partitionType() != Dune::BorderEntity)
68 {
69 Index vIdx = static_cast<Index>(map_.index(*vIt));
70 blackList_.addIndex(vIdx);
71 }
72 }
73 }
74
75 // data handle methods
76 bool contains(int dim, int codim) const
77 { return dim == codim; }
78
79#if DUNE_VERSION_LT(DUNE_GRID, 2, 8)
80 bool fixedsize(int, int) const
81#else
82 bool fixedSize(int, int) const
83#endif
84 { return true; }
85
86 template <class EntityType>
87 size_t size(const EntityType&) const
88 { return 2; }
89
90 template <class MessageBufferImp, class EntityType>
91 void gather(MessageBufferImp& buff, const EntityType& e) const
92 {
93 buff.write(static_cast<int>(gridView_.comm().rank()));
94 buff.write(static_cast<int>(map_.index(e)));
95 }
96
97 template <class MessageBufferImp, class EntityType>
98 void scatter(MessageBufferImp& buff, const EntityType& e, size_t)
99 {
100 BorderIndex bIdx;
101
102 bIdx.localIdx = static_cast<Index>(map_.index(e));
103 {
104 int tmp;
105 buff.read(tmp);
106 bIdx.peerRank = static_cast<ProcessRank>(tmp);
107 }
108 {
109 int tmp;
110 buff.read(tmp);
111 bIdx.peerIdx = static_cast<Index>(tmp);
112 }
113 bIdx.borderDistance = 0;
114
115 borderList_.push_back(bIdx);
116 }
117
118 // Access to the border list.
119 const BorderList& borderList() const
120 { return borderList_; }
121
122 // Access to the black-list indices.
123 const BlackList& blackList() const
124 { return blackList_; }
125
126private:
127 const GridView gridView_;
128 const VertexMapper& map_;
129 BorderList borderList_;
130 BlackList blackList_;
131};
132
133} // namespace Linear
134} // namespace Opm
135
136#endif
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Definition: blacklist.hh:49
void addIndex(Index nativeIdx)
Definition: blacklist.hh:66
Uses communication on the grid to find the initial seed list of indices.
Definition: vertexborderlistfromgrid.hh:52
const BlackList & blackList() const
Definition: vertexborderlistfromgrid.hh:123
void scatter(MessageBufferImp &buff, const EntityType &e, size_t)
Definition: vertexborderlistfromgrid.hh:98
bool fixedSize(int, int) const
Definition: vertexborderlistfromgrid.hh:82
const BorderList & borderList() const
Definition: vertexborderlistfromgrid.hh:119
size_t size(const EntityType &) const
Definition: vertexborderlistfromgrid.hh:87
bool contains(int dim, int codim) const
Definition: vertexborderlistfromgrid.hh:76
VertexBorderListFromGrid(const GridView &gridView, const VertexMapper &map)
Definition: vertexborderlistfromgrid.hh:56
void gather(MessageBufferImp &buff, const EntityType &e) const
Definition: vertexborderlistfromgrid.hh:91
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
static const int dim
Definition: structuredgridvanguard.hh:67
Definition: blackoilboundaryratevector.hh:37
This files provides several data structures for storing tuples of indices of remote and/or local proc...
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