OGS
ComputeSparsityPattern.cpp
Go to the documentation of this file.
1
12
13#include <numeric>
14#include <range/v3/range/conversion.hpp>
15#include <range/v3/view/transform.hpp>
16
18#include "MeshLib/Location.h"
20
21#ifdef USE_PETSC
23
25 NumLib::LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
26{
27 assert(dynamic_cast<MeshLib::NodePartitionedMesh const*>(&mesh));
28 auto const& npmesh =
29 *static_cast<MeshLib::NodePartitionedMesh const*>(&mesh);
30
31 auto const max_nonzeroes = dof_table.getNumberOfGlobalComponents() *
32 npmesh.getMaximumNConnectedNodesToNode();
33
34 // The sparsity pattern is misused here in the sense that it will only
35 // contain a single value.
36 return GlobalSparsityPattern(1, max_nonzeroes);
37}
38#else
39GlobalSparsityPattern computeSparsityPatternNonPETSc(
40 NumLib::LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
41{
42 MeshLib::NodeAdjacencyTable const node_adjacency_table(mesh);
43
44 // A mapping mesh node id -> global indices
45 // It acts as a cache for dof table queries.
46 auto const global_idcs =
48 ranges::views::transform([&](auto&& l)
49 { return dof_table.getGlobalIndices(l); }) |
50 ranges::to<std::vector>();
51
52 GlobalSparsityPattern sparsity_pattern(dof_table.dofSizeWithGhosts());
53
54 // Map adjacent mesh nodes to "adjacent global indices".
55 for (std::size_t n = 0; n < mesh.getNumberOfNodes(); ++n)
56 {
57 auto const& an = node_adjacency_table.getAdjacentNodes(n);
58 auto const n_self_dof = global_idcs[n].size();
59 auto const n_connected_dof =
60 std::accumulate(cbegin(an), cend(an), 0,
61 [&](auto const result, auto const i)
62 { return result + global_idcs[i].size(); });
63 auto const n_dof = n_self_dof + n_connected_dof;
64 for (auto global_index : global_idcs[n])
65 {
66 sparsity_pattern[global_index] = n_dof;
67 }
68 }
69
70 return sparsity_pattern;
71}
72#endif
73
74namespace NumLib
75{
77 LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
78{
79#ifdef USE_PETSC
80 return computeSparsityPatternPETSc(dof_table, mesh);
81#else
82 return computeSparsityPatternNonPETSc(dof_table, mesh);
83#endif
84}
85
86} // namespace NumLib
GlobalSparsityPattern computeSparsityPatternPETSc(NumLib::LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
MathLib::SparsityPattern< GlobalIndexType > GlobalSparsityPattern
Definition of mesh class for partitioned mesh (by node) for parallel computing within the framework o...
std::size_t getNumberOfNodes() const
Get the number of nodes.
Definition Mesh.h:102
std::vector< GlobalIndexType > getGlobalIndices(const MeshLib::Location &l) const
Forwards the respective method from MeshComponentMap.
auto meshLocations(Mesh const &mesh, MeshItemType const item_type)
Definition Mesh.h:238
GlobalSparsityPattern computeSparsityPattern(LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
Computes a sparsity pattern for the given inputs.