27 #ifndef EWOMS_VTK_MULTI_WRITER_HH
28 #define EWOMS_VTK_MULTI_WRITER_HH
36 #include <opm/material/common/Valgrind.hpp>
38 #include <dune/common/fvector.hh>
39 #include <dune/istl/bvector.hh>
40 #include <dune/grid/io/file/vtk/vtkwriter.hh>
60 template <
class Gr
idView,
int vtkFormat>
63 enum { dim = GridView::dimension };
65 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGVertexLayout> VertexMapper;
66 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGElementLayout> ElementMapper;
79 const std::string &simName =
"",
80 std::string multiFileName =
"")
82 , elementMapper_(gridView)
83 , vertexMapper_(gridView)
85 simName_ = (simName.empty()) ?
"sim" : simName;
86 multiFileName_ = multiFileName;
87 if (multiFileName_.empty()) {
88 multiFileName_ = simName_;
89 multiFileName_ +=
".pvd";
94 commRank_ = gridView.comm().rank();
95 commSize_ = gridView.comm().size();
110 {
return curWriterNum_; }
121 elementMapper_.update();
122 vertexMapper_.update();
130 if (!multiFile_.is_open()) {
131 startMultiFile_(multiFileName_);
134 curWriter_ =
new VtkWriter(gridView_, Dune::VTK::conforming);
138 curOutFileName_ = fileName_();
150 managedScalarBuffers_.push_back(buf);
163 for (
int i = 0; i < numOuter; ++ i)
164 (*buf)[i].resize(numInner);
166 managedVectorBuffers_.push_back(buf);
188 sanitizeScalarBuffer_(buf);
190 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
192 FunctionPtr fnPtr(
new VtkFn(name,
197 curWriter_->addVertexData(fnPtr);
217 sanitizeScalarBuffer_(buf);
219 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
221 FunctionPtr fnPtr(
new VtkFn(name,
226 curWriter_->addCellData(fnPtr);
247 sanitizeVectorBuffer_(buf);
249 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
251 FunctionPtr fnPtr(
new VtkFn(name,
256 curWriter_->addVertexData(fnPtr);
264 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
267 for (
size_t colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
268 std::ostringstream oss;
269 oss << name <<
"[" << colIdx <<
"]";
271 FunctionPtr fnPtr(
new VtkFn(oss.str(),
277 curWriter_->addVertexData(fnPtr);
298 sanitizeVectorBuffer_(buf);
300 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
302 FunctionPtr fnPtr(
new VtkFn(name,
307 curWriter_->addCellData(fnPtr);
315 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
318 for (
size_t colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
319 std::ostringstream oss;
320 oss << name <<
"[" << colIdx <<
"]";
322 FunctionPtr fnPtr(
new VtkFn(oss.str(),
328 curWriter_->addCellData(fnPtr);
342 std::string fileName;
344 fileName = curWriter_->write(curOutFileName_.c_str(),
345 static_cast<Dune::VTK::OutputType
>(vtkFormat));
349 multiFile_.precision(16);
350 multiFile_ <<
" <DataSet timestep=\"" << curTime_ <<
"\" file=\""
351 << fileName <<
"\"/>\n";
358 while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
359 delete managedScalarBuffers_.front();
360 managedScalarBuffers_.pop_front();
362 while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
363 delete managedVectorBuffers_.front();
364 managedVectorBuffers_.pop_front();
376 template <
class Restarter>
379 res.serializeSectionBegin(
"VTKMultiWriter");
380 res.serializeStream() << curWriterNum_ <<
"\n";
382 if (commRank_ == 0) {
385 if (multiFile_.is_open()) {
387 filePos = multiFile_.tellp();
388 multiFile_.seekp(0, std::ios::end);
389 fileLen = multiFile_.tellp();
390 multiFile_.seekp(filePos);
393 res.serializeStream() << fileLen <<
" " << filePos <<
"\n";
396 std::ifstream multiFileIn(multiFileName_.c_str());
397 char *tmp =
new char[fileLen];
398 multiFileIn.read(tmp, static_cast<long>(fileLen));
399 res.serializeStream().write(tmp, fileLen);
404 res.serializeSectionEnd();
410 template <
class Restarter>
413 res.deserializeSectionBegin(
"VTKMultiWriter");
414 res.deserializeStream() >> curWriterNum_;
416 if (commRank_ == 0) {
418 std::getline(res.deserializeStream(), dummy);
421 size_t filePos, fileLen;
422 res.deserializeStream() >> fileLen >> filePos;
423 std::getline(res.deserializeStream(), dummy);
424 if (multiFile_.is_open())
428 multiFile_.open(multiFileName_.c_str());
430 char *tmp =
new char[fileLen];
431 res.deserializeStream().read(tmp, fileLen);
432 multiFile_.write(tmp, fileLen);
436 multiFile_.seekp(filePos);
440 std::getline(res.deserializeStream(), tmp);
442 res.deserializeSectionEnd();
446 std::string fileName_()
449 std::ostringstream oss;
450 oss << simName_ <<
"-" << std::setw(5) << std::setfill(
'0')
455 std::string fileSuffix_()
456 {
return (GridView::dimension == 1) ?
"vtp" :
"vtu"; }
458 void startMultiFile_(
const std::string &multiFileName)
461 if (commRank_ == 0) {
463 multiFile_.open(multiFileName.c_str());
464 multiFile_ <<
"<?xml version=\"1.0\"?>\n"
465 "<VTKFile type=\"Collection\"\n"
467 " byte_order=\"LittleEndian\"\n"
468 " compressor=\"vtkZLibDataCompressor\">\n"
473 void finishMultiFile_()
476 if (commRank_ == 0) {
478 std::ofstream::pos_type pos = multiFile_.tellp();
479 multiFile_ <<
" </Collection>\n"
481 multiFile_.seekp(pos);
488 void sanitizeScalarBuffer_(ScalarBuffer &b)
490 for (
unsigned j = 0; j < b.size(); ++j) {
491 Valgrind::CheckDefined(b[j]);
495 if (std::abs(b[j]) < std::numeric_limits<float>::min()) {
501 void sanitizeVectorBuffer_(VectorBuffer &b)
503 size_t nComps = b[0].size();
504 for (
size_t i = 0; i < b.size(); ++i) {
505 for (
size_t j = 0; j < nComps; ++j) {
506 Valgrind::CheckDefined(b[i][j]);
510 if (std::abs(b[i][j]) < std::numeric_limits<float>::min()) {
517 const GridView gridView_;
518 ElementMapper elementMapper_;
519 VertexMapper vertexMapper_;
521 std::string simName_;
522 std::ofstream multiFile_;
523 std::string multiFileName_;
528 VtkWriter *curWriter_;
530 std::string curOutFileName_;
533 std::list<ScalarBuffer *> managedScalarBuffers_;
534 std::list<VectorBuffer *> managedVectorBuffers_;
VectorBuffer * allocateManagedVectorBuffer(int numOuter, int numInner)
Allocate a managed buffer for a vector field.
Definition: vtkmultiwriter.hh:160
void attachVectorVertexData(VectorBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:245
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: vtkmultiwriter.hh:72
std::vector< Vector > VectorBuffer
Definition: baseoutputwriter.hh:48
std::vector< Scalar > ScalarBuffer
Definition: baseoutputwriter.hh:47
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: vtkmultiwriter.hh:73
std::vector< Tensor > TensorBuffer
Definition: baseoutputwriter.hh:49
void serialize(Restarter &res)
Write the multi-writer's state to a restart file.
Definition: vtkmultiwriter.hh:377
Dune::VTKWriter< GridView > VtkWriter
Definition: vtkmultiwriter.hh:76
void attachScalarVertexData(ScalarBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:186
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition: vtkmultiwriter.hh:119
double Scalar
Definition: baseoutputwriter.hh:44
The base class for all output writers.
Definition: baseoutputwriter.hh:41
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkscalarfunction.hh:49
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
BaseOutputWriter::Tensor Tensor
Definition: vtkmultiwriter.hh:71
void attachVectorElementData(VectorBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:296
Definition: baseauxiliarymodule.hh:35
ScalarBuffer * allocateManagedScalarBuffer(int numEntities)
Allocate a managed buffer for a scalar field.
Definition: vtkmultiwriter.hh:147
Provides a vector-valued function using Dune::FieldVectors as elements.
The base class for all output writers.
void beginWrite(double t)
Called whenever a new time step must be written.
Definition: vtkmultiwriter.hh:128
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
Definition: vtktensorfunction.hh:47
void attachScalarElementData(ScalarBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:215
void attachTensorVertexData(TensorBuffer &buf, std::string name)
Add a finished vertex-centered tensor field to the output.
Definition: vtkmultiwriter.hh:262
void attachTensorElementData(TensorBuffer &buf, std::string name)
Add a finished element-centered tensor field to the output.
Definition: vtkmultiwriter.hh:313
BaseOutputWriter::Scalar Scalar
Definition: vtkmultiwriter.hh:69
int curWriterNum() const
Returns the number of the current VTK file.
Definition: vtkmultiwriter.hh:109
~VtkMultiWriter()
Definition: vtkmultiwriter.hh:98
Dune::DynamicVector< double > Vector
Definition: baseoutputwriter.hh:45
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkvectorfunction.hh:49
BaseOutputWriter::Vector Vector
Definition: vtkmultiwriter.hh:70
void deserialize(Restarter &res)
Read the multi-writer's state from a restart file.
Definition: vtkmultiwriter.hh:411
Dune::DynamicMatrix< double > Tensor
Definition: baseoutputwriter.hh:46
Provides a vector-valued function using Dune::FieldVectors as elements.
VtkMultiWriter(const GridView &gridView, const std::string &simName="", std::string multiFileName="")
Definition: vtkmultiwriter.hh:78
void endWrite(bool onlyDiscard=false)
Finalizes the current writer.
Definition: vtkmultiwriter.hh:339
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:61
BaseOutputWriter::TensorBuffer TensorBuffer
Definition: vtkmultiwriter.hh:74