OGS
OGSMesh Class Reference

Detailed Description

Definition at line 20 of file OGSMesh.h.

#include <OGSMesh.h>

Collaboration diagram for OGSMesh:
[legend]

Public Member Functions

 OGSMesh (MeshLib::Mesh &mesh)
 
std::vector< double > getPointCoordinates () const
 
std::pair< std::vector< int >, std::vector< int > > getCells () const
 
std::vector< std::string > getDataArrayNames () const
 
MeshLib::MeshItemType meshItemType (std::string_view const name) const
 
template<typename T >
pybind11::array & dataArray (std::string const &name)
 
pybind11::object dataArray_dispatch (std::string const &name, std::string const &dtype)
 

Private Attributes

MeshLib::Mesh_mesh
 
std::map< std::string, pybind11::array > data_array_mapping
 

Constructor & Destructor Documentation

◆ OGSMesh()

OGSMesh::OGSMesh ( MeshLib::Mesh & mesh)
explicit

Definition at line 30 of file OGSMesh.cpp.

30: _mesh(mesh) {}
MeshLib::Mesh & _mesh
Definition OGSMesh.h:149

Member Function Documentation

◆ dataArray()

template<typename T >
pybind11::array & OGSMesh::dataArray ( std::string const & name)
inline

Definition at line 33 of file OGSMesh.h.

34 {
35 auto const data_array_it = data_array_mapping.find(name);
36 if (data_array_it != data_array_mapping.end())
37 {
38 INFO("found data array '{}' with address: {}", name,
39 fmt::ptr(&(data_array_it->second)));
40 return data_array_it->second;
41 }
42
43 auto& mesh_properties = _mesh.getProperties();
44
45 // check if PropertyVector with specified data type T exists
46 auto* property = mesh_properties.getPropertyVector<T>(name);
47 if (property == nullptr)
48 {
49 OGS_FATAL("Couldn't access data array '{}'.", name);
50 }
51
52 // Capsule ties lifetime of the numpy array to the PropertyVector
53 auto capsule = pybind11::capsule(property, [](void* /*ignored*/) {});
54
55 auto const n_components = property->getNumberOfGlobalComponents();
56 if (n_components == 1)
57 {
58 auto const& [it, success] = data_array_mapping.insert(
59 {name,
60 pybind11::array(
61 pybind11::buffer_info(
62 property->data(), // data pointer
63 sizeof(T), // item size
64 pybind11::format_descriptor<T>::format(), // dtype
65 // format
66 // string
67 1, // ndim
68 {property->size()}, // shape
69 {sizeof(T)} // strides
70 ),
71 capsule)});
72 if (!success)
73 {
74 OGS_FATAL("Could not insert data array '{}' into internal map.",
75 name);
76 }
77 INFO("insert data array '{}' with address: {}", name,
78 fmt::ptr(&(it->second)));
79 return it->second;
80 }
81 else
82 {
83 // 2D case: (n_items, n_components)
84 auto const& [it, success] = data_array_mapping.insert(
85 {name,
86 pybind11::array(
87 pybind11::buffer_info(
88 property->data(), // data pointer
89 sizeof(T), // item size
90 pybind11::format_descriptor<T>::format(), // dtype
91 // format
92 // string
93 2, // ndim
94 std::vector<pybind11::ssize_t>{
95 static_cast<pybind11::ssize_t>(property->size() /
96 n_components),
97 static_cast<pybind11::ssize_t>(
98 n_components)}, // shape
99 std::vector<pybind11::ssize_t>{
100 static_cast<pybind11::ssize_t>(sizeof(T) *
101 n_components),
102 static_cast<pybind11::ssize_t>(sizeof(T))}
103 // strides
104 ),
105 capsule)});
106 if (!success)
107 {
108 OGS_FATAL("Could not insert data array '{}' into internal map.",
109 name);
110 }
111 INFO("insert data array '{}' with address: {}", name,
112 fmt::ptr(&(it->second)));
113 return it->second;
114 }
115 }
#define OGS_FATAL(...)
Definition Error.h:26
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
Properties & getProperties()
Definition Mesh.h:136
PropertyVector< T > const * getPropertyVector(std::string_view name) const
std::map< std::string, pybind11::array > data_array_mapping
Definition OGSMesh.h:150

References INFO(), and OGS_FATAL.

◆ dataArray_dispatch()

pybind11::object OGSMesh::dataArray_dispatch ( std::string const & name,
std::string const & dtype )
inline

Definition at line 117 of file OGSMesh.h.

119 {
120 if (dtype == "double")
121 {
122 return dataArray<double>(name);
123 }
124 if (dtype == "float")
125 {
126 return dataArray<float>(name);
127 }
128 if (dtype == "int")
129 {
130 return dataArray<int>(name);
131 }
132 if (dtype == "int64")
133 {
134 return dataArray<int64_t>(name);
135 }
136 if (dtype == "std::size_t")
137 {
138 return dataArray<std::size_t>(name);
139 }
140 if (dtype == "char")
141 {
142 return dataArray<char>(name);
143 }
144
145 throw std::runtime_error("Unsupported dtype: " + dtype);
146 }
pybind11::array & dataArray(std::string const &name)
Definition OGSMesh.h:33

Referenced by PYBIND11_MODULE().

◆ getCells()

std::pair< std::vector< int >, std::vector< int > > OGSMesh::getCells ( ) const

Definition at line 38 of file OGSMesh.cpp.

39{
40 auto const& elements = _mesh.getElements();
41 std::vector<int> cells;
42 std::vector<int> cell_types;
43 for (auto const* element : elements)
44 {
45 ranges::copy(element->nodes() | MeshLib::views::ids,
46 std::back_inserter(cells));
47 cell_types.push_back(OGSToVtkCellType(element->getCellType()));
48 }
49 return {cells, cell_types};
50}
int OGSToVtkCellType(MeshLib::CellType ogs)
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:111
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
Definition Mesh.h:227

References _mesh, MeshLib::Mesh::getElements(), MeshLib::views::ids, and OGSToVtkCellType().

Referenced by PYBIND11_MODULE().

◆ getDataArrayNames()

std::vector< std::string > OGSMesh::getDataArrayNames ( ) const

Definition at line 52 of file OGSMesh.cpp.

53{
55}
std::vector< std::string > getPropertyVectorNames() const

References _mesh, MeshLib::Mesh::getProperties(), and MeshLib::Properties::getPropertyVectorNames().

Referenced by PYBIND11_MODULE().

◆ getPointCoordinates()

std::vector< double > OGSMesh::getPointCoordinates ( ) const

Definition at line 32 of file OGSMesh.cpp.

33{
34 return ranges::to<std::vector>(_mesh.getNodes() | MeshLib::views::coords |
35 ranges::views::join);
36}
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:108
constexpr ranges::views::view_closure coords
Definition Mesh.h:234

References _mesh, MeshLib::views::coords, and MeshLib::Mesh::getNodes().

Referenced by PYBIND11_MODULE().

◆ meshItemType()

MeshLib::MeshItemType OGSMesh::meshItemType ( std::string_view const name) const

Definition at line 57 of file OGSMesh.cpp.

58{
59 auto const& properties = _mesh.getProperties();
60
61 auto const& found = BaseLib::findElementOrError(
62 properties,
63 [&name](auto const& p) { return p.first == name; },
64 [&name]()
65 {
66 OGS_FATAL("A property with the name '{}' doesn't exist in the mesh",
67 name);
68 });
69 return found.second->getMeshItemType();
70}
ranges::range_reference_t< Range > findElementOrError(Range &range, std::predicate< ranges::range_reference_t< Range > > auto &&predicate, std::invocable auto error_callback)
Definition Algorithm.h:81

References _mesh, BaseLib::findElementOrError(), MeshLib::Mesh::getProperties(), and OGS_FATAL.

Referenced by PYBIND11_MODULE().

Member Data Documentation

◆ _mesh

MeshLib::Mesh& OGSMesh::_mesh
private

Definition at line 149 of file OGSMesh.h.

Referenced by getCells(), getDataArrayNames(), getPointCoordinates(), and meshItemType().

◆ data_array_mapping

std::map<std::string, pybind11::array> OGSMesh::data_array_mapping
private

Definition at line 150 of file OGSMesh.h.


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