GridFormat 0.2.1
I/O-Library for grid-like data structures
Loading...
Searching...
No Matches
pvtu_writer.hpp
Go to the documentation of this file.
1// SPDX-FileCopyrightText: 2022-2023 Dennis Gläser <dennis.glaeser@iws.uni-stuttgart.de>
2// SPDX-License-Identifier: MIT
8#ifndef GRIDFORMAT_VTK_PVTU_WRITER_HPP_
9#define GRIDFORMAT_VTK_PVTU_WRITER_HPP_
10
11#include <ostream>
12#include <string>
13#include <fstream>
14#include <filesystem>
15
16#include <gridformat/common/exceptions.hpp>
17#include <gridformat/common/lvalue_reference.hpp>
18#include <gridformat/parallel/communication.hpp>
19#include <gridformat/parallel/helpers.hpp>
20
21#include <gridformat/grid/grid.hpp>
22#include <gridformat/xml/element.hpp>
25
26namespace GridFormat {
27
32template<Concepts::UnstructuredGrid Grid,
33 Concepts::Communicator Communicator>
34class PVTUWriter : public VTK::XMLWriterBase<Grid, PVTUWriter<Grid, Communicator>> {
36
37 public:
38 explicit PVTUWriter(LValueReferenceOf<const Grid> grid,
39 Communicator comm,
40 VTK::XMLOptions xml_opts = {})
41 : ParentType(grid.get(), ".pvtu", false, xml_opts)
42 , _comm(comm)
43 {}
44
45 const Communicator& communicator() const {
46 return _comm;
47 }
48
49 private:
50 Communicator _comm;
51
52 PVTUWriter _with(VTK::XMLOptions xml_opts) const override {
53 return PVTUWriter{this->grid(), _comm, std::move(xml_opts)};
54 }
55
56 void _write(std::ostream&) const override {
57 throw InvalidState(
58 "PVTUWriter does not support direct export into stream. "
59 "Use overload with filename instead!"
60 );
61 }
62
63 virtual void _write(const std::string& filename_with_ext) const override {
64 _write_piece(filename_with_ext);
65 Parallel::barrier(_comm); // ensure all pieces finished successfully
66 if (Parallel::rank(_comm) == 0)
67 _write_pvtu_file(filename_with_ext);
68 Parallel::barrier(_comm); // ensure .pvtu file is written before returning
69 }
70
71 void _write_piece(const std::string& par_filename) const {
72 VTUWriter writer{this->grid(), this->_xml_opts};
73 this->copy_fields(writer);
74 writer.write(PVTK::piece_basefilename(par_filename, Parallel::rank(_comm)));
75 }
76
77 void _write_pvtu_file(const std::string& filename_with_ext) const {
78 std::ofstream file_stream(filename_with_ext, std::ios::out);
79
80 XMLElement pvtk_xml("VTKFile");
81 pvtk_xml.set_attribute("type", "PUnstructuredGrid");
82
83 XMLElement& grid = pvtk_xml.add_child("PUnstructuredGrid");
84 XMLElement& ppoint_data = grid.add_child("PPointData");
85 XMLElement& pcell_data = grid.add_child("PCellData");
86 std::visit([&] (const auto& encoder) {
87 std::visit([&] (const auto& data_format) {
88 PVTK::PDataArrayHelper pdata_helper{encoder, data_format, ppoint_data};
89 PVTK::PDataArrayHelper cdata_helper{encoder, data_format, pcell_data};
90 std::ranges::for_each(this->_point_field_names(), [&] (const std::string& name) {
91 pdata_helper.add(name, this->_get_point_field(name));
92 });
93 std::ranges::for_each(this->_cell_field_names(), [&] (const std::string& name) {
94 cdata_helper.add(name, this->_get_cell_field(name));
95 });
96 }, this->_xml_settings.data_format);
97 }, this->_xml_settings.encoder);
98
99 XMLElement& point_array = grid.add_child("PPoints").add_child("PDataArray");
100 point_array.set_attribute("NumberOfComponents", "3");
101 std::visit([&] <typename T> (const Precision<T>& prec) {
102 point_array.set_attribute("type", VTK::attribute_name(prec));
103 }, this->_xml_settings.coordinate_precision);
104
105 std::ranges::for_each(Parallel::ranks(_comm), [&] (int rank) {
106 grid.add_child("Piece").set_attribute("Source", std::filesystem::path{
107 PVTK::piece_basefilename(filename_with_ext, rank) + ".vtu"
108 }.filename());
109 });
110
111 this->_set_default_active_fields(pvtk_xml.get_child("PUnstructuredGrid"));
112 write_xml_with_version_header(pvtk_xml, file_stream, Indentation{{.width = 2}});
113 }
114};
115
116template<typename G, Concepts::Communicator C>
117PVTUWriter(G&&, const C&, VTK::XMLOptions = {}) -> PVTUWriter<std::remove_cvref_t<G>, C>;
118
119} // namespace GridFormat
120
121#endif // GRIDFORMAT_VTK_PVTU_WRITER_HPP_
Helper to add a PDataArray child to an xml element.
Definition: parallel.hpp:45
Writer for parallel .pvtu files.
Definition: pvtu_writer.hpp:34
Base class for VTK-XML Writer implementations.
Definition: xml.hpp:190
Writer for .vtu file format.
Definition: vtu_writer.hpp:30
Helper function for writing parallel VTK files.
Options for VTK-XML files for setting the desired encoding, data format and compression.
Definition: xml.hpp:99