AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
Communication.hpp
1 #pragma once
2 
3 #include <optional>
4 
5 #include <petsc.h>
6 
7 #include <dune/common/parallel/indexset.hh>
8 #include <dune/common/parallel/localindex.hh>
9 #include <dune/common/parallel/plocalindex.hh>
10 #include <dune/common/parallel/remoteindices.hh>
11 #include <dune/common/typeutilities.hh>
12 
13 #include <amdis/functions/GlobalIdSet.hpp>
14 #include <amdis/linearalgebra/Communication.hpp>
15 #include <amdis/linearalgebra/DOFMapping.hpp>
16 #include <amdis/linearalgebra/ParallelIndexSet.hpp>
17 
18 namespace AMDiS
19 {
20  template <class GlobalId, class SizeType>
22  {
24  public:
25 #if HAVE_MPI
26  using LocalIndex = Dune::ParallelLocalIndex<DefaultAttributeSet::Type>;
27  using IndexSet = Dune::ParallelIndexSet<GlobalId, LocalIndex, 512>;
28  using RemoteIndices = Dune::RemoteIndices<IndexSet>;
29 #else
30  struct IndexSet
31  {
32  constexpr SizeType size() const { return size_; }
33  SizeType size_ = 0;
34  };
35 
36  struct RemoteIndices {};
37 #endif
39 
40  public:
41  template <class Basis,
42  Dune::disableCopyMove<Self, Basis> = 0>
43  DistributedCommunication(Basis const& basis)
44  : remoteIndices_(std::make_unique<RemoteIndices>())
45  {
46  update(basis);
47  }
48 
49  DistributedCommunication(Self const& other) = delete;
50  DistributedCommunication(Self&& other) = default;
51 
53  IndexSet const& indexSet() const { return *indexSet_; }
54  IndexSet& indexSet() { return *indexSet_; }
55 
57  RemoteIndices const& remoteIndices() const { return *remoteIndices_; }
58  RemoteIndices& remoteIndices() { return *remoteIndices_; }
59 
61  DofMap const& dofMap() const { return dofMap_; }
62  DofMap& dofMap() { return dofMap_; }
63 
65  template <class Basis>
66  void update(Basis const& basis)
67  {
68  indexSet_ = std::make_unique<IndexSet>();
69 
70 #if HAVE_MPI
71  buildParallelIndexSet(basis, *indexSet_);
72  remoteIndices_->setIndexSets(*indexSet_, *indexSet_, Environment::comm());
73  remoteIndices_->template rebuild<true>();
74 #else
75  (*indexSet_).size_ = basis.dimension();
76 #endif
77 
78  dofMap_.update(*this);
79  }
80 
81  protected:
82  std::unique_ptr<IndexSet> indexSet_ = nullptr;
83  std::unique_ptr<RemoteIndices> remoteIndices_;
84  DofMap dofMap_;
85  };
86 
87  template <class B>
89  = DistributedCommunication<GlobalIdType_t<B>, typename B::size_type>;
90 
91 
92  template <class G, class L>
94  {
96 
97  template <class Basis>
98  static Communication create(Basis const& basis, std::string const& prefix = "")
99  {
100  return {basis};
101  }
102  };
103 
104 } // end namespace AMDiS
Definition: Communication.hpp:30
IndexSet const & indexSet() const
Return the ParallelIndexSet.
Definition: Communication.hpp:53
RemoteIndices const & remoteIndices() const
Return the RemoteIndices.
Definition: Communication.hpp:57
Fallback for ParallelDofMapping in case there is only one mpi core.
Definition: DOFMapping.hpp:21
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
Implementation of a creator pattern for Communication types.
Definition: Communication.hpp:47
Definition: Communication.hpp:21
void update(Basis const &basis)
Update the indexSet, remoteIndices and dofmapping.
Definition: Communication.hpp:66
void buildParallelIndexSet(Basis const &basis, PIS &parallelIndexSet)
Fills a parallelIndexSet with indices from a basis.
Definition: ParallelIndexSet.hpp:24
DofMap const & dofMap() const
Return the DofMapping.
Definition: Communication.hpp:61
Definition: Communication.hpp:36
static Dune::MPIHelper::MPICommunicator comm()
Return the MPI_Comm object (or a fake communicator)
Definition: Environment.hpp:86