AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
Communication.hpp
1 #pragma once
2 
3 #include <memory>
4 #include <optional>
5 #include <string>
6 
7 #include <dune/istl/owneroverlapcopy.hh>
8 #include <dune/istl/solvercategory.hh>
9 #include <dune/istl/paamg/pinfo.hh>
10 
11 #include <amdis/Environment.hpp>
12 #include <amdis/functions/GlobalIdSet.hpp>
13 #include <amdis/linearalgebra/Communication.hpp>
14 #include <amdis/linearalgebra/ParallelIndexSet.hpp>
15 
16 namespace AMDiS
17 {
19  template <class G, class L>
21  {
22  using Impl = Dune::Amg::SequentialInformation;
23 
24  public:
26 
27  SequentialISTLCommunication() = default;
28 
29  template <class Basis>
30  SequentialISTLCommunication(Basis const&, typename Dune::SolverCategory::Category cat)
31  {
32  assert(cat == Dune::SolverCategory::sequential);
33  }
34 
35  typename Dune::SolverCategory::Category category() const
36  {
37  return Dune::SolverCategory::sequential;
38  }
39 
40  Impl const& get() const
41  {
42  return impl_;
43  }
44 
45  Sequential const& sequential() const
46  {
47  return *this;
48  }
49 
50  template <class Basis>
51  void update(Basis const&) { /* do nothing */ }
52 
53  private:
54  Impl impl_;
55  };
56 
57 
58 #if HAVE_MPI
59  template <class GlobalId, class LocalIndex>
61  class ISTLCommunication
62  {
63  using Impl = Dune::OwnerOverlapCopyCommunication<GlobalId, LocalIndex>;
64 
65  public:
67  using IndexSet = typename Impl::ParallelIndexSet;
68  using RemoteIndices = typename Impl::RemoteIndices;
69 
70  public:
71  template <class Basis>
72  ISTLCommunication(Basis const& basis, typename Dune::SolverCategory::Category cat)
73  : cat_(cat)
74  {
75  update(basis);
76  }
77 
78  IndexSet const& indexSet() const
79  {
80  assert(bool(impl_));
81  return impl_->indexSet();
82  }
83 
84  RemoteIndices const& remoteIndices() const
85  {
86  assert(bool(impl_));
87  return impl_->remoteIndices();
88  }
89 
90  typename Dune::SolverCategory::Category category() const
91  {
92  return cat_;
93  }
94 
95  Impl const& get() const
96  {
97  assert(bool(impl_));
98  return *impl_;
99  }
100 
101  Sequential sequential() const
102  {
103  return Sequential{};
104  }
105 
106  template <class Basis>
107  void update(Basis const& basis)
108  {
109  impl_ = std::make_unique<Impl>(Environment::comm(), cat_);
110 
111  auto& pis = impl_->indexSet();
112  buildParallelIndexSet(basis, pis);
113 
114  auto& ris = impl_->remoteIndices();
115  ris.setIndexSets(pis, pis, impl_->communicator());
116  ris.template rebuild<true>();
117  }
118 
119  private:
120  typename Dune::SolverCategory::Category cat_;
121  std::unique_ptr<Impl> impl_ = nullptr;
122  };
123 
124 #else // !HAVE_MPI
125 
126  template <class G, class L>
128 
129 #endif
130 
131  template <class B>
132  using ISTLCommunication_t
133  = ISTLCommunication<GlobalIdType_t<B>, typename B::size_type>;
134 
135 
136  template <class G, class L>
138  {
140 
149  template <class Basis>
150  static Communication create(Basis const& basis, std::string const& prefix = "");
151  };
152 
153 } // end namespace AMDiS
154 
155 #include <amdis/linearalgebra/istl/Communication.inc.hpp>
Dummy implementation for ISTL-specific communication when no MPI is found.
Definition: Communication.hpp:20
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
auto cat(Dune::TypeTree::HybridTreePath< S... > const &tp0, Dune::TypeTree::HybridTreePath< T... > const &tp1)
Concatenate two treepaths.
Definition: TreePath.hpp:209
Implementation of a creator pattern for Communication types.
Definition: Communication.hpp:47
void buildParallelIndexSet(Basis const &basis, PIS &parallelIndexSet)
Fills a parallelIndexSet with indices from a basis.
Definition: ParallelIndexSet.hpp:24
static Dune::MPIHelper::MPICommunicator comm()
Return the MPI_Comm object (or a fake communicator)
Definition: Environment.hpp:86