OGS
DOFTableUtil.h
Go to the documentation of this file.
1
11#pragma once
12
13#include <Eigen/Core>
14
18
19namespace NumLib
20{
23double getNonGhostNodalValue(GlobalVector const& x, MeshLib::Mesh const& mesh,
24 NumLib::LocalToGlobalIndexMap const& dof_table,
25 std::size_t const node_id,
26 std::size_t const global_component_id);
27
30double getNodalValue(GlobalVector const& x, MeshLib::Mesh const& mesh,
31 NumLib::LocalToGlobalIndexMap const& dof_table,
32 std::size_t const node_id,
33 std::size_t const global_component_id);
34
37std::vector<GlobalIndexType> getIndices(
38 std::size_t const mesh_item_id,
39 NumLib::LocalToGlobalIndexMap const& dof_table);
40
44 std::size_t const id,
45 NumLib::LocalToGlobalIndexMap const& dof_table,
46 std::vector<GlobalIndexType>& indices);
47
51double norm(GlobalVector const& x, unsigned const global_component,
52 MathLib::VecNormType norm_type,
53 LocalToGlobalIndexMap const& dof_table, MeshLib::Mesh const& mesh);
54
60template <typename Functor>
62 GlobalVector const& input_vector, int const variable_id,
63 NumLib::LocalToGlobalIndexMap const& local_to_global_index_map,
64 MeshLib::PropertyVector<double>& output_vector, Functor map_function)
65{
67
68 // We fill the output with zeroes, because filling with NaN
69 // "breaks" visualization, e.g. of heat or mass flow rate
70 // with Taylor-Hood elements. I.e., all lower order
71 // properties would have NaN on the higher order nodes.
72 std::fill(begin(output_vector), end(output_vector), 0);
73
74 int const n_components =
75 local_to_global_index_map.getNumberOfVariableComponents(variable_id);
76 for (int component = 0; component < n_components; ++component)
77 {
78 auto const& mesh_subset =
79 local_to_global_index_map.getMeshSubset(variable_id, component);
80 auto const mesh_id = mesh_subset.getMeshID();
81 for (auto const& node : mesh_subset.getNodes())
82 {
83 auto const node_id = node->getID();
85 node_id);
86 auto const input_index = local_to_global_index_map.getGlobalIndex(
87 l, variable_id, component);
88 double const value = input_vector[input_index];
89
90 output_vector.getComponent(node_id, component) =
91 map_function(value);
92 }
93 }
94}
95
98template <typename Functor>
100 GlobalVector const& input_vector_on_bulk_mesh, int const variable_id,
101 NumLib::LocalToGlobalIndexMap const& local_to_global_index_map,
102 MeshLib::PropertyVector<double>& output_vector_on_submesh,
103 std::vector<std::size_t> const& map_submesh_node_id_to_bulk_mesh_node_id,
104 Functor map_function)
105{
106 MathLib::LinAlg::setLocalAccessibleVector(input_vector_on_bulk_mesh);
107
108 // We fill the output with zeroes, because filling with NaN
109 // "breaks" visualization, e.g. of heat or mass flow rate
110 // with Taylor-Hood elements. I.e., all lower order
111 // properties would have NaN on the higher order nodes.
112 std::fill(begin(output_vector_on_submesh), end(output_vector_on_submesh),
113 0);
114
115 int const n_components =
116 local_to_global_index_map.getNumberOfVariableComponents(variable_id);
117 for (int component = 0; component < n_components; ++component)
118 {
119 auto const& mesh_subset =
120 local_to_global_index_map.getMeshSubset(variable_id, component);
121 auto const mesh_id = mesh_subset.getMeshID();
122
123 for (std::size_t submesh_node_id = 0;
124 submesh_node_id < map_submesh_node_id_to_bulk_mesh_node_id.size();
125 ++submesh_node_id)
126 {
127 std::size_t const bulk_node_id =
128 map_submesh_node_id_to_bulk_mesh_node_id[submesh_node_id];
130 bulk_node_id);
131 auto const input_index = local_to_global_index_map.getGlobalIndex(
132 l, variable_id, component);
133
134 // On higher-order meshes some d.o.f.s might not be defined on all
135 // nodes. We silently ignore the d.o.f.s we cannot find.
136 [[unlikely]] if (input_index == NumLib::MeshComponentMap::nop)
137 {
138 continue;
139 }
140
141 double const value = input_vector_on_bulk_mesh[input_index];
142
143 output_vector_on_submesh.getComponent(submesh_node_id, component) =
144 map_function(value);
145 }
146 }
147}
148
149Eigen::VectorXd getLocalX(
150 std::size_t const mesh_item_id,
151 std::vector<NumLib::LocalToGlobalIndexMap const*> const& dof_tables,
152 std::vector<GlobalVector*> const& x);
153} // namespace NumLib
Global vector based on Eigen vector.
Definition EigenVector.h:25
std::size_t getMeshID() const
return this mesh ID
Definition MeshSubset.h:76
PROP_VAL_TYPE & getComponent(std::size_t tuple_index, int component)
Returns the value for the given component stored in the given tuple.
MathLib::RowColumnIndices< GlobalIndexType > RowColumnIndices
GlobalIndexType getGlobalIndex(MeshLib::Location const &l, int const variable_id, int const component_id) const
int getNumberOfVariableComponents(int variable_id) const
MeshLib::MeshSubset const & getMeshSubset(int const variable_id, int const component_id) const
static constexpr NUMLIB_EXPORT GlobalIndexType const nop
void setLocalAccessibleVector(PETScVector const &x)
Definition LinAlg.cpp:27
double getNodalValue(GlobalVector const &x, MeshLib::Mesh const &mesh, NumLib::LocalToGlobalIndexMap const &dof_table, std::size_t const node_id, std::size_t const global_component_id)
void transformVariableFromGlobalVector(GlobalVector const &input_vector, int const variable_id, NumLib::LocalToGlobalIndexMap const &local_to_global_index_map, MeshLib::PropertyVector< double > &output_vector, Functor map_function)
double norm(GlobalVector const &x, unsigned const global_component, MathLib::VecNormType norm_type, LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
double getNonGhostNodalValue(GlobalVector const &x, MeshLib::Mesh const &mesh, NumLib::LocalToGlobalIndexMap const &dof_table, std::size_t const node_id, std::size_t const global_component_id)
std::vector< GlobalIndexType > getIndices(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
NumLib::LocalToGlobalIndexMap::RowColumnIndices getRowColumnIndices(std::size_t const id, NumLib::LocalToGlobalIndexMap const &dof_table, std::vector< GlobalIndexType > &indices)
Eigen::VectorXd getLocalX(std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x)