Skip to content
Snippets Groups Projects
Commit 5b651490 authored by Christoph Lehmann's avatar Christoph Lehmann
Browse files

[NL] misuse sparsity pattern for PETSc

parent 40cb1dd7
No related branches found
No related tags found
No related merge requests found
...@@ -12,20 +12,36 @@ ...@@ -12,20 +12,36 @@
#include "LocalToGlobalIndexMap.h" #include "LocalToGlobalIndexMap.h"
#include "MeshLib/NodeAdjacencyTable.h" #include "MeshLib/NodeAdjacencyTable.h"
namespace NumLib #ifdef USE_PETSC
#include "MeshLib/NodePartitionedMesh.h"
GlobalSparsityPattern computeSparsityPatternPETSc(
NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
MeshLib::Mesh const& mesh)
{ {
GlobalSparsityPattern computeSparsityPattern( assert(dynamic_cast<MeshLib::NodePartitionedMesh const*>(&mesh));
LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh) auto const& npmesh =
*static_cast<MeshLib::NodePartitionedMesh const*>(&mesh);
auto const max_nonzeroes = npmesh.getMaximumNConnectedNodesToNode();
// The sparsity pattern is misused here in the sense that it will only
// contain a single value.
return GlobalSparsityPattern(1, max_nonzeroes);
}
#else
GlobalSparsityPattern computeSparsityPatternNonPETSc(
NumLib::LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
{ {
MeshLib::NodeAdjacencyTable node_adjacency_table; MeshLib::NodeAdjacencyTable node_adjacency_table;
node_adjacency_table.createTable(mesh.getNodes()); node_adjacency_table.createTable(mesh.getNodes());
// A mapping mesh node id -> global indices // A mapping mesh node id -> global indices
// It acts as a cache for dof table queries. // It acts as a cache for dof table queries.
std::vector<std::vector<GlobalIndexType> > global_idcs; std::vector<std::vector<GlobalIndexType>> global_idcs;
global_idcs.reserve(mesh.getNNodes()); global_idcs.reserve(mesh.getNNodes());
for (std::size_t n=0; n<mesh.getNNodes(); ++n) for (std::size_t n = 0; n < mesh.getNNodes(); ++n)
{ {
MeshLib::Location l(mesh.getID(), MeshLib::MeshItemType::Node, n); MeshLib::Location l(mesh.getID(), MeshLib::MeshItemType::Node, n);
global_idcs.push_back(dof_table.getGlobalIndices(l)); global_idcs.push_back(dof_table.getGlobalIndices(l));
...@@ -34,15 +50,18 @@ GlobalSparsityPattern computeSparsityPattern( ...@@ -34,15 +50,18 @@ GlobalSparsityPattern computeSparsityPattern(
GlobalSparsityPattern sparsity_pattern(dof_table.dofSizeWithGhosts()); GlobalSparsityPattern sparsity_pattern(dof_table.dofSizeWithGhosts());
// Map adjacent mesh nodes to "adjacent global indices". // Map adjacent mesh nodes to "adjacent global indices".
for (std::size_t n=0; n<mesh.getNNodes(); ++n) for (std::size_t n = 0; n < mesh.getNNodes(); ++n)
{ {
auto const& node_ids = node_adjacency_table.getAdjacentNodes(n); auto const& node_ids = node_adjacency_table.getAdjacentNodes(n);
for (auto an : node_ids) { for (auto an : node_ids)
{
auto const& row_ids = global_idcs[an]; auto const& row_ids = global_idcs[an];
auto const num_components = row_ids.size(); auto const num_components = row_ids.size();
for (auto r : row_ids) { for (auto r : row_ids)
{
// Each component leads to an entry in the row. // Each component leads to an entry in the row.
// For the sparsity pattern only the number of entries are needed. // For the sparsity pattern only the number of entries are
// needed.
sparsity_pattern[r] += num_components; sparsity_pattern[r] += num_components;
} }
} }
...@@ -50,5 +69,18 @@ GlobalSparsityPattern computeSparsityPattern( ...@@ -50,5 +69,18 @@ GlobalSparsityPattern computeSparsityPattern(
return sparsity_pattern; return sparsity_pattern;
} }
#endif
namespace NumLib
{
GlobalSparsityPattern computeSparsityPattern(
LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh)
{
#ifdef USE_PETSC
return computeSparsityPatternPETSc(dof_table, mesh);
#else
return computeSparsityPatternNonPETSc(dof_table, mesh);
#endif
}
} }
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment