10 #include <dune/common/parallel/remoteindices.hh> 11 #include <amdis/common/parallel/Communicator.hpp> 14 #include <amdis/Environment.hpp> 15 #include <amdis/linearalgebra/AttributeSet.hpp> 20 template <
class IS,
class GI = std::
size_t>
26 using size_type = std::size_t;
27 using DofIndex = size_type;
28 using LocalIndex = size_type;
29 using GlobalIndex = GI;
34 template <
class Communication>
85 assert(
false &&
"There are no ghost indices in sequential dofmappings");
92 assert(
false &&
"There are no ghost indices in sequential dofmappings");
97 GlobalIndex
global(LocalIndex
const& n)
const 137 template <
class Communication>
140 localSize_ = c.indexSet().size();
141 globalSize_ = c.indexSet().size();
142 indices_.resize(globalSize_);
143 std::iota(indices_.begin(), indices_.end(), size_type(0));
146 void debug()
const {}
149 size_type localSize_ = 0;
150 size_type globalSize_ = 0;
151 std::vector<GlobalIndex> indices_;
156 template <
class PIS,
class GI = std::
size_t>
158 class ParallelDofMapping
160 using ParallelIndexSet = PIS;
161 using Attribute =
typename AttributeSet<PIS>::type;
162 using RemoteIndices = Dune::RemoteIndices<ParallelIndexSet>;
165 using size_type = std::size_t;
166 using DofIndex = size_type;
167 using LocalIndex = size_type;
168 using GlobalIndex = GI;
171 ParallelDofMapping() =
default;
173 template <
class Communication>
174 ParallelDofMapping(Communication& c)
186 std::vector<size_type>
const&
localSizes()
const 206 return globalIndices_;
218 return ghostIndices_;
224 auto it = std::find(ghostIndices_.begin(), ghostIndices_.end(), n);
225 assert(it != ghostIndices_.end());
226 return std::distance(ghostIndices_.begin(), it);
230 LocalIndex
dofToGhost(DofIndex
const& n)
const 233 assert(n < ghostLocalIndices_.size());
234 assert(ghostLocalIndices_[n] < ghostSize_);
236 return ghostLocalIndices_[n];
240 GlobalIndex
global(DofIndex
const& n)
const 242 assert(n < globalIndices_.size());
243 return globalIndices_[n];
254 LocalIndex
dofToLocal(DofIndex
const& n)
const 256 assert(n < globalIndices_.size());
261 bool owner(DofIndex
const& n)
const 263 assert(n < owner_.size());
274 bool globalOwner(
int p, GlobalIndex
const& n)
const 277 return n >= starts_[p] && n < starts_[p+1];
281 template <
class Communication>
282 void update(Communication& c);
295 globalIndices_.clear();
296 ghostIndices_.clear();
297 ghostLocalIndices_.clear();
302 std::vector<size_type> sizes_;
303 std::vector<GlobalIndex> starts_;
304 size_type localSize_;
305 size_type globalSize_;
306 size_type ghostSize_;
308 std::vector<GlobalIndex> globalIndices_;
309 std::vector<GlobalIndex> ghostIndices_;
310 std::vector<LocalIndex> ghostLocalIndices_;
311 std::vector<bool> owner_;
313 const Mpi::Tag tag_{7513};
316 template <
class PIS,
class GI>
317 using DofMapping = ParallelDofMapping<PIS,GI>;
319 template <
class IS,
class GI>
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