OGS 6.2.1-97-g73d1aeda3
LocalToGlobalIndexMap.h
Go to the documentation of this file.
1 
10 #pragma once
11 
12 #ifndef NDEBUG
13 #include <iosfwd>
14 #endif // NDEBUG
15 
16 #include <vector>
17 
18 #include <Eigen/Dense>
19 
21 
22 #include "MeshComponentMap.h"
23 
24 namespace MeshLib
25 {
26 class Node;
27 class Element;
28 } // namespace MeshLib
29 
30 namespace NumLib
31 {
32 
43 {
44  // Enables using std::make_unique with private constructors from within
45  // member functions of LocalToGlobalIndexMap. Cf.
46  // http://seanmiddleditch.com/enabling-make_unique-with-private-constructors/
48  {
49  explicit ConstructorTag() = default;
50  };
51 
52 public:
55 
56 public:
62  LocalToGlobalIndexMap(std::vector<MeshLib::MeshSubset>&& mesh_subsets,
63  NumLib::ComponentOrder const order);
64 
73  LocalToGlobalIndexMap(std::vector<MeshLib::MeshSubset>&& mesh_subsets,
74  std::vector<int> const& vec_var_n_components,
75  NumLib::ComponentOrder const order);
76 
87  std::vector<MeshLib::MeshSubset>&& mesh_subsets,
88  std::vector<int> const& vec_var_n_components,
89  std::vector<std::vector<MeshLib::Element*> const*> const&
90  vec_var_elements,
91  NumLib::ComponentOrder const order);
92 
96  LocalToGlobalIndexMap* deriveBoundaryConstrainedMap(
97  int const variable_id,
98  std::vector<int> const& component_ids,
99  MeshLib::MeshSubset&& new_mesh_subset) const;
100 
105  std::unique_ptr<LocalToGlobalIndexMap> deriveBoundaryConstrainedMap(
106  MeshLib::MeshSubset&& new_mesh_subset) const;
107 
110  std::size_t dofSizeWithGhosts() const;
111 
115  std::size_t dofSizeWithoutGhosts() const;
116 
117  std::size_t size() const;
118 
119  int getNumberOfVariables() const;
120 
121  int getNumberOfVariableComponents(int variable_id) const;
122 
123  int getNumberOfComponents() const;
124 
125  RowColumnIndices operator()(std::size_t const mesh_item_id,
126  const int component_id) const;
127 
128  std::size_t getNumberOfElementDOF(std::size_t const mesh_item_id) const;
129 
130  std::size_t getNumberOfElementComponents(std::size_t const mesh_item_id) const;
131 
132  std::vector<int> getElementVariableIDs(
133  std::size_t const mesh_item_id) const;
134 
135  GlobalIndexType getGlobalIndex(MeshLib::Location const& l,
136  int const variable_id,
137  int const component_id) const;
138 
139  GlobalIndexType getGlobalIndex(MeshLib::Location const& l,
140  int const global_component_id) const;
141 
143  std::vector<GlobalIndexType> getGlobalIndices(
144  const MeshLib::Location& l) const;
145 
147  std::vector<GlobalIndexType> const& getGhostIndices() const;
148 
151  GlobalIndexType getLocalIndex(MeshLib::Location const& l,
152  std::size_t const comp_id,
153  std::size_t const range_begin,
154  std::size_t const range_end) const;
155 
156  MeshLib::MeshSubset const& getMeshSubset(int const variable_id,
157  int const component_id) const;
158 
159  MeshLib::MeshSubset const& getMeshSubset(
160  int const global_component_id) const;
161 
164  int getGlobalComponent(int const variable_id, int const component_id) const;
165 
172  explicit LocalToGlobalIndexMap(
173  std::vector<MeshLib::MeshSubset>&& mesh_subsets,
174  std::vector<int> const& global_component_ids,
175  std::vector<int> const& variable_component_offsets,
176  std::vector<MeshLib::Element*> const& elements,
177  NumLib::MeshComponentMap&& mesh_component_map,
178  ConstructorTag /*unused*/);
179 
180 private:
181  template <typename ElementIterator>
182  void findGlobalIndices(ElementIterator first, ElementIterator last,
183  std::vector<MeshLib::Node*> const& nodes,
184  std::size_t const mesh_id, const int comp_id,
185  const int comp_id_write);
186 
187  template <typename ElementIterator>
188  void findGlobalIndicesWithElementID(
189  ElementIterator first, ElementIterator last,
190  std::vector<MeshLib::Node*> const& nodes, std::size_t const mesh_id,
191  const int comp_id, const int comp_id_write);
192 
194  std::vector<MeshLib::MeshSubset> _mesh_subsets;
196 
197  using Table = Eigen::Matrix<LineIndex, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>;
198 
202 
204  Table const& _columns = _rows;
205 
206  std::vector<int> const _variable_component_offsets;
207 #ifndef NDEBUG
208  friend std::ostream& operator<<(std::ostream& os, LocalToGlobalIndexMap const& map);
210 #endif // NDEBUG
211 
212 };
213 
214 } // namespace NumLib
RowColumnIndices::LineIndex LineIndex
Multidirectional mapping between mesh entities and degrees of freedom.
ComponentOrder
Ordering of components in global matrix/vector.
std::ostream & operator<<(std::ostream &os, Element const &e)
Definition: Element.cpp:216
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
Interface for heuristic search length strategy.
Definition: ProjectData.h:30
Eigen::Matrix< LineIndex, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > Table
GlobalMatrix::IndexType GlobalIndexType
A subset of nodes on a single mesh.
Definition: MeshSubset.h:26
std::vector< int > const _variable_component_offsets
typename std::vector< IDX_TYPE > LineIndex