OGS 6.2.1-76-gbb689931b
NumLib::LocalToGlobalIndexMap Class Referencefinal

Detailed Description

Row and column indices in global linear algebra objects for each mesh item.

The row and column indices in global matrix and rhs vector for example, required for addition of local contributions from every mesh item (like node or cell) to global objects.

The number of rows should be equal to the number of mesh items and the number of columns should be equal to the number of the components on that mesh item.

Definition at line 42 of file LocalToGlobalIndexMap.h.

#include <LocalToGlobalIndexMap.h>

Collaboration diagram for NumLib::LocalToGlobalIndexMap:

Classes

struct  ConstructorTag
 

Public Types

using RowColumnIndices = MathLib::RowColumnIndices< GlobalIndexType >
 
using LineIndex = RowColumnIndices::LineIndex
 

Public Member Functions

 LocalToGlobalIndexMap (std::vector< MeshLib::MeshSubset > &&mesh_subsets, NumLib::ComponentOrder const order)
 
 LocalToGlobalIndexMap (std::vector< MeshLib::MeshSubset > &&mesh_subsets, std::vector< int > const &vec_var_n_components, NumLib::ComponentOrder const order)
 
 LocalToGlobalIndexMap (std::vector< MeshLib::MeshSubset > &&mesh_subsets, std::vector< int > const &vec_var_n_components, std::vector< std::vector< MeshLib::Element *> const *> const &vec_var_elements, NumLib::ComponentOrder const order)
 
LocalToGlobalIndexMapderiveBoundaryConstrainedMap (int const variable_id, std::vector< int > const &component_ids, MeshLib::MeshSubset &&new_mesh_subset) const
 
std::unique_ptr< LocalToGlobalIndexMapderiveBoundaryConstrainedMap (MeshLib::MeshSubset &&new_mesh_subset) const
 
std::size_t dofSizeWithGhosts () const
 
std::size_t dofSizeWithoutGhosts () const
 
std::size_t size () const
 
int getNumberOfVariables () const
 
int getNumberOfVariableComponents (int variable_id) const
 
int getNumberOfComponents () const
 
RowColumnIndices operator() (std::size_t const mesh_item_id, const int component_id) const
 
std::size_t getNumberOfElementDOF (std::size_t const mesh_item_id) const
 
std::size_t getNumberOfElementComponents (std::size_t const mesh_item_id) const
 
std::vector< int > getElementVariableIDs (std::size_t const mesh_item_id) const
 
GlobalIndexType getGlobalIndex (MeshLib::Location const &l, int const variable_id, int const component_id) const
 
GlobalIndexType getGlobalIndex (MeshLib::Location const &l, int const global_component_id) const
 
std::vector< GlobalIndexTypegetGlobalIndices (const MeshLib::Location &l) const
 Forwards the respective method from MeshComponentMap. More...
 
std::vector< GlobalIndexType > const & getGhostIndices () const
 Get ghost indices, forwarded from MeshComponentMap. More...
 
GlobalIndexType getLocalIndex (MeshLib::Location const &l, std::size_t const comp_id, std::size_t const range_begin, std::size_t const range_end) const
 
MeshLib::MeshSubset const & getMeshSubset (int const variable_id, int const component_id) const
 
MeshLib::MeshSubset const & getMeshSubset (int const global_component_id) const
 
int getGlobalComponent (int const variable_id, int const component_id) const
 
 LocalToGlobalIndexMap (std::vector< MeshLib::MeshSubset > &&mesh_subsets, std::vector< int > const &global_component_ids, std::vector< int > const &variable_component_offsets, std::vector< MeshLib::Element *> const &elements, NumLib::MeshComponentMap &&mesh_component_map, ConstructorTag)
 

Private Types

using Table = Eigen::Matrix< LineIndex, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >
 

Private Member Functions

template<typename ElementIterator >
void findGlobalIndices (ElementIterator first, ElementIterator last, std::vector< MeshLib::Node *> const &nodes, std::size_t const mesh_id, const int comp_id, const int comp_id_write)
 
template<typename ElementIterator >
void findGlobalIndicesWithElementID (ElementIterator first, ElementIterator last, std::vector< MeshLib::Node *> const &nodes, std::size_t const mesh_id, const int comp_id, const int comp_id_write)
 

Private Attributes

std::vector< MeshLib::MeshSubset_mesh_subsets
 A vector of mesh subsets for each process variables' components. More...
 
NumLib::MeshComponentMap _mesh_component_map
 
Table _rows
 
Table const & _columns = _rows
 
std::vector< int > const _variable_component_offsets
 

Friends

std::ostream & operator<< (std::ostream &os, LocalToGlobalIndexMap const &map)
 Prints first rows of the table, every line, and the mesh component map. More...
 

Member Typedef Documentation

◆ LineIndex

◆ RowColumnIndices

◆ Table

using NumLib::LocalToGlobalIndexMap::Table = Eigen::Matrix<LineIndex, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
private

Definition at line 196 of file LocalToGlobalIndexMap.h.

Constructor & Destructor Documentation

◆ LocalToGlobalIndexMap() [1/4]

NumLib::LocalToGlobalIndexMap::LocalToGlobalIndexMap ( std::vector< MeshLib::MeshSubset > &&  mesh_subsets,
NumLib::ComponentOrder const  order 
)

Creates a MeshComponentMap internally and stores the global indices for each mesh element of the given mesh_subsets.

Attention
This constructor assumes the number of the given mesh subsets is equal to the number of variables, i.e. every variable has a single component.

Definition at line 115 of file LocalToGlobalIndexMap.cpp.

Referenced by deriveBoundaryConstrainedMap().

118  : LocalToGlobalIndexMap(std::move(mesh_subsets),
119  std::vector<int>(mesh_subsets.size(), 1), order)
120 {
121 }
LocalToGlobalIndexMap(std::vector< MeshLib::MeshSubset > &&mesh_subsets, NumLib::ComponentOrder const order)

◆ LocalToGlobalIndexMap() [2/4]

NumLib::LocalToGlobalIndexMap::LocalToGlobalIndexMap ( std::vector< MeshLib::MeshSubset > &&  mesh_subsets,
std::vector< int > const &  vec_var_n_components,
NumLib::ComponentOrder const  order 
)

Creates a MeshComponentMap internally and stores the global indices for each mesh element of the given mesh_subsets.

Parameters
mesh_subsetsa vector of components
vec_var_n_componentsa vector of the number of variable components. The size of the vector should be equal to the number of variables. Sum of the entries should be equal to the size of the mesh_subsets.
ordertype of ordering values in a vector

Definition at line 123 of file LocalToGlobalIndexMap.cpp.

References _mesh_subsets, findGlobalIndices(), and getGlobalComponent().

127  : _mesh_subsets(std::move(mesh_subsets)),
129  _variable_component_offsets(to_cumulative(vec_var_n_components))
130 {
131  // For each element of that MeshSubset save a line of global indices.
132  for (int variable_id = 0; variable_id < static_cast<int>(vec_var_n_components.size());
133  ++variable_id)
134  {
135  for (int component_id = 0; component_id < static_cast<int>(vec_var_n_components[variable_id]);
136  ++component_id)
137  {
138  auto const global_component_id =
139  getGlobalComponent(variable_id, component_id);
140 
141  auto const& ms = _mesh_subsets[global_component_id];
142  std::size_t const mesh_id = ms.getMeshID();
143 
144  findGlobalIndices(ms.elementsBegin(), ms.elementsEnd(),
145  ms.getNodes(), mesh_id, global_component_id,
146  global_component_id);
147  }
148  }
149 }
int getGlobalComponent(int const variable_id, int const component_id) const
void findGlobalIndices(ElementIterator first, ElementIterator last, std::vector< MeshLib::Node *> const &nodes, std::size_t const mesh_id, const int comp_id, const int comp_id_write)
std::vector< T > to_cumulative(std::vector< T > const &vec)
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
std::vector< int > const _variable_component_offsets

◆ LocalToGlobalIndexMap() [3/4]

NumLib::LocalToGlobalIndexMap::LocalToGlobalIndexMap ( std::vector< MeshLib::MeshSubset > &&  mesh_subsets,
std::vector< int > const &  vec_var_n_components,
std::vector< std::vector< MeshLib::Element *> const *> const &  vec_var_elements,
NumLib::ComponentOrder const  order 
)

Creates a MeshComponentMap internally and stores the global indices for the given mesh elements

Parameters
mesh_subsetsa vector of components
vec_var_n_componentsa vector of the number of variable components. The size of the vector should be equal to the number of variables. Sum of the entries should be equal to the size of the mesh_subsets.
vec_var_elementsa vector of active mesh elements for each variable.
ordertype of ordering values in a vector

Definition at line 151 of file LocalToGlobalIndexMap.cpp.

References _mesh_subsets, _rows, findGlobalIndicesWithElementID(), and getGlobalComponent().

156  : _mesh_subsets(std::move(mesh_subsets)),
158  _variable_component_offsets(to_cumulative(vec_var_n_components))
159 {
160  assert(vec_var_n_components.size() == vec_var_elements.size());
161 
162  // For each element of that MeshSubset save a line of global indices.
163 
164  // _rows should be resized based on an element ID
165  std::size_t max_elem_id = 0;
166  for (std::vector<MeshLib::Element*>const* eles : vec_var_elements)
167  {
168  for (auto e : *eles)
169  {
170  max_elem_id = std::max(max_elem_id, e->getID());
171  }
172  }
173  _rows.resize(max_elem_id + 1, _mesh_subsets.size());
174 
175  for (int variable_id = 0; variable_id < static_cast<int>(vec_var_n_components.size());
176  ++variable_id)
177  {
178  std::vector<MeshLib::Element*> const& var_elements = *vec_var_elements[variable_id];
179  for (int component_id = 0; component_id < static_cast<int>(vec_var_n_components[variable_id]);
180  ++component_id)
181  {
182  auto const global_component_id =
183  getGlobalComponent(variable_id, component_id);
184 
185  auto const& ms = _mesh_subsets[global_component_id];
186  std::size_t const mesh_id = ms.getMeshID();
187 
189  var_elements.cbegin(), var_elements.cend(), ms.getNodes(),
190  mesh_id, global_component_id, global_component_id);
191  }
192  }
193 }
int getGlobalComponent(int const variable_id, int const component_id) const
void findGlobalIndicesWithElementID(ElementIterator first, ElementIterator last, std::vector< MeshLib::Node *> const &nodes, std::size_t const mesh_id, const int comp_id, const int comp_id_write)
std::vector< T > to_cumulative(std::vector< T > const &vec)
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
std::vector< int > const _variable_component_offsets

◆ LocalToGlobalIndexMap() [4/4]

NumLib::LocalToGlobalIndexMap::LocalToGlobalIndexMap ( std::vector< MeshLib::MeshSubset > &&  mesh_subsets,
std::vector< int > const &  global_component_ids,
std::vector< int > const &  variable_component_offsets,
std::vector< MeshLib::Element *> const &  elements,
NumLib::MeshComponentMap &&  mesh_component_map,
LocalToGlobalIndexMap::ConstructorTag   
)
explicit

Private constructor (ensured by ConstructorTag) used by internally created local-to-global index maps. The mesh_component_map is passed as argument instead of being created by the constructor.

Attention
The passed mesh_component_map is in undefined state after this construtor.

Definition at line 195 of file LocalToGlobalIndexMap.cpp.

References _mesh_subsets, findGlobalIndices(), and OGS_FATAL.

202  : _mesh_subsets(std::move(mesh_subsets)),
203  _mesh_component_map(std::move(mesh_component_map)),
204  _variable_component_offsets(variable_component_offsets)
205 {
206  // Each subset in the mesh_subsets represents a single component.
207  if (_mesh_subsets.size() != global_component_ids.size())
208  {
209  OGS_FATAL(
210  "Number of mesh subsets is not equal to number of components. "
211  "There are %d mesh subsets and %d components.",
212  _mesh_subsets.size(), global_component_ids.size());
213  }
214 
215  for (int i = 0; i < static_cast<int>(global_component_ids.size()); ++i)
216  {
217  auto const& ms = _mesh_subsets[i];
218 
219  // For all MeshSubset in mesh_subsets and each element of that
220  // MeshSubset save a line of global indices.
221  std::size_t const mesh_id = ms.getMeshID();
222 
223  findGlobalIndices(elements.cbegin(), elements.cend(), ms.getNodes(),
224  mesh_id, global_component_ids[i], i);
225  }
226 }
void findGlobalIndices(ElementIterator first, ElementIterator last, std::vector< MeshLib::Node *> const &nodes, std::size_t const mesh_id, const int comp_id, const int comp_id_write)
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::vector< int > const _variable_component_offsets

Member Function Documentation

◆ deriveBoundaryConstrainedMap() [1/2]

LocalToGlobalIndexMap * NumLib::LocalToGlobalIndexMap::deriveBoundaryConstrainedMap ( int const  variable_id,
std::vector< int > const &  component_ids,
MeshLib::MeshSubset &&  new_mesh_subset 
) const

Derive a LocalToGlobalIndexMap constrained to the mesh subset and mesh subset's elements. A new mesh component map will be constructed using the passed mesh_subset for the given variable and component ids.

Definition at line 228 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, _mesh_subsets, _variable_component_offsets, getGlobalComponent(), NumLib::MeshComponentMap::getSubset(), LocalToGlobalIndexMap(), and OGS_FATAL.

Referenced by ProcessLib::DirichletBoundaryConditionWithinTimeInterval::config(), ProcessLib::ConstraintDirichletBoundaryCondition::ConstraintDirichletBoundaryCondition(), ProcessLib::createSourceTerm(), ProcessLib::createVariableDependentNeumannBoundaryCondition(), ProcessLib::DirichletBoundaryCondition::DirichletBoundaryCondition(), ProcessLib::Output::doOutputAlways(), ProcessLib::GenericNaturalBoundaryCondition< BoundaryConditionData, LocalAssemblerImplementation >::GenericNaturalBoundaryCondition(), ProcessLib::NormalTractionBoundaryCondition::NormalTractionBoundaryCondition< LocalAssemblerImplementation >::NormalTractionBoundaryCondition(), and ProcessLib::PythonBoundaryCondition::PythonBoundaryCondition().

232 {
233  DBUG("Construct reduced local to global index map.");
234 
235  if (component_ids.empty())
236  {
237  OGS_FATAL("Expected non-empty vector of component ids.");
238  }
239 
240  // Elements of the new_mesh_subset's mesh.
241  std::vector<MeshLib::Element*> const& elements =
242  new_mesh_subset.getMesh().getElements();
243 
244  // Create a subset of the current mesh component map.
245  std::vector<int> global_component_ids;
246 
247  for (auto component_id : component_ids)
248  {
249  global_component_ids.push_back(
250  getGlobalComponent(variable_id, component_id));
251  }
252 
253  auto mesh_component_map = _mesh_component_map.getSubset(
254  _mesh_subsets, new_mesh_subset, global_component_ids);
255 
256  // Create copies of the new_mesh_subset for each of the global components.
257  // The last component is moved after the for-loop.
258  std::vector<MeshLib::MeshSubset> all_mesh_subsets;
259  for (int i = 0; i < static_cast<int>(global_component_ids.size()) - 1; ++i)
260  {
261  all_mesh_subsets.emplace_back(new_mesh_subset);
262  }
263  all_mesh_subsets.emplace_back(std::move(new_mesh_subset));
264 
265  return new LocalToGlobalIndexMap(
266  std::move(all_mesh_subsets), global_component_ids,
267  _variable_component_offsets, elements, std::move(mesh_component_map),
268  ConstructorTag{});
269 }
int getGlobalComponent(int const variable_id, int const component_id) const
LocalToGlobalIndexMap(std::vector< MeshLib::MeshSubset > &&mesh_subsets, NumLib::ComponentOrder const order)
MeshComponentMap getSubset(std::vector< MeshLib::MeshSubset > const &bulk_mesh_subsets, MeshLib::MeshSubset const &new_mesh_subset, std::vector< int > const &new_global_component_ids) const
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
#define OGS_FATAL(fmt,...)
Definition: Error.h:63
std::vector< int > const _variable_component_offsets
Mesh const & getMesh() const
Definition: MeshSubset.h:106

◆ deriveBoundaryConstrainedMap() [2/2]

std::unique_ptr< LocalToGlobalIndexMap > NumLib::LocalToGlobalIndexMap::deriveBoundaryConstrainedMap ( MeshLib::MeshSubset &&  new_mesh_subset) const

Derive a LocalToGlobalIndexMap constrained to the mesh subset and mesh subset's elements. A new mesh component map will be constructed using the passed mesh_subset for all variables and components of the current LocalToGlobalIndexMap.

Definition at line 272 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, _mesh_subsets, _variable_component_offsets, getNumberOfComponents(), and NumLib::MeshComponentMap::getSubset().

274 {
275  DBUG("Construct reduced local to global index map.");
276 
277  // Create a subset of the current mesh component map.
278  std::vector<int> global_component_ids;
279 
280  for (int i = 0; i < getNumberOfComponents(); ++i)
281  {
282  global_component_ids.push_back(i);
283  }
284 
285  // Elements of the new_mesh_subset's mesh.
286  std::vector<MeshLib::Element*> const& elements =
287  new_mesh_subset.getMesh().getElements();
288 
289  auto mesh_component_map = _mesh_component_map.getSubset(
290  _mesh_subsets, new_mesh_subset, global_component_ids);
291 
292  // Create copies of the new_mesh_subset for each of the global components.
293  // The last component is moved after the for-loop.
294  std::vector<MeshLib::MeshSubset> all_mesh_subsets;
295  for (int i = 0; i < static_cast<int>(global_component_ids.size()) - 1; ++i)
296  {
297  all_mesh_subsets.emplace_back(new_mesh_subset);
298  }
299  all_mesh_subsets.emplace_back(std::move(new_mesh_subset));
300 
301  return std::make_unique<LocalToGlobalIndexMap>(
302  std::move(all_mesh_subsets), global_component_ids,
303  _variable_component_offsets, elements, std::move(mesh_component_map),
304  ConstructorTag{});
305 }
MeshComponentMap getSubset(std::vector< MeshLib::MeshSubset > const &bulk_mesh_subsets, MeshLib::MeshSubset const &new_mesh_subset, std::vector< int > const &new_global_component_ids) const
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition: Mesh.h:108
std::vector< int > const _variable_component_offsets
Mesh const & getMesh() const
Definition: MeshSubset.h:106

◆ dofSizeWithGhosts()

std::size_t NumLib::LocalToGlobalIndexMap::dofSizeWithGhosts ( ) const

Returns total number of degrees of freedom including those located in the ghost nodes.

Definition at line 307 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::dofSizeWithGhosts().

Referenced by computeSparsityPatternNonPETSc().

308 {
310 }
std::size_t dofSizeWithGhosts() const
The number of dofs including the those located in the ghost nodes.
NumLib::MeshComponentMap _mesh_component_map

◆ dofSizeWithoutGhosts()

std::size_t NumLib::LocalToGlobalIndexMap::dofSizeWithoutGhosts ( ) const

Returns total number of local degrees of freedom of the present rank, which does not count the unknowns associated with ghost nodes (for DDC with node-wise mesh partitioning).

Definition at line 312 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::dofSizeWithoutGhosts().

Referenced by ProcessLib::GlobalVectorFromNamedFunction::call(), and NumLib::LocalLinearLeastSquaresExtrapolator::extrapolate().

313 {
315 }
std::size_t dofSizeWithoutGhosts() const
NumLib::MeshComponentMap _mesh_component_map

◆ findGlobalIndices()

template<typename ElementIterator >
void NumLib::LocalToGlobalIndexMap::findGlobalIndices ( ElementIterator  first,
ElementIterator  last,
std::vector< MeshLib::Node *> const &  nodes,
std::size_t const  mesh_id,
const int  comp_id,
const int  comp_id_write 
)
private

Definition at line 74 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, _mesh_subsets, _rows, NumLib::MeshComponentMap::getGlobalIndex(), and MeshLib::Node.

Referenced by LocalToGlobalIndexMap().

78 {
79  _rows.resize(std::distance(first, last), _mesh_subsets.size());
80 
81  std::unordered_set<MeshLib::Node*> const set_nodes(nodes.begin(), nodes.end());
82 
83  // For each element find the global indices for node/element
84  // components.
85  std::size_t elem_id = 0;
86  for (ElementIterator e = first; e != last; ++e, ++elem_id)
87  {
88  LineIndex indices;
89  indices.reserve((*e)->getNumberOfNodes());
90 
91  for (auto* n = (*e)->getNodes();
92  n < (*e)->getNodes() + (*e)->getNumberOfNodes(); ++n)
93  {
94  // Check if the element's node is in the given list of nodes.
95  if (set_nodes.find(*n) == set_nodes.end())
96  {
97  continue;
98  }
100  mesh_id, MeshLib::MeshItemType::Node, (*n)->getID());
101  auto const global_index =
103  if (global_index == std::numeric_limits<GlobalIndexType>::max())
104  {
105  continue;
106  }
107  indices.push_back(global_index);
108  }
109 
110  indices.shrink_to_fit();
111  _rows(elem_id, comp_id_write) = std::move(indices);
112  }
113 }
RowColumnIndices::LineIndex LineIndex
NumLib::MeshComponentMap _mesh_component_map
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.
GlobalIndexType getGlobalIndex(Location const &l, int const comp_id) const

◆ findGlobalIndicesWithElementID()

template<typename ElementIterator >
void NumLib::LocalToGlobalIndexMap::findGlobalIndicesWithElementID ( ElementIterator  first,
ElementIterator  last,
std::vector< MeshLib::Node *> const &  nodes,
std::size_t const  mesh_id,
const int  comp_id,
const int  comp_id_write 
)
private

Definition at line 41 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, _rows, NumLib::MeshComponentMap::getGlobalIndex(), and MeshLib::Node.

Referenced by LocalToGlobalIndexMap().

45 {
46  std::unordered_set<MeshLib::Node*> const set_nodes(nodes.begin(), nodes.end());
47 
48  // For each element find the global indices for node/element
49  // components.
50  for (ElementIterator e = first; e != last; ++e)
51  {
52  LineIndex indices;
53  indices.reserve((*e)->getNumberOfNodes());
54 
55  for (auto* n = (*e)->getNodes();
56  n < (*e)->getNodes()+(*e)->getNumberOfNodes(); ++n)
57  {
58  // Check if the element's node is in the given list of nodes.
59  if (set_nodes.find(*n) == set_nodes.end())
60  {
61  continue;
62  }
64  mesh_id, MeshLib::MeshItemType::Node, (*n)->getID());
65  indices.push_back(_mesh_component_map.getGlobalIndex(l, comp_id));
66  }
67 
68  indices.shrink_to_fit();
69  _rows((*e)->getID(), comp_id_write) = std::move(indices);
70  }
71 }
RowColumnIndices::LineIndex LineIndex
NumLib::MeshComponentMap _mesh_component_map
GlobalIndexType getGlobalIndex(Location const &l, int const comp_id) const

◆ getElementVariableIDs()

std::vector< int > NumLib::LocalToGlobalIndexMap::getElementVariableIDs ( std::size_t const  mesh_item_id) const

Definition at line 373 of file LocalToGlobalIndexMap.cpp.

References _rows, getGlobalComponent(), getNumberOfVariableComponents(), and getNumberOfVariables().

Referenced by ProcessLib::LIE::HydroMechanics::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), and ProcessLib::LIE::SmallDeformation::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()().

375 {
376  std::vector<int> vec;
377  for (int i = 0; i < getNumberOfVariables(); i++)
378  {
379  for (int j=0; j<getNumberOfVariableComponents(i); j++)
380  {
381  auto comp_id = getGlobalComponent(i, j);
382  if (!_rows(mesh_item_id, comp_id).empty())
383  {
384  vec.push_back(i);
385  }
386  }
387  }
388  std::sort(vec.begin(), vec.end());
389  vec.erase(std::unique(vec.begin(), vec.end()), vec.end());
390 
391  return vec;
392 }
int getGlobalComponent(int const variable_id, int const component_id) const
int getNumberOfVariableComponents(int variable_id) const

◆ getGhostIndices()

std::vector< GlobalIndexType > const & NumLib::LocalToGlobalIndexMap::getGhostIndices ( ) const

Get ghost indices, forwarded from MeshComponentMap.

Definition at line 417 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::getGhostIndices().

Referenced by ProcessLib::GlobalVectorFromNamedFunction::call(), and NumLib::LocalLinearLeastSquaresExtrapolator::extrapolate().

419 {
421 }
std::vector< GlobalIndexType > const & getGhostIndices() const
Get ghost indices (for DDC).
NumLib::MeshComponentMap _mesh_component_map

◆ getGlobalComponent()

int NumLib::LocalToGlobalIndexMap::getGlobalComponent ( int const  variable_id,
int const  component_id 
) const

◆ getGlobalIndex() [1/2]

GlobalIndexType NumLib::LocalToGlobalIndexMap::getGlobalIndex ( MeshLib::Location const &  l,
int const  variable_id,
int const  component_id 
) const

Definition at line 394 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, anonymous_namespace{Density100MPa.cpp}::c, getGlobalComponent(), and NumLib::MeshComponentMap::getGlobalIndex().

Referenced by ProcessLib::SourceTerms::Python::PythonSourceTermLocalAssembler< ShapeFunction, IntegrationMethod, GlobalDim >::assemble(), ProcessLib::PythonBoundaryCondition::getEssentialBCValues(), ProcessLib::getEssentialBCValuesLocal(), NumLib::getNodalValue(), NumLib::getNonGhostNodalValue(), ProcessLib::LIE::HydroMechanics::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::LIE::SmallDeformation::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::PhaseFieldIrreversibleDamageOracleBoundaryCondition::preTimestep(), and NumLib::transformVariableFromGlobalVector().

398 {
399  auto const c = getGlobalComponent(variable_id, component_id);
401 }
int getGlobalComponent(int const variable_id, int const component_id) const
NumLib::MeshComponentMap _mesh_component_map
GlobalIndexType getGlobalIndex(Location const &l, int const comp_id) const

◆ getGlobalIndex() [2/2]

GlobalIndexType NumLib::LocalToGlobalIndexMap::getGlobalIndex ( MeshLib::Location const &  l,
int const  global_component_id 
) const

Definition at line 403 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::getGlobalIndex().

405 {
406  return _mesh_component_map.getGlobalIndex(l, global_component_id);
407 }
NumLib::MeshComponentMap _mesh_component_map
GlobalIndexType getGlobalIndex(Location const &l, int const comp_id) const

◆ getGlobalIndices()

std::vector< GlobalIndexType > NumLib::LocalToGlobalIndexMap::getGlobalIndices ( const MeshLib::Location l) const

Forwards the respective method from MeshComponentMap.

Definition at line 410 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::getGlobalIndices().

Referenced by computeSparsityPatternNonPETSc().

412 {
414 }
NumLib::MeshComponentMap _mesh_component_map
std::vector< GlobalIndexType > getGlobalIndices(const Location &l) const

◆ getLocalIndex()

GlobalIndexType NumLib::LocalToGlobalIndexMap::getLocalIndex ( MeshLib::Location const &  l,
std::size_t const  comp_id,
std::size_t const  range_begin,
std::size_t const  range_end 
) const

Computes the index in a local (for DDC) vector for a given location and component; forwarded from MeshComponentMap.

Definition at line 425 of file LocalToGlobalIndexMap.cpp.

References _mesh_component_map, and NumLib::MeshComponentMap::getLocalIndex().

Referenced by ProcessLib::processOutputData().

428 {
429  return _mesh_component_map.getLocalIndex(l, comp_id, range_begin,
430  range_end);
431 }
GlobalIndexType getLocalIndex(Location const &l, int const comp_id, std::size_t const range_begin, std::size_t const range_end) const
NumLib::MeshComponentMap _mesh_component_map

◆ getMeshSubset() [1/2]

◆ getMeshSubset() [2/2]

MeshLib::MeshSubset const & NumLib::LocalToGlobalIndexMap::getMeshSubset ( int const  global_component_id) const

Definition at line 439 of file LocalToGlobalIndexMap.cpp.

References _mesh_subsets.

441 {
442  return _mesh_subsets[global_component_id];
443 }
std::vector< MeshLib::MeshSubset > _mesh_subsets
A vector of mesh subsets for each process variables&#39; components.

◆ getNumberOfComponents()

◆ getNumberOfElementComponents()

std::size_t NumLib::LocalToGlobalIndexMap::getNumberOfElementComponents ( std::size_t const  mesh_item_id) const

Definition at line 360 of file LocalToGlobalIndexMap.cpp.

References _rows, and anonymous_namespace{Density100MPa.cpp}::c.

Referenced by ProcessLib::LIE::SmallDeformation::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()().

361 {
362  std::size_t n = 0;
363  for (Table::Index c = 0; c < _rows.cols(); ++c)
364  {
365  if (!_rows(mesh_item_id, c).empty())
366  {
367  n++;
368  }
369  }
370  return n;
371 }

◆ getNumberOfElementDOF()

std::size_t NumLib::LocalToGlobalIndexMap::getNumberOfElementDOF ( std::size_t const  mesh_item_id) const

Definition at line 347 of file LocalToGlobalIndexMap.cpp.

References _rows, and anonymous_namespace{Density100MPa.cpp}::c.

Referenced by ProcessLib::ThermoHydroMechanics::LocalDataInitializer< LocalAssemblerInterface, ThermoHydroMechanicsLocalAssembler, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::RichardsMechanics::LocalDataInitializer< LocalAssemblerInterface, RichardsMechanicsLocalAssembler, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::HydroMechanics::LocalDataInitializer< LocalAssemblerInterface, HydroMechanicsLocalAssembler, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::LIE::HydroMechanics::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::LocalDataInitializer< LocalAssemblerInterface, SmallDeformationLocalAssembler, GlobalDim, ConstructorArgs >::operator()(), and ProcessLib::LIE::SmallDeformation::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()().

348 {
349  std::size_t ndof = 0;
350 
351  for (Table::Index c = 0; c < _rows.cols(); ++c)
352  {
353  ndof += _rows(mesh_item_id, c).size();
354  }
355 
356  return ndof;
357 }

◆ getNumberOfVariableComponents()

int NumLib::LocalToGlobalIndexMap::getNumberOfVariableComponents ( int  variable_id) const

Definition at line 322 of file LocalToGlobalIndexMap.cpp.

References _variable_component_offsets, and getNumberOfVariables().

Referenced by ProcessLib::SourceTerms::Python::PythonSourceTermLocalAssembler< ShapeFunction, IntegrationMethod, GlobalDim >::assemble(), ProcessLib::checkParametersOfDirichletBoundaryCondition(), ProcessLib::ConstraintDirichletBoundaryCondition::ConstraintDirichletBoundaryCondition(), ProcessLib::createPythonBoundaryCondition(), ProcessLib::createSourceTerm(), ProcessLib::GenericNaturalBoundaryCondition< BoundaryConditionData, LocalAssemblerImplementation >::GenericNaturalBoundaryCondition(), getElementVariableIDs(), ProcessLib::NormalTractionBoundaryCondition::NormalTractionBoundaryCondition< LocalAssemblerImplementation >::NormalTractionBoundaryCondition(), ProcessLib::LIE::HydroMechanics::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::LIE::SmallDeformation::LocalDataInitializer< LocalAssemblerInterface, LocalAssemblerDataMatrix, LocalAssemblerDataMatrixNearFracture, LocalAssemblerDataFracture, GlobalDim, ConstructorArgs >::operator()(), ProcessLib::PhaseFieldIrreversibleDamageOracleBoundaryCondition::PhaseFieldIrreversibleDamageOracleBoundaryCondition(), and NumLib::transformVariableFromGlobalVector().

323 {
324  assert(variable_id < getNumberOfVariables());
325  return _variable_component_offsets[variable_id + 1] -
326  _variable_component_offsets[variable_id];
327 }
std::vector< int > const _variable_component_offsets

◆ getNumberOfVariables()

◆ operator()()

LocalToGlobalIndexMap::RowColumnIndices NumLib::LocalToGlobalIndexMap::operator() ( std::size_t const  mesh_item_id,
const int  component_id 
) const

Definition at line 339 of file LocalToGlobalIndexMap.cpp.

References _columns, and _rows.

341 {
342  return RowColumnIndices(_rows(mesh_item_id, component_id),
343  _columns(mesh_item_id, component_id));
344 }
MathLib::RowColumnIndices< GlobalIndexType > RowColumnIndices

◆ size()

std::size_t NumLib::LocalToGlobalIndexMap::size ( ) const

Friends And Related Function Documentation

◆ operator<<

std::ostream& operator<< ( std::ostream &  os,
LocalToGlobalIndexMap const &  map 
)
friend

Prints first rows of the table, every line, and the mesh component map.

Definition at line 446 of file LocalToGlobalIndexMap.cpp.

447 {
448  std::size_t const max_lines = 10;
449  std::size_t lines_printed = 0;
450 
451  os << "Rows of the local to global index map; " << map._rows.size()
452  << " rows\n";
453  for (std::size_t e=0; e<map.size(); ++e)
454  {
455  os << "== e " << e << " ==\n";
456  for (int c = 0; c < map.getNumberOfComponents(); ++c)
457  {
458  auto const& line = map._rows(e, c);
459 
460  os << "c" << c << " { ";
461  std::copy(line.cbegin(), line.cend(),
462  std::ostream_iterator<std::size_t>(os, " "));
463  os << " }\n";
464  }
465 
466  if (lines_printed++ > max_lines)
467  {
468  os << "...\n";
469  break;
470  }
471  }
472 
473  os << "Mesh component map:\n" << map._mesh_component_map;
474  return os;
475 }
void copy(MatrixOrVector const &x, MatrixOrVector &y)
Copies x to y.
Definition: LinAlg.h:36

Member Data Documentation

◆ _columns

Table const& NumLib::LocalToGlobalIndexMap::_columns = _rows
private
See also
_rows

Definition at line 203 of file LocalToGlobalIndexMap.h.

Referenced by operator()().

◆ _mesh_component_map

◆ _mesh_subsets

std::vector<MeshLib::MeshSubset> NumLib::LocalToGlobalIndexMap::_mesh_subsets
private

A vector of mesh subsets for each process variables' components.

Definition at line 193 of file LocalToGlobalIndexMap.h.

Referenced by deriveBoundaryConstrainedMap(), findGlobalIndices(), getMeshSubset(), getNumberOfComponents(), and LocalToGlobalIndexMap().

◆ _rows

Table NumLib::LocalToGlobalIndexMap::_rows
private

Table contains for each element (first index) and each component (second index) a vector (LineIndex) of indices in the global stiffness matrix or vector

Definition at line 200 of file LocalToGlobalIndexMap.h.

Referenced by findGlobalIndices(), findGlobalIndicesWithElementID(), getElementVariableIDs(), getNumberOfElementComponents(), getNumberOfElementDOF(), LocalToGlobalIndexMap(), operator()(), NumLib::operator<<(), and size().

◆ _variable_component_offsets

std::vector<int> const NumLib::LocalToGlobalIndexMap::_variable_component_offsets
private

The documentation for this class was generated from the following files: