Loading [MathJax]/extensions/tex2jax.js
OGS
transformData.cpp
Go to the documentation of this file.
1
10#include "transformData.h"
11
12#include <algorithm>
13#include <array>
14#include <optional>
15#include <range/v3/action/push_back.hpp>
16#include <range/v3/view/transform.hpp>
17#include <string>
18
19#include "BaseLib/cpp23.h"
20#if defined(USE_PETSC)
21#include "BaseLib/MPI.h"
22#endif
23#include "InfoLib/GitInfo.h"
25#include "MeshLib/Mesh.h"
26#include "MeshLib/MeshEnums.h"
27#include "MeshLib/Node.h"
29#include "partition.h"
30
31using namespace BaseLib;
32namespace MeshLib::IO
33{
35{
36 // https://xdmf.org/index.php/XDMF_Model_and_Format#Topology, Section
37 // Arbitrary
39 unsigned int number_of_nodes;
40};
41
42static constexpr auto elemOGSTypeToXDMFType()
43{
45 elem_type{};
73 ParentDataType::WEDGE, 6}; // parallel triangle wedge
82 return elem_type;
83}
84
86
87constexpr auto cellTypeOGS2XDMF(MeshLib::CellType const& cell_type)
88{
89 return elem_type_ogs2xdmf[to_underlying(cell_type)];
90}
91
92std::optional<XdmfHdfData> transformAttribute(
93 std::pair<std::string, PropertyVectorBase*> const& property_pair,
94 unsigned int const n_files, unsigned int const chunk_size_bytes)
95{
96 // 3 data that will be captured and written by lambda f below
98 std::size_t num_of_tuples = 0;
99 void const* data_ptr = 0;
100
101 // lambda f : Collects properties from the PropertyVectorBase. It captures
102 // (and overwrites) data that can only be collected via the typed property.
103 // It has boolean return type to allow kind of pipe using || operator.
104 auto f = [&data_type, &num_of_tuples, &data_ptr,
105 &property_pair](auto basic_type) -> bool
106 {
107 auto const property_base = property_pair.second;
108 auto const typed_property =
109 dynamic_cast<PropertyVector<decltype(basic_type)> const*>(
110 property_base);
111 if (typed_property == nullptr)
112 {
113 return false;
114 }
115 // overwrite captured data
116 num_of_tuples = typed_property->getNumberOfTuples();
117 data_ptr = typed_property->data();
118
119 if constexpr (std::is_same_v<double, decltype(basic_type)>)
120 {
121 // The standard 64-bit IEEE 754 floating-point type
122 // (double-precision) has a 53 bit fractional part (52 bits written,
123 // one implied)
124 static_assert((std::numeric_limits<double>::digits == 53),
125 "Double has 52 bits fractional part");
127 }
128 else if constexpr (std::is_same_v<float, decltype(basic_type)>)
129 {
130 // The standard 32-bit IEEE 754 floating-point type
131 // (single-precision) has a 24 bit fractional part (23 bits written,
132 // one implied)
133 static_assert((std::numeric_limits<float>::digits == 24),
134 "Float has 23 bits fractional part");
136 }
137 else if constexpr (std::is_same_v<int32_t, decltype(basic_type)>)
138 {
139 data_type = MeshPropertyDataType::int32;
140 }
141 else if constexpr (std::is_same_v<uint32_t, decltype(basic_type)>)
142 {
144 }
145 else if constexpr (std::is_same_v<int64_t, decltype(basic_type)>)
146 {
147 data_type = MeshPropertyDataType::int64;
148 }
149 else if constexpr (std::is_same_v<uint64_t, decltype(basic_type)>)
150 {
152 }
153 else if constexpr (std::is_same_v<int8_t, decltype(basic_type)>)
154 {
155 data_type = MeshPropertyDataType::int8;
156 }
157 else if constexpr (std::is_same_v<uint8_t, decltype(basic_type)>)
158 {
159 data_type = MeshPropertyDataType::uint8;
160 }
161 else if constexpr (std::is_same_v<char, decltype(basic_type)>)
162 {
164 }
165 else if constexpr (std::is_same_v<unsigned char, decltype(basic_type)>)
166 {
167 static_assert((std::numeric_limits<unsigned char>::digits == 8),
168 "Unsigned char has 8 bits");
169 data_type = MeshPropertyDataType::uchar;
170 }
171 else if constexpr (std::is_same_v<unsigned long, decltype(basic_type)>)
172 {
173 if (sizeof(unsigned long) == 8 &&
174 std::numeric_limits<unsigned long>::digits == 64)
175 {
177 }
178 else if (sizeof(unsigned long) == 4 &&
179 std::numeric_limits<unsigned long>::digits == 32)
180 {
182 }
183 else
184 {
185 return false;
186 }
187 }
188 else if constexpr (std::is_same_v<std::size_t, decltype(basic_type)>)
189 {
190 if (sizeof(std::size_t) == 8 &&
191 std::numeric_limits<std::size_t>::digits == 64)
192 {
194 }
195 else if (sizeof(std::size_t) == 4 &&
196 std::numeric_limits<std::size_t>::digits == 32)
197 {
199 }
200 else
201 {
202 return false;
203 }
204 }
205 else
206 {
207 return false;
208 }
209 return true;
210 };
211
212 f(double{}) || f(float{}) || f(int{}) || f(long{}) || f(unsigned{}) ||
213 f(long{}) || f(static_cast<unsigned long>(0)) || f(std::size_t{}) ||
214 f(char{}) || f(static_cast<unsigned char>(0));
215
216 if (data_type == MeshPropertyDataType::unknown)
217 {
218 return std::nullopt;
219 }
220
221 auto const& property_base = property_pair.second;
222 auto const& global_components =
223 property_base->getNumberOfGlobalComponents();
224 // TODO (tm) property_pair vector::getNumberOfGlobalComponents should return
225 // unsigned value. Then explicit cast from signed to unsigned int and
226 // assert can be removed here. Implicit cast to long long is fine and
227 // can be kept
228 assert(global_components >= 0);
229 auto const ui_global_components =
230 static_cast<unsigned int>(global_components);
231
232 MeshLib::MeshItemType const mesh_item_type =
233 property_base->getMeshItemType();
234
235 std::string const& name = property_base->getPropertyName();
236
237 HdfData hdf = {data_ptr, num_of_tuples, ui_global_components, name,
238 data_type, n_files, chunk_size_bytes};
239
240 XdmfData xdmf{num_of_tuples, ui_global_components, data_type,
241 name, mesh_item_type, 0,
242 n_files, std::nullopt};
243
244 return XdmfHdfData{std::move(hdf), std::move(xdmf)};
245}
246
247std::vector<XdmfHdfData> transformAttributes(
248 MeshLib::Mesh const& mesh, unsigned int const n_files,
249 unsigned int const chunk_size_bytes)
250{
251 MeshLib::Properties const& properties = mesh.getProperties();
252
253 // \TODO (tm) use c++20 ranges
254 // a = p | filter (first!=OGS_VERSION) | filter null_opt | transformAttr |
255 std::vector<XdmfHdfData> attributes;
256 for (auto const& [name, property_base] : properties)
257 {
259 {
260 continue;
261 }
262
263 if (!property_base->is_for_output)
264 {
265 continue;
266 }
267
268 if (auto const attribute = transformAttribute(
269
270 std::pair(std::string(name), property_base), n_files,
271 chunk_size_bytes))
272
273 {
274 attributes.push_back(attribute.value());
275 }
276 else
277 {
278 WARN("Could not create attribute meta of {:s}.", name);
279 }
280 }
281 return attributes;
282}
283
284std::vector<double> transformToXDMFGeometry(MeshLib::Mesh const& mesh)
285{
286 std::vector<MeshLib::Node*> const& nodes = mesh.getNodes();
287
288 int const point_size = 3;
289 std::vector<double> values;
290 values.reserve(nodes.size() * point_size);
291 for (auto const& n : nodes)
292 {
293 const double* x = n->data();
294 values.insert(values.cend(), x, x + point_size);
295 }
296
297 return values;
298}
299
300XdmfHdfData transformGeometry(MeshLib::Mesh const& mesh, double const* data_ptr,
301 unsigned int const n_files,
302 unsigned int const chunk_size_bytes)
303{
304 std::string const name = "geometry";
305 std::vector<MeshLib::Node*> const& nodes = mesh.getNodes();
306
307 int const point_size = 3;
308 auto const& partition_dim = nodes.size();
309
310 HdfData const hdf = {data_ptr,
311 partition_dim,
312 point_size,
313 name,
315 n_files,
316 chunk_size_bytes};
317 XdmfData const xdmf = {
318 partition_dim, point_size, MeshPropertyDataType::float64,
319 name, std::nullopt, 2,
320 n_files, std::nullopt};
321
322 return XdmfHdfData{std::move(hdf), std::move(xdmf)};
323}
324
326{
327 auto const& elements = mesh.getElements();
328
329 if (elements.empty())
330 {
332 }
333 auto const ogs_cell_type = elements[0]->getCellType();
334
335 if (std::any_of(elements.begin(), elements.end(),
336 [&](auto const& cell)
337 { return cell->getCellType() != ogs_cell_type; }))
338 {
340 }
341
342 return cellTypeOGS2XDMF(ogs_cell_type).id;
343}
344
345#if !defined(USE_PETSC)
347{
348 return getLocalTopologyType(mesh);
349}
350#else
352{
353 auto const local_topology_type = getLocalTopologyType(mesh);
354 // get the topology_type from other partitions
355 BaseLib::MPI::Mpi const mpi;
356
357 std::vector<int> const topology_types =
358 BaseLib::MPI::allgather(static_cast<int>(local_topology_type), mpi);
359
360 auto const common_topology_type =
361 std::all_of(topology_types.begin(), topology_types.end(),
362 [&local_topology_type](int other)
363 { return other == static_cast<int>(local_topology_type); })
364 ? local_topology_type
366 return common_topology_type;
367}
368#endif // USE_PETSC
369
370std::pair<std::vector<std::size_t>, ParentDataType> transformToXDMFTopology(
371 MeshLib::Mesh const& mesh, std::size_t const offset)
372{
373 std::vector<MeshLib::Element*> const& elements = mesh.getElements();
374 std::vector<std::size_t> values;
375
376 auto const push_cellnode_ids_to_vector =
377 [&values, &offset](auto const& cell)
378 {
379 values |= ranges::actions::push_back(
380 cell->nodes() | MeshLib::views::ids |
381 ranges::views::transform([&offset](auto const node_id)
382 { return node_id + offset; }));
383 };
384
385 auto const topology_type = getTopologyType(mesh);
386 if (topology_type == ParentDataType::MIXED)
387 {
388 values.reserve(elements.size() * 2); // each cell has at least two
389 // numbers
390 for (auto const& cell : elements)
391 {
392 auto const ogs_cell_type = cell->getCellType();
393 auto const xdmf_cell_id = cellTypeOGS2XDMF(ogs_cell_type).id;
394 values.push_back(to_underlying(xdmf_cell_id));
395 if (ogs_cell_type == MeshLib::CellType::LINE2)
396 {
397 values.push_back(2);
398 }
399 // LINE3 maps to ParentDataType::EDGE_3 -> no special care required
400 push_cellnode_ids_to_vector(cell);
401 }
402 }
403 else
404 {
405 values.reserve(elements.size() * elements[0]->getNumberOfNodes());
406 for (auto const& cell : elements)
407 {
408 push_cellnode_ids_to_vector(cell);
409 }
410 }
411
412 return {values, topology_type};
413}
414
415XdmfHdfData transformTopology(std::vector<std::size_t> const& values,
416 ParentDataType const parent_data_type,
417 unsigned int const n_files,
418 unsigned int const chunk_size_bytes)
419{
420 std::string const name = "topology";
421 auto const tuple_size = ParentDataType2String(parent_data_type).second;
422 HdfData const hdf = {values.data(),
423 values.size() / tuple_size,
424 tuple_size,
425 name,
427 n_files,
428 chunk_size_bytes};
429 XdmfData const xdmf{values.size() / tuple_size,
430 tuple_size,
432 name,
433 std::nullopt,
434 3,
435 n_files,
436 parent_data_type};
437
438 return XdmfHdfData{std::move(hdf), std::move(xdmf)};
439}
440} // namespace MeshLib::IO
Definition of the Element class.
Git information.
void WARN(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:40
Definition of mesh-related Enumerations.
std::pair< std::string, std::size_t > ParentDataType2String(ParentDataType p)
Enum for all propertyVector data types and XDMF ParentDataTypes.
MeshPropertyDataType
Definition of the Mesh class.
Definition of the Node class.
std::vector< Node * > const & getNodes() const
Get the nodes-vector for the mesh.
Definition Mesh.h:108
std::vector< Element * > const & getElements() const
Get the element-vector for the mesh.
Definition Mesh.h:111
Properties & getProperties()
Definition Mesh.h:136
Property manager on mesh items. Class Properties manages scalar, vector or matrix properties....
Definition Properties.h:33
std::size_t getNumberOfTuples() const
static std::vector< T > allgather(T const &value, Mpi const &mpi)
Definition MPI.h:104
constexpr auto to_underlying(E e) noexcept
Converts an enumeration to its underlying type.
Definition cpp23.h:29
const std::string OGS_VERSION
Definition GitInfo.cpp:20
ParentDataType getLocalTopologyType(MeshLib::Mesh const &mesh)
constexpr auto cellTypeOGS2XDMF(MeshLib::CellType const &cell_type)
ParentDataType getTopologyType(MeshLib::Mesh const &mesh)
std::optional< XdmfHdfData > transformAttribute(std::pair< std::string, PropertyVectorBase * > const &property_pair, unsigned int const n_files, unsigned int const chunk_size_bytes)
std::pair< std::vector< std::size_t >, ParentDataType > transformToXDMFTopology(MeshLib::Mesh const &mesh, std::size_t const offset)
Copies all cells into a new vector. Contiguous data used for writing. The topology is specific to xdm...
std::vector< double > transformToXDMFGeometry(MeshLib::Mesh const &mesh)
Copies all node points into a new vector. Contiguous data used for writing. Conform with XDMF standar...
XdmfHdfData transformGeometry(MeshLib::Mesh const &mesh, double const *data_ptr, unsigned int const n_files, unsigned int const chunk_size_bytes)
Create meta data for geometry used for hdf5 and xdmf.
static constexpr auto elemOGSTypeToXDMFType()
std::vector< XdmfHdfData > transformAttributes(MeshLib::Mesh const &mesh, unsigned int const n_files, unsigned int const chunk_size_bytes)
Create meta data for attributes used for hdf5 and xdmf.
XdmfHdfData transformTopology(std::vector< std::size_t > const &values, ParentDataType const parent_data_type, unsigned int const n_files, unsigned int const chunk_size_bytes)
Create meta data for topology used for HDF5 and XDMF.
constexpr auto elem_type_ogs2xdmf
constexpr ranges::views::view_closure ids
For an element of a range view return its id.
Definition Mesh.h:227
CellType
Types of mesh elements supported by OpenGeoSys.
Definition MeshEnums.h:64
Dispatches functions specific to execution platform (w/o MPI)
Transforms OGS Mesh into vectorized data.