dune-grid  2.5.1
vtksequencewriterbase.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_GRID_IO_FILE_VTK_VTKSEQUENCEWRITERBASE_HH
5 #define DUNE_GRID_IO_FILE_VTK_VTKSEQUENCEWRITERBASE_HH
6 
7 #include <vector>
8 #include <iostream>
9 #include <sstream>
10 #include <fstream>
11 #include <iomanip>
12 #include <memory>
13 
15 #include <dune/common/path.hh>
16 
18 
19 namespace Dune {
20 
30  template<class GridView>
32  {
33  std::shared_ptr<VTKWriter<GridView> > vtkWriter_;
34  std::vector<double> timesteps_;
35  std::string name_,path_,extendpath_;
36  int rank_;
37  int size_;
38  public:
45  explicit VTKSequenceWriterBase( std::shared_ptr<VTKWriter<GridView> > vtkWriter,
46  const std::string& name,
47  const std::string& path,
48  const std::string& extendpath,
49  int rank,
50  int size)
51  : vtkWriter_(vtkWriter),
52  name_(name), path_(path),
53  extendpath_(extendpath),
54  rank_(rank),
55  size_(size)
56  {}
57 
59 
61  void addCellData (const std::shared_ptr<const typename VTKWriter<GridView>::VTKFunction> &p)
62  {
63  vtkWriter_->addCellData(p);
64  }
65 
71  template<class V >
72  void addCellData (const V &v, const std::string &name, int ncomps=1)
73  {
74  vtkWriter_->addCellData(v, name, ncomps);
75  }
76 
78  void addVertexData (const std::shared_ptr<const typename VTKWriter<GridView>::VTKFunction> &p)
79  {
80  vtkWriter_->addVertexData(p);
81  }
82 
88  template<class V >
89  void addVertexData (const V &v, const std::string &name, int ncomps=1)
90  {
91  vtkWriter_->addVertexData(v, name, ncomps);
92  }
93 
94 
100  void write (double time, VTK::OutputType type = VTK::ascii)
101  {
102  /* remember current time step */
103  unsigned int count = timesteps_.size();
104  timesteps_.push_back(time);
105 
106  /* write VTK file */
107  if(size_==1)
108  vtkWriter_->write(concatPaths(path_,seqName(count)),type);
109  else
110  vtkWriter_->pwrite(seqName(count), path_,extendpath_,type);
111 
112  /* write pvd file ... only on rank 0 */
113  if (rank_==0) {
114  std::ofstream pvdFile;
115  pvdFile.exceptions(std::ios_base::badbit | std::ios_base::failbit |
116  std::ios_base::eofbit);
117  std::string pvdname = name_ + ".pvd";
118  pvdFile.open(pvdname.c_str());
119  pvdFile << "<?xml version=\"1.0\"?> \n"
120  << "<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"" << VTK::getEndiannessString() << "\"> \n"
121  << "<Collection> \n";
122  for (unsigned int i=0; i<=count; i++)
123  {
124  // filename
125  std::string piecepath;
126  std::string fullname;
127  if(size_==1) {
128  piecepath = path_;
129  fullname = vtkWriter_->getSerialPieceName(seqName(i), piecepath);
130  }
131  else {
132  piecepath = concatPaths(path_, extendpath_);
133  fullname = vtkWriter_->getParallelHeaderName(seqName(i), piecepath, size_);
134  }
135  pvdFile << "<DataSet timestep=\"" << timesteps_[i]
136  << "\" group=\"\" part=\"0\" name=\"\" file=\""
137  << fullname << "\"/> \n";
138  }
139  pvdFile << "</Collection> \n"
140  << "</VTKFile> \n" << std::flush;
141  pvdFile.close();
142  }
143  }
144  private:
145 
146  // create sequence name
147  std::string seqName(unsigned int count) const
148  {
149  std::stringstream n;
150  n.fill('0');
151  n << name_ << "-" << std::setw(5) << count;
152  return n.str();
153  }
154  };
155 
156 } // end namespace Dune
157 
158 #endif
Writer for the ouput of grid functions in the vtk format.Writes arbitrary grid functions (living on c...
Definition: vtkwriter.hh:87
A base class for grid functions with any return type and dimension.
Definition: function.hh:38
void addCellData(const V &v, const std::string &name, int ncomps=1)
Adds a field of cell data to the VTK file.
Definition: vtksequencewriterbase.hh:72
void addCellData(const std::shared_ptr< const typename VTKWriter< GridView >::VTKFunction > &p)
Adds a field of cell data to the VTK file.
Definition: vtksequencewriterbase.hh:61
OutputType
How the bulk data should be stored in the file.
Definition: common.hh:40
Common stuff for the VTKWriter.
VTKSequenceWriterBase(std::shared_ptr< VTKWriter< GridView > > vtkWriter, const std::string &name, const std::string &path, const std::string &extendpath, int rank, int size)
Set up the VTKSequenceWriterBase class.
Definition: vtksequencewriterbase.hh:45
std::string getEndiannessString()
determine endianness of this C++ implementation
Definition: common.hh:278
void write(double time, VTK::OutputType type=VTK::ascii)
Writes VTK data for the given time,.
Definition: vtksequencewriterbase.hh:100
Base class to write pvd-files which contains a list of all collected vtk-files.
Definition: vtksequencewriterbase.hh:31
void addVertexData(const V &v, const std::string &name, int ncomps=1)
Adds a field of vertex data to the VTK file.
Definition: vtksequencewriterbase.hh:89
Include standard header files.
Definition: agrid.hh:59
Output to the file is in ascii.
Definition: common.hh:42
~VTKSequenceWriterBase()
Definition: vtksequencewriterbase.hh:58
Provides file i/o for the visualization toolkit.
void addVertexData(const std::shared_ptr< const typename VTKWriter< GridView >::VTKFunction > &p)
Adds a field of vertex data to the VTK file.
Definition: vtksequencewriterbase.hh:78