vtkmultiwriter.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
28 #ifndef EWOMS_VTK_MULTI_WRITER_HH
29 #define EWOMS_VTK_MULTI_WRITER_HH
30 
31 #include "vtkscalarfunction.hh"
32 #include "vtkvectorfunction.hh"
33 #include "vtktensorfunction.hh"
34 
36 
37 #include <opm/common/Valgrind.hpp>
38 #include <opm/common/Unused.hpp>
39 
40 #include <dune/common/fvector.hh>
41 #include <dune/istl/bvector.hh>
42 #include <dune/grid/io/file/vtk/vtkwriter.hh>
43 
44 #if HAVE_MPI
45 #include <mpi.h>
46 #endif
47 
48 #include <list>
49 #include <string>
50 #include <limits>
51 #include <sstream>
52 #include <fstream>
53 
54 namespace Ewoms {
62 template <class GridView, int vtkFormat>
64 {
65  enum { dim = GridView::dimension };
66 
67  typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGVertexLayout> VertexMapper;
68  typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGElementLayout> ElementMapper;
69 
70 public:
71  typedef BaseOutputWriter::Scalar Scalar;
72  typedef BaseOutputWriter::Vector Vector;
73  typedef BaseOutputWriter::Tensor Tensor;
74  typedef BaseOutputWriter::ScalarBuffer ScalarBuffer;
75  typedef BaseOutputWriter::VectorBuffer VectorBuffer;
76  typedef BaseOutputWriter::TensorBuffer TensorBuffer;
77 
78  typedef Dune::VTKWriter<GridView> VtkWriter;
79 #if DUNE_VERSION_NEWER(DUNE_GRID, 2,5)
80  typedef std::shared_ptr< Dune::VTKFunction< GridView > > FunctionPtr;
81 #else
82  typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
83 #endif
84 
85  VtkMultiWriter(const GridView& gridView,
86  const std::string& simName = "",
87  std::string multiFileName = "")
88  : gridView_(gridView)
89  , elementMapper_(gridView)
90  , vertexMapper_(gridView)
91  {
92  simName_ = (simName.empty()) ? "sim" : simName;
93  multiFileName_ = multiFileName;
94  if (multiFileName_.empty()) {
95  multiFileName_ = simName_;
96  multiFileName_ += ".pvd";
97  }
98 
99  curWriterNum_ = 0;
100 
101  commRank_ = gridView.comm().rank();
102  commSize_ = gridView.comm().size();
103  }
104 
105  ~VtkMultiWriter()
106  {
107  finishMultiFile_();
108 
109  if (commRank_ == 0)
110  multiFile_.close();
111  }
112 
116  int curWriterNum() const
117  { return curWriterNum_; }
118 
126  void gridChanged()
127  {
128  elementMapper_.update();
129  vertexMapper_.update();
130  }
131 
135  void beginWrite(double t)
136  {
137  if (!multiFile_.is_open()) {
138  startMultiFile_(multiFileName_);
139  }
140 
141  curTime_ = t;
142  curOutFileName_ = fileName_();
143 
144  curWriter_ = new VtkWriter(gridView_, Dune::VTK::conforming);
145  ++curWriterNum_;
146  }
147 
154  ScalarBuffer *allocateManagedScalarBuffer(size_t numEntities)
155  {
156  ScalarBuffer *buf = new ScalarBuffer(numEntities);
157  managedScalarBuffers_.push_back(buf);
158  return buf;
159  }
160 
167  VectorBuffer *allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
168  {
169  VectorBuffer *buf = new VectorBuffer(numOuter);
170  for (size_t i = 0; i < numOuter; ++ i)
171  (*buf)[i].resize(numInner);
172 
173  managedVectorBuffers_.push_back(buf);
174  return buf;
175  }
176 
193  void attachScalarVertexData(ScalarBuffer& buf, std::string name)
194  {
195  sanitizeScalarBuffer_(buf);
196 
198  FunctionPtr fnPtr(new VtkFn(name,
199  gridView_,
200  vertexMapper_,
201  buf,
202  /*codim=*/dim));
203  curWriter_->addVertexData(fnPtr);
204  }
205 
221  void attachScalarElementData(ScalarBuffer& buf, std::string name)
222  {
223  sanitizeScalarBuffer_(buf);
224 
226  FunctionPtr fnPtr(new VtkFn(name,
227  gridView_,
228  elementMapper_,
229  buf,
230  /*codim=*/0));
231  curWriter_->addCellData(fnPtr);
232  }
233 
250  void attachVectorVertexData(VectorBuffer& buf, std::string name)
251  {
252  sanitizeVectorBuffer_(buf);
253 
255  FunctionPtr fnPtr(new VtkFn(name,
256  gridView_,
257  vertexMapper_,
258  buf,
259  /*codim=*/dim));
260  curWriter_->addVertexData(fnPtr);
261  }
262 
266  void attachTensorVertexData(TensorBuffer& buf, std::string name)
267  {
269 
270  for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
271  std::ostringstream oss;
272  oss << name << "[" << colIdx << "]";
273 
274  FunctionPtr fnPtr(new VtkFn(oss.str(),
275  gridView_,
276  vertexMapper_,
277  buf,
278  /*codim=*/dim,
279  colIdx));
280  curWriter_->addVertexData(fnPtr);
281  }
282  }
283 
299  void attachVectorElementData(VectorBuffer& buf, std::string name)
300  {
301  sanitizeVectorBuffer_(buf);
302 
304  FunctionPtr fnPtr(new VtkFn(name,
305  gridView_,
306  elementMapper_,
307  buf,
308  /*codim=*/0));
309  curWriter_->addCellData(fnPtr);
310  }
311 
315  void attachTensorElementData(TensorBuffer& buf, std::string name)
316  {
318 
319  for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
320  std::ostringstream oss;
321  oss << name << "[" << colIdx << "]";
322 
323  FunctionPtr fnPtr(new VtkFn(oss.str(),
324  gridView_,
325  elementMapper_,
326  buf,
327  /*codim=*/0,
328  colIdx));
329  curWriter_->addCellData(fnPtr);
330  }
331  }
332 
340  void endWrite(bool onlyDiscard = false)
341  {
342  if (!onlyDiscard) {
343  std::string fileName;
344  // write the actual data as vtu or vtp (plus the pieces file in the parallel case)
345  fileName = curWriter_->write(/*name=*/curOutFileName_.c_str(),
346  static_cast<Dune::VTK::OutputType>(vtkFormat));
347 
348  // determine name to write into the multi-file for the
349  // current time step
350  multiFile_.precision(16);
351  multiFile_ << " <DataSet timestep=\"" << curTime_ << "\" file=\""
352  << fileName << "\"/>\n";
353  }
354  else
355  --curWriterNum_;
356 
357  // discard managed objects and the current VTK writer
358  delete curWriter_;
359  while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
360  delete managedScalarBuffers_.front();
361  managedScalarBuffers_.pop_front();
362  }
363  while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
364  delete managedVectorBuffers_.front();
365  managedVectorBuffers_.pop_front();
366  }
367 
368  // temporarily write the closing XML mumbo-jumbo to the mashup
369  // file so that the data set can be loaded even if the
370  // simulation is aborted (or not yet finished)
371  finishMultiFile_();
372  }
373 
377  template <class Restarter>
378  void serialize(Restarter& res)
379  {
380  res.serializeSectionBegin("VTKMultiWriter");
381  res.serializeStream() << curWriterNum_ << "\n";
382 
383  if (commRank_ == 0) {
384  std::streamsize fileLen = 0;
385  std::streamoff filePos = 0;
386  if (multiFile_.is_open()) {
387  // write the meta file into the restart file
388  filePos = multiFile_.tellp();
389  multiFile_.seekp(0, std::ios::end);
390  fileLen = multiFile_.tellp();
391  multiFile_.seekp(filePos);
392  }
393 
394  res.serializeStream() << fileLen << " " << filePos << "\n";
395 
396  if (fileLen > 0) {
397  std::ifstream multiFileIn(multiFileName_.c_str());
398  char *tmp = new char[fileLen];
399  multiFileIn.read(tmp, static_cast<long>(fileLen));
400  res.serializeStream().write(tmp, fileLen);
401  delete[] tmp;
402  }
403  }
404 
405  res.serializeSectionEnd();
406  }
407 
411  template <class Restarter>
412  void deserialize(Restarter& res)
413  {
414  res.deserializeSectionBegin("VTKMultiWriter");
415  res.deserializeStream() >> curWriterNum_;
416 
417  if (commRank_ == 0) {
418  std::string dummy;
419  std::getline(res.deserializeStream(), dummy);
420 
421  // recreate the meta file from the restart file
422  std::streamoff filePos;
423  std::streamsize fileLen;
424  res.deserializeStream() >> fileLen >> filePos;
425  std::getline(res.deserializeStream(), dummy);
426  if (multiFile_.is_open())
427  multiFile_.close();
428 
429  if (fileLen > 0) {
430  multiFile_.open(multiFileName_.c_str());
431 
432  char *tmp = new char[fileLen];
433  res.deserializeStream().read(tmp, fileLen);
434  multiFile_.write(tmp, fileLen);
435  delete[] tmp;
436  }
437 
438  multiFile_.seekp(filePos);
439  }
440  else {
441  std::string tmp;
442  std::getline(res.deserializeStream(), tmp);
443  }
444  res.deserializeSectionEnd();
445  }
446 
447 private:
448  std::string fileName_()
449  {
450  // use a new file name for each time step
451  std::ostringstream oss;
452  oss << simName_ << "-" << std::setw(5) << std::setfill('0')
453  << curWriterNum_;
454  return oss.str();
455  }
456 
457  std::string fileSuffix_()
458  { return (GridView::dimension == 1) ? "vtp" : "vtu"; }
459 
460  void startMultiFile_(const std::string& multiFileName)
461  {
462  // only the first process writes to the multi-file
463  if (commRank_ == 0) {
464  // generate one meta vtk-file holding the individual time steps
465  multiFile_.open(multiFileName.c_str());
466  multiFile_ << "<?xml version=\"1.0\"?>\n"
467  "<VTKFile type=\"Collection\"\n"
468  " version=\"0.1\"\n"
469  " byte_order=\"LittleEndian\"\n"
470  " compressor=\"vtkZLibDataCompressor\">\n"
471  " <Collection>\n";
472  }
473  }
474 
475  void finishMultiFile_()
476  {
477  // only the first process writes to the multi-file
478  if (commRank_ == 0) {
479  // make sure that we always have a working meta file
480  std::ofstream::pos_type pos = multiFile_.tellp();
481  multiFile_ << " </Collection>\n"
482  "</VTKFile>\n";
483  multiFile_.seekp(pos);
484  multiFile_.flush();
485  }
486  }
487 
488  // make sure the field is well defined if running under valgrind
489  // and make sure that all values can be displayed by paraview
490  void sanitizeScalarBuffer_(ScalarBuffer& b OPM_UNUSED)
491  {
492  // nothing to do: this is done by VtkScalarFunction
493  }
494 
495  void sanitizeVectorBuffer_(VectorBuffer& b OPM_UNUSED)
496  {
497  // nothing to do: this is done by VtkVectorFunction
498  }
499 
500  const GridView gridView_;
501  ElementMapper elementMapper_;
502  VertexMapper vertexMapper_;
503 
504  std::string simName_;
505  std::ofstream multiFile_;
506  std::string multiFileName_;
507 
508  int commSize_; // number of processes in the communicator
509  int commRank_; // rank of the current process in the communicator
510 
511  VtkWriter *curWriter_;
512  double curTime_;
513  std::string curOutFileName_;
514  int curWriterNum_;
515 
516  std::list<ScalarBuffer *> managedScalarBuffers_;
517  std::list<VectorBuffer *> managedVectorBuffers_;
518 };
519 } // namespace Ewoms
520 
521 #endif
The base class for all output writers.
Definition: baseoutputwriter.hh:43
void serialize(Restarter &res)
Write the multi-writer&#39;s state to a restart file.
Definition: vtkmultiwriter.hh:378
void beginWrite(double t)
Called whenever a new time step must be written.
Definition: vtkmultiwriter.hh:135
ScalarBuffer * allocateManagedScalarBuffer(size_t numEntities)
Allocate a managed buffer for a scalar field.
Definition: vtkmultiwriter.hh:154
Definition: baseauxiliarymodule.hh:37
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
Provides a vector-valued function using Dune::FieldVectors as elements.
VectorBuffer * allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
Allocate a managed buffer for a vector field.
Definition: vtkmultiwriter.hh:167
void attachScalarVertexData(ScalarBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:193
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
void endWrite(bool onlyDiscard=false)
Finalizes the current writer.
Definition: vtkmultiwriter.hh:340
void attachVectorElementData(VectorBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:299
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkvectorfunction.hh:51
Provides a vector-valued function using Dune::FieldVectors as elements.
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkscalarfunction.hh:53
void deserialize(Restarter &res)
Read the multi-writer&#39;s state from a restart file.
Definition: vtkmultiwriter.hh:412
void attachVectorVertexData(VectorBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:250
The base class for all output writers.
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
Definition: vtktensorfunction.hh:49
int curWriterNum() const
Returns the number of the current VTK file.
Definition: vtkmultiwriter.hh:116
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition: vtkmultiwriter.hh:126
void attachTensorVertexData(TensorBuffer &buf, std::string name)
Add a finished vertex-centered tensor field to the output.
Definition: vtkmultiwriter.hh:266
void attachTensorElementData(TensorBuffer &buf, std::string name)
Add a finished element-centered tensor field to the output.
Definition: vtkmultiwriter.hh:315
void attachScalarElementData(ScalarBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:221