vtkwriterinterface.hh 5.54 KB
Newer Older
1
2
3
4
#pragma once

#include <iosfwd>
#include <map>
5
6
#include <string>
#include <vector>
7

8
9
10
11
#ifdef HAVE_ZLIB
#include <zlib.h>
#endif

12
13
14
15
16
#include <dune/common/std/optional.hh>
#include <dune/vtk/filewriter.hh>
#include <dune/vtk/vtkfunction.hh>
#include <dune/vtk/vtktypes.hh>

17
namespace Dune
18
{
19
20
21
22
  // forward declaration
  template <class VtkWriter>
  class VtkTimeseriesWriter;

23
24
25
26
  template <class VtkWriter>
  class PvdWriter;

  /// Interface for file writers for the Vtk XML file formats
27
28
29
30
  template <class GridView, class DataCollector>
  class VtkWriterInterface
      : public FileWriter
  {
31
    template <class> friend class VtkTimeseriesWriter;
32
    template <class> friend class PvdWriter;
33

34
35
36
  protected:
    static constexpr int dimension = GridView::dimension;

37
    using VtkFunction = Dune::VtkFunction<GridView>;
38
39
40
41
42
43
44
45
46
    using pos_type = typename std::ostream::pos_type;

    enum PositionTypes {
      POINT_DATA,
      CELL_DATA
    };

  public:
    /// Constructor, stores the gridView
47
48
49
    VtkWriterInterface (GridView const& gridView,
                        Vtk::FormatTypes format = Vtk::BINARY,
                        Vtk::DataTypes datatype = Vtk::FLOAT32)
50
      : dataCollector_(gridView)
51
52
      , format_(format)
      , datatype_(datatype)
53
    {
54
55
56
57
58
59
#ifndef HAVE_ZLIB
      if (format_ == Vtk::COMPRESSED) {
        std::cout << "Dune is compiled without compression. Falling back to BINARY VTK output!\n";
        format_ = Vtk::BINARY;
      }
#endif
60
61
    }

62
    /// Write the attached data to the file
63
    virtual void write (std::string const& fn) const override;
64

Praetorius, Simon's avatar
Praetorius, Simon committed
65
66
67
    /// Attach point data to the writer, \see VtkFunction for possible arguments
    template <class Function, class... Args>
    VtkWriterInterface& addPointData (Function const& fct, Args&&... args)
68
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
69
      pointData_.emplace_back(fct, std::forward<Args>(args)...);
70
71
72
      return *this;
    }

Praetorius, Simon's avatar
Praetorius, Simon committed
73
74
75
    /// Attach cell data to the writer, \see VtkFunction for possible arguments
    template <class Function, class... Args>
    VtkWriterInterface& addCellData (Function const& fct, Args&&... args)
76
    {
Praetorius, Simon's avatar
Praetorius, Simon committed
77
      cellData_.emplace_back(fct, std::forward<Args>(args)...);
78
79
80
      return *this;
    }

81
  private:
82
    /// Write a serial VTK file in Unstructured format
83
    virtual void writeSerialFile (std::ofstream& out) const = 0;
84
85
86
87

    /// Write a parallel VTK file `pfilename.pvtu` in Unstructured format,
    /// with `size` the number of pieces and serial files given by `pfilename_p[i].vtu`
    /// for [i] in [0,...,size).
88
    virtual void writeParallelFile (std::ofstream& out, std::string const& pfilename, int size) const = 0;
89
90
91
92

    /// Return the file extension of the serial file (not including the dot)
    virtual std::string fileExtension () const = 0;

93
94
95
    /// Write points and cells in raw/compressed format to output stream
    virtual void writeGridAppended (std::ofstream& out, std::vector<std::uint64_t>& blocks) const = 0;

96
  protected:
97
98
99
100
101
    // Write the point or cell values given by the grid function `fct` to the
    // output stream `out`. In case of binary format, stores the streampos of XML
    // attributes "offset" in the vector `offsets`.
    void writeData (std::ofstream& out,
                    std::vector<pos_type>& offsets,
Praetorius, Simon's avatar
Praetorius, Simon committed
102
                    VtkFunction const& fct,
103
104
                    PositionTypes type,
                    Std::optional<std::size_t> timestep = {}) const;
105

106
107
    // Write points-data and cell-data in raw/compressed format to output stream
    void writeDataAppended (std::ofstream& out, std::vector<std::uint64_t>& blocks) const;
108

109
110
111
112
    // Write the coordinates of the vertices to the output stream `out`. In case
    // of binary format, stores the streampos of XML attributes "offset" in the
    // vector `offsets`.
    void writePoints (std::ofstream& out,
113
114
                      std::vector<pos_type>& offsets,
                      Std::optional<std::size_t> timestep = {}) const;
115

116
117
    // Write Appended section and fillin offset values to XML attributes
    void writeAppended (std::ofstream& out, std::vector<pos_type> const& offsets) const;
118
119
120
121

    // Write the `values` in blocks (possibly compressed) to the output
    // stream `out`. Return the written block size.
    template <class T>
122
    std::uint64_t writeValuesAppended (std::ofstream& out, std::vector<T> const& values) const;
123

124
    /// Return PointData/CellData attributes for the name of the first scalar/vector/tensor DataArray
125
    std::string getNames (std::vector<VtkFunction> const& data) const;
126
127
128
129
130
131
132
133

    // Returns endianness
    std::string getEndian () const
    {
      short i = 1;
      return (reinterpret_cast<char*>(&i)[1] == 1 ? "BigEndian" : "LittleEndian");
    }

134
    // provide accessor to \ref fileExtension virtual method
135
136
137
138
139
    std::string getFileExtension () const
    {
      return fileExtension();
    }

140
141
142
143
144
145
146
147
148
149
    Vtk::FormatTypes getFormat () const
    {
      return format_;
    }

    Vtk::DataTypes getDatatype () const
    {
      return datatype_;
    }

150
151
152
153
154
155
156
  protected:
    mutable DataCollector dataCollector_;

    Vtk::FormatTypes format_;
    Vtk::DataTypes datatype_;

    // attached data
Praetorius, Simon's avatar
Praetorius, Simon committed
157
158
    std::vector<VtkFunction> pointData_;
    std::vector<VtkFunction> cellData_;
159
160
161
162
163

    std::size_t const block_size = 1024*32;
    int compression_level = -1; // in [0,9], -1 ... use default value
  };

164
165
166
167
168
169
170
171

  template <class Writer>
  struct IsVtkWriter
  {
    template <class GV, class DC>
    static std::uint16_t test(VtkWriterInterface<GV,DC> const&);
    static std::uint8_t  test(...); // fall-back overload

172
    static constexpr bool value = sizeof(test(std::declval<Writer>())) > sizeof(std::uint8_t);
173
174
  };

175
} // end namespace Dune
176
177

#include "vtkwriterinterface.impl.hh"