GridFormat 0.2.1
I/O-Library for grid-like data structures
Loading...
Searching...
No Matches
pvtr_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_PVTR_WRITER_HPP_
9#define GRIDFORMAT_VTK_PVTR_WRITER_HPP_
10
11#include <ostream>
12#include <string>
13#include <fstream>
14#include <algorithm>
15#include <filesystem>
16#include <array>
17#include <tuple>
18
20#include <gridformat/common/ranges.hpp>
21#include <gridformat/common/exceptions.hpp>
22#include <gridformat/common/lvalue_reference.hpp>
23
24#include <gridformat/parallel/communication.hpp>
25#include <gridformat/parallel/helpers.hpp>
26
27#include <gridformat/grid/grid.hpp>
28#include <gridformat/xml/element.hpp>
31
32namespace GridFormat {
33
38template<Concepts::RectilinearGrid Grid,
39 Concepts::Communicator Communicator>
40class PVTRWriter : public VTK::XMLWriterBase<Grid, PVTRWriter<Grid, Communicator>> {
42 using CT = CoordinateType<Grid>;
43
44 static constexpr std::size_t space_dim = 3;
45 static constexpr std::size_t dim = dimension<Grid>;
46 static constexpr int root_rank = 0;
47
48 public:
49 explicit PVTRWriter(LValueReferenceOf<const Grid> grid,
50 Communicator comm,
51 VTK::XMLOptions xml_opts = {})
52 : ParentType(grid.get(), ".pvtr", true, xml_opts)
53 , _comm(comm)
54 {}
55
56 const Communicator& communicator() const {
57 return _comm;
58 }
59
60 private:
61 Communicator _comm;
62
63 PVTRWriter _with(VTK::XMLOptions xml_opts) const override {
64 return PVTRWriter{this->grid(), _comm, std::move(xml_opts)};
65 }
66
67 void _write(std::ostream&) const override {
68 throw InvalidState(
69 "PVTRWriter does not support direct export into stream. "
70 "Use overload with filename instead!"
71 );
72 }
73
74 virtual void _write(const std::string& filename_with_ext) const override {
75 const auto& local_extents = extents(this->grid());
76 const auto [origin, is_negative_axis] = _get_origin_and_orientations();
77
79 const auto all_origins = Parallel::gather(_comm, origin, root_rank);
80 const auto all_extents = Parallel::gather(_comm, local_extents, root_rank);
81 const auto [exts_begin, exts_end, whole_extent, _] = helper.compute_extents_and_origin(
82 all_origins,
83 all_extents,
84 is_negative_axis
85 );
86
87 const auto my_whole_extent = Parallel::broadcast(_comm, whole_extent, root_rank);
88 const auto my_extent_offset = Parallel::scatter(_comm, Ranges::flat(exts_begin), root_rank);
89
90 _write_piece(filename_with_ext, Ranges::to_array<dim>(my_extent_offset), {my_whole_extent});
91 Parallel::barrier(_comm); // ensure all pieces finished successfully
92 if (Parallel::rank(_comm) == 0)
93 _write_pvtr_file(filename_with_ext, my_whole_extent, exts_begin, exts_end);
94 Parallel::barrier(_comm); // ensure .pvtr file is written before returning
95 }
96
97 auto _get_origin_and_orientations() const {
98 std::array<CT, dim> origin;
99 std::array<bool, dim> is_negative_axis;
100 for (unsigned dir = 0; dir < dim; ++dir) {
101 std::array<CT, 2> ordinates_01{0, 0};
102 const auto& dir_ordinates = ordinates(this->grid(), dir);
103 std::ranges::copy(std::views::take(dir_ordinates, 2), ordinates_01.begin());
104 origin[dir] = ordinates_01[0];
105 is_negative_axis[dir] = ordinates_01[1] - ordinates_01[0] < CT{0};
106 }
107 return std::make_tuple(origin, is_negative_axis);
108 }
109
110 void _write_piece(const std::string& par_filename,
111 const std::array<std::size_t, dim>& offset,
112 typename VTRWriter<Grid>::Domain domain) const {
113 auto writer = VTRWriter{this->grid(), this->_xml_opts}
114 .as_piece_for(std::move(domain))
115 .with_offset(offset);
116 this->copy_fields(writer);
117 writer.write(PVTK::piece_basefilename(par_filename, Parallel::rank(_comm)));
118 }
119
120 void _write_pvtr_file(const std::string& filename_with_ext,
121 const std::array<std::size_t, dim>& extents,
122 const std::vector<std::array<std::size_t, dim>>& proc_extents_begin,
123 const std::vector<std::array<std::size_t, dim>>& proc_extents_end) const {
124 std::ofstream file_stream(filename_with_ext, std::ios::out);
125
126 XMLElement pvtk_xml("VTKFile");
127 pvtk_xml.set_attribute("type", "PRectilinearGrid");
128
129 XMLElement& grid = pvtk_xml.add_child("PRectilinearGrid");
130 grid.set_attribute("WholeExtent", VTK::CommonDetail::extents_string(extents));
131
132 XMLElement& ppoint_data = grid.add_child("PPointData");
133 XMLElement& pcell_data = grid.add_child("PCellData");
134 XMLElement& pcoords = grid.add_child("PCoordinates");
135 std::visit([&] (const auto& encoder) {
136 std::visit([&] (const auto& data_format) {
137 PVTK::PDataArrayHelper pdata_helper{encoder, data_format, ppoint_data};
138 PVTK::PDataArrayHelper cdata_helper{encoder, data_format, pcell_data};
139 std::ranges::for_each(this->_point_field_names(), [&] (const std::string& name) {
140 pdata_helper.add(name, this->_get_point_field(name));
141 });
142 std::ranges::for_each(this->_cell_field_names(), [&] (const std::string& name) {
143 cdata_helper.add(name, this->_get_cell_field(name));
144 });
145
146 std::visit([&] <typename T> (const Precision<T>& prec) {
147 for (unsigned int i = 0; i < space_dim; ++i) {
148 XMLElement& pdata_array = pcoords.add_child("PDataArray");
149 pdata_array.set_attribute("NumberOfComponents", "1");
150 pdata_array.set_attribute("Name", "X_" + std::to_string(i));
151 pdata_array.set_attribute("format", VTK::data_format_name(encoder, data_format));
152 pdata_array.set_attribute("type", VTK::attribute_name(prec));
153 }
154 }, this->_xml_settings.coordinate_precision);
155 }, this->_xml_settings.data_format);
156 }, this->_xml_settings.encoder);
157
158 std::ranges::for_each(Parallel::ranks(_comm), [&] (int rank) {
159 std::array<std::size_t, dim> extents_begin;
160 std::array<std::size_t, dim> extents_end;
161 for (unsigned dir = 0; dir < dim; ++dir) {
162 extents_begin[dir] = proc_extents_begin[rank][dir];
163 extents_end[dir] = proc_extents_end[rank][dir];
164 }
165 auto& piece = grid.add_child("Piece");
166 piece.set_attribute("Extent", VTK::CommonDetail::extents_string(extents_begin, extents_end));
167 piece.set_attribute("Source", std::filesystem::path{
168 PVTK::piece_basefilename(filename_with_ext, rank) + ".vtr"
169 }.filename());
170 });
171
172 this->_set_default_active_fields(pvtk_xml.get_child("PRectilinearGrid"));
173 write_xml_with_version_header(pvtk_xml, file_stream, Indentation{{.width = 2}});
174 }
175};
176
177template<typename G, Concepts::Communicator C>
178PVTRWriter(G&&, const C&, VTK::XMLOptions = {}) -> PVTRWriter<std::remove_cvref_t<G>, C>;
179
180namespace Traits {
181
182template<typename... Args>
183struct WritesConnectivity<PVTRWriter<Args...>> : public std::false_type {};
184
185} // namespace Traits
186} // namespace GridFormat
187
188#endif // GRIDFORMAT_VTK_PVTR_WRITER_HPP_
Helper to add a PDataArray child to an xml element.
Definition: parallel.hpp:45
Writer for parallel .pvtr files.
Definition: pvtr_writer.hpp:40
Base class for VTK-XML Writer implementations.
Definition: xml.hpp:190
Writer for .vtr file format.
Definition: vtr_writer.hpp:35
Helper function for writing parallel VTK files.
Can be specialized by writers in case the file format does not contain connectivity information.
Definition: writer.hpp:36
Options for VTK-XML files for setting the desired encoding, data format and compression.
Definition: xml.hpp:99
Definition: vtr_writer.hpp:43