28 #ifndef EWOMS_VTK_MULTI_WRITER_HH 29 #define EWOMS_VTK_MULTI_WRITER_HH 37 #include <opm/common/Valgrind.hpp> 38 #include <opm/common/Unused.hpp> 40 #include <dune/common/fvector.hh> 41 #include <dune/istl/bvector.hh> 42 #include <dune/grid/io/file/vtk/vtkwriter.hh> 62 template <
class Gr
idView,
int vtkFormat>
65 enum { dim = GridView::dimension };
67 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGVertexLayout> VertexMapper;
68 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView, Dune::MCMGElementLayout> ElementMapper;
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;
78 typedef Dune::VTKWriter<GridView> VtkWriter;
79 #if DUNE_VERSION_NEWER(DUNE_GRID, 2,5) 80 typedef std::shared_ptr< Dune::VTKFunction< GridView > > FunctionPtr;
82 typedef typename VtkWriter::VTKFunctionPtr FunctionPtr;
86 const std::string& simName =
"",
87 std::string multiFileName =
"")
89 , elementMapper_(gridView)
90 , vertexMapper_(gridView)
92 simName_ = (simName.empty()) ?
"sim" : simName;
93 multiFileName_ = multiFileName;
94 if (multiFileName_.empty()) {
95 multiFileName_ = simName_;
96 multiFileName_ +=
".pvd";
101 commRank_ = gridView.comm().rank();
102 commSize_ = gridView.comm().size();
117 {
return curWriterNum_; }
128 elementMapper_.update();
129 vertexMapper_.update();
137 if (!multiFile_.is_open()) {
138 startMultiFile_(multiFileName_);
142 curOutFileName_ = fileName_();
144 curWriter_ =
new VtkWriter(gridView_, Dune::VTK::conforming);
156 ScalarBuffer *buf =
new ScalarBuffer(numEntities);
157 managedScalarBuffers_.push_back(buf);
169 VectorBuffer *buf =
new VectorBuffer(numOuter);
170 for (
size_t i = 0; i < numOuter; ++ i)
171 (*buf)[i].resize(numInner);
173 managedVectorBuffers_.push_back(buf);
195 sanitizeScalarBuffer_(buf);
198 FunctionPtr fnPtr(
new VtkFn(name,
203 curWriter_->addVertexData(fnPtr);
223 sanitizeScalarBuffer_(buf);
226 FunctionPtr fnPtr(
new VtkFn(name,
231 curWriter_->addCellData(fnPtr);
252 sanitizeVectorBuffer_(buf);
255 FunctionPtr fnPtr(
new VtkFn(name,
260 curWriter_->addVertexData(fnPtr);
270 for (
unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
271 std::ostringstream oss;
272 oss << name <<
"[" << colIdx <<
"]";
274 FunctionPtr fnPtr(
new VtkFn(oss.str(),
280 curWriter_->addVertexData(fnPtr);
301 sanitizeVectorBuffer_(buf);
304 FunctionPtr fnPtr(
new VtkFn(name,
309 curWriter_->addCellData(fnPtr);
319 for (
unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
320 std::ostringstream oss;
321 oss << name <<
"[" << colIdx <<
"]";
323 FunctionPtr fnPtr(
new VtkFn(oss.str(),
329 curWriter_->addCellData(fnPtr);
343 std::string fileName;
345 fileName = curWriter_->write(curOutFileName_.c_str(),
346 static_cast<Dune::VTK::OutputType
>(vtkFormat));
350 multiFile_.precision(16);
351 multiFile_ <<
" <DataSet timestep=\"" << curTime_ <<
"\" file=\"" 352 << fileName <<
"\"/>\n";
359 while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
360 delete managedScalarBuffers_.front();
361 managedScalarBuffers_.pop_front();
363 while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
364 delete managedVectorBuffers_.front();
365 managedVectorBuffers_.pop_front();
377 template <
class Restarter>
380 res.serializeSectionBegin(
"VTKMultiWriter");
381 res.serializeStream() << curWriterNum_ <<
"\n";
383 if (commRank_ == 0) {
384 std::streamsize fileLen = 0;
385 std::streamoff filePos = 0;
386 if (multiFile_.is_open()) {
388 filePos = multiFile_.tellp();
389 multiFile_.seekp(0, std::ios::end);
390 fileLen = multiFile_.tellp();
391 multiFile_.seekp(filePos);
394 res.serializeStream() << fileLen <<
" " << filePos <<
"\n";
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);
405 res.serializeSectionEnd();
411 template <
class Restarter>
414 res.deserializeSectionBegin(
"VTKMultiWriter");
415 res.deserializeStream() >> curWriterNum_;
417 if (commRank_ == 0) {
419 std::getline(res.deserializeStream(), dummy);
422 std::streamoff filePos;
423 std::streamsize fileLen;
424 res.deserializeStream() >> fileLen >> filePos;
425 std::getline(res.deserializeStream(), dummy);
426 if (multiFile_.is_open())
430 multiFile_.open(multiFileName_.c_str());
432 char *tmp =
new char[fileLen];
433 res.deserializeStream().read(tmp, fileLen);
434 multiFile_.write(tmp, fileLen);
438 multiFile_.seekp(filePos);
442 std::getline(res.deserializeStream(), tmp);
444 res.deserializeSectionEnd();
448 std::string fileName_()
451 std::ostringstream oss;
452 oss << simName_ <<
"-" << std::setw(5) << std::setfill(
'0')
457 std::string fileSuffix_()
458 {
return (GridView::dimension == 1) ?
"vtp" :
"vtu"; }
460 void startMultiFile_(
const std::string& multiFileName)
463 if (commRank_ == 0) {
465 multiFile_.open(multiFileName.c_str());
466 multiFile_ <<
"<?xml version=\"1.0\"?>\n" 467 "<VTKFile type=\"Collection\"\n" 469 " byte_order=\"LittleEndian\"\n" 470 " compressor=\"vtkZLibDataCompressor\">\n" 475 void finishMultiFile_()
478 if (commRank_ == 0) {
480 std::ofstream::pos_type pos = multiFile_.tellp();
481 multiFile_ <<
" </Collection>\n" 483 multiFile_.seekp(pos);
490 void sanitizeScalarBuffer_(ScalarBuffer& b OPM_UNUSED)
495 void sanitizeVectorBuffer_(VectorBuffer& b OPM_UNUSED)
500 const GridView gridView_;
501 ElementMapper elementMapper_;
502 VertexMapper vertexMapper_;
504 std::string simName_;
505 std::ofstream multiFile_;
506 std::string multiFileName_;
511 VtkWriter *curWriter_;
513 std::string curOutFileName_;
516 std::list<ScalarBuffer *> managedScalarBuffers_;
517 std::list<VectorBuffer *> managedVectorBuffers_;
The base class for all output writers.
Definition: baseoutputwriter.hh:43
void serialize(Restarter &res)
Write the multi-writer'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'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