AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
DOFMapping.hpp
1 #pragma once
2 
3 #include <array>
4 #include <algorithm>
5 #include <iterator>
6 #include <numeric>
7 #include <vector>
8 
9 #if HAVE_MPI
10 #include <dune/common/parallel/remoteindices.hh>
11 #include <amdis/common/parallel/Communicator.hpp>
12 #endif
13 
14 #include <amdis/Environment.hpp>
15 #include <amdis/linearalgebra/AttributeSet.hpp>
16 
17 namespace AMDiS
18 {
20  template <class IS, class GI = std::size_t>
22  {
23  using IndexSet = IS;
24 
25  public:
26  using size_type = std::size_t;
27  using DofIndex = size_type;
28  using LocalIndex = size_type;
29  using GlobalIndex = GI;
30 
31  public:
32  SequentialDofMapping() = default;
33 
34  template <class Communication>
35  SequentialDofMapping(Communication& c)
36  {
37  update(c);
38  }
39 
41  size_type localSize() const
42  {
43  return localSize_;
44  }
45 
47  std::array<size_type,1> localSizes() const
48  {
49  return {localSize_};
50  }
51 
53  size_type globalSize() const
54  {
55  return globalSize_;
56  }
57 
59  std::array<GlobalIndex,1> globalStarts() const
60  {
61  return {0u};
62  }
63 
65  std::vector<GlobalIndex> const& globalIndices() const
66  {
67  return indices_;
68  }
69 
71  GlobalIndex ghostSize() const
72  {
73  return 0;
74  }
75 
77  std::array<GlobalIndex,0> ghostIndices() const
78  {
79  return {};
80  }
81 
83  LocalIndex globalToGhost(GlobalIndex const& n) const
84  {
85  assert(false && "There are no ghost indices in sequential dofmappings");
86  return 0;
87  }
88 
90  LocalIndex dofToGhost(DofIndex const& n) const
91  {
92  assert(false && "There are no ghost indices in sequential dofmappings");
93  return 0;
94  }
95 
97  GlobalIndex global(LocalIndex const& n) const
98  {
99  return n;
100  }
101 
103  LocalIndex globalToLocal(GlobalIndex const& n) const
104  {
105  return n;
106  }
107 
109  LocalIndex dofToLocal(DofIndex const& n) const
110  {
111  return n;
112  }
113 
114 
116  bool owner(DofIndex const& n) const
117  {
118  assert(n < localSize());
119  return true;
120  }
121 
123  bool globalOwner(GlobalIndex const& n) const
124  {
125  return globalOwner(0, n);
126  }
127 
129  bool globalOwner(int p, GlobalIndex const& n) const
130  {
131  assert(p == 0);
132  assert(n < globalSize());
133  return true;
134  }
135 
137  template <class Communication>
138  void update(Communication& c)
139  {
140  localSize_ = c.indexSet().size();
141  globalSize_ = c.indexSet().size();
142  indices_.resize(globalSize_);
143  std::iota(indices_.begin(), indices_.end(), size_type(0));
144  }
145 
146  void debug() const {}
147 
148  private:
149  size_type localSize_ = 0;
150  size_type globalSize_ = 0;
151  std::vector<GlobalIndex> indices_;
152  };
153 
154 
155 #if HAVE_MPI
156  template <class PIS, class GI = std::size_t>
158  class ParallelDofMapping
159  {
160  using ParallelIndexSet = PIS;
161  using Attribute = typename AttributeSet<PIS>::type;
162  using RemoteIndices = Dune::RemoteIndices<ParallelIndexSet>;
163 
164  public:
165  using size_type = std::size_t;
166  using DofIndex = size_type;
167  using LocalIndex = size_type;
168  using GlobalIndex = GI;
169 
170  public:
171  ParallelDofMapping() = default;
172 
173  template <class Communication>
174  ParallelDofMapping(Communication& c)
175  {
176  update(c);
177  }
178 
180  size_type localSize() const
181  {
182  return localSize_;
183  }
184 
186  std::vector<size_type> const& localSizes() const
187  {
188  return sizes_;
189  }
190 
192  size_type globalSize() const
193  {
194  return globalSize_;
195  }
196 
198  std::vector<GlobalIndex> const& globalStarts() const
199  {
200  return starts_;
201  }
202 
204  std::vector<GlobalIndex> const& globalIndices() const
205  {
206  return globalIndices_;
207  }
208 
210  size_type ghostSize() const
211  {
212  return ghostSize_;
213  }
214 
216  std::vector<GlobalIndex> const& ghostIndices() const
217  {
218  return ghostIndices_;
219  }
220 
222  LocalIndex globalToGhost(GlobalIndex const& n) const
223  {
224  auto it = std::find(ghostIndices_.begin(), ghostIndices_.end(), n);
225  assert(it != ghostIndices_.end());
226  return std::distance(ghostIndices_.begin(), it);
227  }
228 
230  LocalIndex dofToGhost(DofIndex const& n) const
231  {
232  assert(!owner(n));
233  assert(n < ghostLocalIndices_.size());
234  assert(ghostLocalIndices_[n] < ghostSize_);
235 
236  return ghostLocalIndices_[n];
237  }
238 
240  GlobalIndex global(DofIndex const& n) const
241  {
242  assert(n < globalIndices_.size());
243  return globalIndices_[n];
244  }
245 
246 
248  LocalIndex globalToLocal(GlobalIndex const& n) const
249  {
250  return n - starts_[Environment::mpiRank()];
251  }
252 
254  LocalIndex dofToLocal(DofIndex const& n) const
255  {
256  assert(n < globalIndices_.size());
257  return globalToLocal(globalIndices_[n]);
258  }
259 
261  bool owner(DofIndex const& n) const
262  {
263  assert(n < owner_.size());
264  return owner_[n];
265  }
266 
268  bool globalOwner(GlobalIndex const& n) const
269  {
270  return globalOwner(Environment::mpiRank(), n);
271  }
272 
274  bool globalOwner(int p, GlobalIndex const& n) const
275  {
276  assert(p < Environment::mpiSize());
277  return n >= starts_[p] && n < starts_[p+1];
278  }
279 
281  template <class Communication>
282  void update(Communication& c);
283 
284  void debug() const;
285 
286  private:
287  void reset()
288  {
289  sizes_.clear();
290  starts_.clear();
291  localSize_ = 0;
292  globalSize_ = 0;
293  ghostSize_ = 0;
294 
295  globalIndices_.clear();
296  ghostIndices_.clear();
297  ghostLocalIndices_.clear();
298  owner_.clear();
299  }
300 
301  private:
302  std::vector<size_type> sizes_;
303  std::vector<GlobalIndex> starts_;
304  size_type localSize_;
305  size_type globalSize_;
306  size_type ghostSize_;
307 
308  std::vector<GlobalIndex> globalIndices_; // indexed by LocalIndex
309  std::vector<GlobalIndex> ghostIndices_;
310  std::vector<LocalIndex> ghostLocalIndices_;
311  std::vector<bool> owner_;
312 
313  const Mpi::Tag tag_{7513};
314  };
315 
316  template <class PIS,class GI>
317  using DofMapping = ParallelDofMapping<PIS,GI>;
318 #else
319  template <class IS, class GI>
321 #endif
322 
323 } // end namespace AMDiS
324 
325 #include <amdis/linearalgebra/DOFMapping.inc.hpp>
bool globalOwner(int p, GlobalIndex const &n) const
Global index n is owned by processor p.
Definition: DOFMapping.hpp:129
void update(Communication &c)
Update the local to global mapping. Must be called before mapping local to global.
Definition: DOFMapping.hpp:138
Fallback for ParallelDofMapping in case there is only one mpi core.
Definition: DOFMapping.hpp:21
size_type localSize() const
How many DOFs are owned by my processor?
Definition: DOFMapping.hpp:41
LocalIndex globalToLocal(GlobalIndex const &n) const
Map global index to consecutive local owner index.
Definition: DOFMapping.hpp:103
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
std::array< size_type, 1 > localSizes() const
Return the sequence of number of local indices for all processors.
Definition: DOFMapping.hpp:47
static int mpiSize()
Return the MPI_Size of the group created by Dune::MPIHelper.
Definition: Environment.hpp:74
std::array< GlobalIndex, 1 > globalStarts() const
Return the sequence of starting points of the global indices for all processors.
Definition: DOFMapping.hpp:59
LocalIndex dofToLocal(DofIndex const &n) const
Map DOF index to consecutive local owner index.
Definition: DOFMapping.hpp:109
std::vector< GlobalIndex > const & globalIndices() const
Return the vector of global indices.
Definition: DOFMapping.hpp:65
LocalIndex globalToGhost(GlobalIndex const &n) const
Map global index to local ghost index.
Definition: DOFMapping.hpp:83
LocalIndex dofToGhost(DofIndex const &n) const
Map DOF index to local ghost index.
Definition: DOFMapping.hpp:90
GlobalIndex ghostSize() const
Return number of ghost indices.
Definition: DOFMapping.hpp:71
size_type globalSize() const
The total number of global DOFs.
Definition: DOFMapping.hpp:53
std::array< GlobalIndex, 0 > ghostIndices() const
Return the vector of ghost indices.
Definition: DOFMapping.hpp:77
bool globalOwner(GlobalIndex const &n) const
Global index n is owned by this processor.
Definition: DOFMapping.hpp:123
bool owner(DofIndex const &n) const
DOF index n is owned by this processor.
Definition: DOFMapping.hpp:116
GlobalIndex global(LocalIndex const &n) const
Global index of local index n.
Definition: DOFMapping.hpp:97
static int mpiRank()
Return the MPI_Rank of the current processor.
Definition: Environment.hpp:68