dune-vtk 2.8
Loading...
Searching...
No Matches
vtkwriterinterface.hh
Go to the documentation of this file.
1#pragma once
2
3#include <iosfwd>
4#include <map>
5#include <memory>
6#include <optional>
7#include <string>
8#include <vector>
9
10#include <dune/common/parallel/mpihelper.hh>
12#include <dune/vtk/function.hh>
13#include <dune/vtk/types.hh>
14
15namespace Dune
16{
18
22 template <class GV, class DC>
24 : public Vtk::FileWriter
25 {
26 template <class> friend class TimeseriesWriter;
27 template <class> friend class PvdWriter;
28
29 public:
30 using GridView = GV;
31 using DataCollector = DC;
32
33 protected:
35 using pos_type = typename std::ostream::pos_type;
36
40 };
41
42 public:
44
57 VtkWriterInterface (GridView const& gridView,
61 : VtkWriterInterface(std::make_shared<DataCollector>(gridView), format, datatype, headertype)
62 {}
63
69 : VtkWriterInterface(stackobject_to_shared_ptr(dataCollector), format, datatype, headertype)
70 {}
71
73 VtkWriterInterface (std::shared_ptr<DataCollector> dataCollector,
77 : dataCollector_(std::move(dataCollector))
78 {
79 setFormat(format);
80 setDatatype(datatype);
81 setHeadertype(headertype);
82 }
83
84
86
92 virtual std::string write (std::string const& fn, std::optional<std::string> dir = {}) const override;
93
95
105 template <class Function, class... Args>
106 VtkWriterInterface& addPointData (Function&& fct, Args&&... args)
107 {
108 pointData_.emplace_back(std::forward<Function>(fct), std::forward<Args>(args)...,
110 return *this;
111 }
112
114
123 template <class Function, class... Args>
124 VtkWriterInterface& addCellData (Function&& fct, Args&&... args)
125 {
126 cellData_.emplace_back(std::forward<Function>(fct), std::forward<Args>(args)...,
128 return *this;
129 }
130
131
132 // Sets the VTK file format
134 {
135 format_ = format;
136
138#if HAVE_VTK_ZLIB
140#else
141 std::cout << "Dune is compiled without compression. Falling back to BINARY VTK output!\n";
143#endif
144 } else {
146 }
147
148 }
149
152 {
153 datatype_ = datatype;
154 }
155
158 {
159 headertype_ = datatype;
160 }
161
164 void setCompressor (Vtk::CompressorTypes compressor, int level = -1)
165 {
166 compressor_ = compressor;
167 compression_level = level;
168 VTK_ASSERT(level >= -1 && level <= 9);
169
172 }
173
174 private:
176 virtual void writeSerialFile (std::ofstream& out) const = 0;
177
181 virtual void writeParallelFile (std::ofstream& out, std::string const& pfilename, int size) const = 0;
182
184 virtual std::string fileExtension () const = 0;
185
187 virtual void writeGridAppended (std::ofstream& out, std::vector<std::uint64_t>& blocks) const = 0;
188
189 protected:
190 // Write the point or cell values given by the grid function `fct` to the
191 // output stream `out`. In case of binary format, append the streampos of XML
192 // attributes "offset" to the vector `offsets`.
193 void writeData (std::ofstream& out,
194 std::vector<pos_type>& offsets,
195 VtkFunction const& fct,
196 PositionTypes type,
197 std::optional<std::size_t> timestep = {}) const;
198
199 // Write point-data and cell-data in raw/compressed format to output stream
200 void writeDataAppended (std::ofstream& out, std::vector<std::uint64_t>& blocks) const;
201
202 // Write the coordinates of the vertices to the output stream `out`. In case
203 // of binary format, appends the streampos of XML attributes "offset" to the
204 // vector `offsets`.
205 void writePoints (std::ofstream& out,
206 std::vector<pos_type>& offsets,
207 std::optional<std::size_t> timestep = {}) const;
208
209 // Write Appended section and fillin offset values to XML attributes
210 void writeAppended (std::ofstream& out, std::vector<pos_type> const& offsets) const;
211
212 // Write the `values` in blocks (possibly compressed) to the output
213 // stream `out`. Return the written block size.
214 template <class HeaderType, class FloatType>
215 std::uint64_t writeValuesAppended (std::ofstream& out, std::vector<FloatType> const& values) const;
216
217 // Write the `values` in a space and newline separated list of ascii representations.
218 // The precision is controlled by the datatype and numerical_limits::digits10.
219 template <class T>
220 void writeValuesAscii (std::ofstream& out, std::vector<T> const& values) const;
221
222 // Write the XML file header of a VTK file `<VTKFile ...>`
223 void writeHeader (std::ofstream& out, std::string const& type) const;
224
226 std::string getNames (std::vector<VtkFunction> const& data) const;
227
228 // Returns endianness
229 std::string getEndian () const
230 {
231 short i = 1;
232 return (reinterpret_cast<char*>(&i)[1] == 1 ? "BigEndian" : "LittleEndian");
233 }
234
235 // provide accessor to \ref fileExtension virtual method
236 std::string getFileExtension () const
237 {
238 return fileExtension();
239 }
240
241 // Returns the VTK file format initialized in the constructor
243 {
244 return format_;
245 }
246
247 // Returns the global datatype used for coordinates and other global float values
249 {
250 return datatype_;
251 }
252
253 // Return the global MPI communicator.
254 auto comm () const
255 {
256 return MPIHelper::getCollectiveCommunication();
257 }
258
259 protected:
260 std::shared_ptr<DataCollector> dataCollector_;
261
266
267 // attached data
268 std::vector<VtkFunction> pointData_;
269 std::vector<VtkFunction> cellData_;
270
271 std::size_t const block_size = 1024*32;
272 int compression_level = -1; // in [0,9], -1 ... use default value
273 };
274
275
276 template <class Writer>
278 {
279 template <class GV, class DC>
280 static std::uint16_t test(VtkWriterInterface<GV,DC> const&);
281 static std::uint8_t test(...); // fall-back overload
282
283 static constexpr bool value = sizeof(test(std::declval<Writer>())) > sizeof(std::uint8_t);
284 };
285
286} // end namespace Dune
287
#define VTK_ASSERT(cond)
check if condition cond holds; otherwise, throw a VtkError.
Definition: errors.hh:29
Definition: writer.hh:13
FormatTypes
Type used for representing the output format.
Definition: types.hh:21
CompressorTypes
Definition: types.hh:149
DataTypes
Definition: types.hh:52
Definition: filewriter.hh:11
Wrapper class for functions allowing local evaluations.
Definition: function.hh:28
File-Writer for ParaView .pvd files.
Definition: pvdwriter.hh:18
Interface for file writers for the Vtk XML file formats.
Definition: vtkwriterinterface.hh:25
VtkWriterInterface(DataCollector &dataCollector, Vtk::FormatTypes format=Vtk::FormatTypes::BINARY, Vtk::DataTypes datatype=Vtk::DataTypes::FLOAT32, Vtk::DataTypes headertype=Vtk::DataTypes::UINT32)
Constructor, wraps the passed DataCollector in a non-destroying shared_ptr.
Definition: vtkwriterinterface.hh:65
friend class TimeseriesWriter
Definition: vtkwriterinterface.hh:26
std::shared_ptr< DataCollector > dataCollector_
Definition: vtkwriterinterface.hh:260
std::string getFileExtension() const
Definition: vtkwriterinterface.hh:236
Vtk::FormatTypes format_
Definition: vtkwriterinterface.hh:262
int compression_level
Definition: vtkwriterinterface.hh:272
Dune::Vtk::Function< GridView > VtkFunction
Definition: vtkwriterinterface.hh:34
std::vector< VtkFunction > pointData_
Definition: vtkwriterinterface.hh:268
void writeAppended(std::ofstream &out, std::vector< pos_type > const &offsets) const
Definition: vtkwriterinterface.impl.hh:156
std::vector< VtkFunction > cellData_
Definition: vtkwriterinterface.hh:269
void setFormat(Vtk::FormatTypes format)
Definition: vtkwriterinterface.hh:133
Vtk::CompressorTypes compressor_
Definition: vtkwriterinterface.hh:265
void writeValuesAscii(std::ofstream &out, std::vector< T > const &values) const
Definition: vtkwriterinterface.impl.hh:192
void setHeadertype(Vtk::DataTypes datatype)
Sets the integer type used in binary data headers.
Definition: vtkwriterinterface.hh:157
auto comm() const
Definition: vtkwriterinterface.hh:254
Vtk::DataTypes getDatatype() const
Definition: vtkwriterinterface.hh:248
GV GridView
Definition: vtkwriterinterface.hh:30
std::string getNames(std::vector< VtkFunction > const &data) const
Return PointData/CellData attributes for the name of the first scalar/vector/tensor DataArray.
Definition: vtkwriterinterface.impl.hh:353
VtkWriterInterface(std::shared_ptr< DataCollector > dataCollector, Vtk::FormatTypes format=Vtk::FormatTypes::BINARY, Vtk::DataTypes datatype=Vtk::DataTypes::FLOAT32, Vtk::DataTypes headertype=Vtk::DataTypes::UINT32)
Constructor, stores the passed DataCollector.
Definition: vtkwriterinterface.hh:73
VtkWriterInterface & addPointData(Function &&fct, Args &&... args)
Attach point data to the writer.
Definition: vtkwriterinterface.hh:106
Vtk::DataTypes headertype_
Definition: vtkwriterinterface.hh:264
std::string getEndian() const
Definition: vtkwriterinterface.hh:229
Vtk::DataTypes datatype_
Definition: vtkwriterinterface.hh:263
typename std::ostream::pos_type pos_type
Definition: vtkwriterinterface.hh:35
void setCompressor(Vtk::CompressorTypes compressor, int level=-1)
Definition: vtkwriterinterface.hh:164
void writeHeader(std::ofstream &out, std::string const &type) const
Definition: vtkwriterinterface.impl.hh:204
VtkWriterInterface & addCellData(Function &&fct, Args &&... args)
Attach cell data to the writer.
Definition: vtkwriterinterface.hh:124
void writeDataAppended(std::ofstream &out, std::vector< std::uint64_t > &blocks) const
Definition: vtkwriterinterface.impl.hh:132
void setDatatype(Vtk::DataTypes datatype)
Sets the global datatype used for coordinates and other global float values.
Definition: vtkwriterinterface.hh:151
DC DataCollector
Definition: vtkwriterinterface.hh:31
VtkWriterInterface(GridView const &gridView, Vtk::FormatTypes format=Vtk::FormatTypes::BINARY, Vtk::DataTypes datatype=Vtk::DataTypes::FLOAT32, Vtk::DataTypes headertype=Vtk::DataTypes::UINT32)
Constructor, passes the gridView to the DataCollector.
Definition: vtkwriterinterface.hh:57
void writeData(std::ofstream &out, std::vector< pos_type > &offsets, VtkFunction const &fct, PositionTypes type, std::optional< std::size_t > timestep={}) const
Definition: vtkwriterinterface.impl.hh:82
void writePoints(std::ofstream &out, std::vector< pos_type > &offsets, std::optional< std::size_t > timestep={}) const
Definition: vtkwriterinterface.impl.hh:109
virtual std::string write(std::string const &fn, std::optional< std::string > dir={}) const override
Write the attached data to the file.
Definition: vtkwriterinterface.impl.hh:28
Vtk::FormatTypes getFormat() const
Definition: vtkwriterinterface.hh:242
std::size_t const block_size
Definition: vtkwriterinterface.hh:271
std::uint64_t writeValuesAppended(std::ofstream &out, std::vector< FloatType > const &values) const
Definition: vtkwriterinterface.impl.hh:286
PositionTypes
Definition: vtkwriterinterface.hh:37
@ POINT_DATA
Definition: vtkwriterinterface.hh:38
@ CELL_DATA
Definition: vtkwriterinterface.hh:39
Definition: vtkwriterinterface.hh:278
static std::uint8_t test(...)
static std::uint16_t test(VtkWriterInterface< GV, DC > const &)
static constexpr bool value
Definition: vtkwriterinterface.hh:283