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
37
38#include <opm/material/common/Valgrind.hpp>
39
40#include <dune/common/fvector.hh>
41#include <dune/common/version.hh>
42#include <dune/istl/bvector.hh>
43#include <dune/grid/io/file/vtk/vtkwriter.hh>
44
45#if HAVE_MPI
46#include <mpi.h>
47#endif
48
49#include <filesystem>
50#include <list>
51#include <string>
52#include <limits>
53#include <sstream>
54#include <fstream>
55
56namespace Opm {
64template <class GridView, int vtkFormat>
66{
67 class WriteDataTasklet : public TaskletInterface
68 {
69 public:
70 WriteDataTasklet(VtkMultiWriter& multiWriter)
71 : multiWriter_(multiWriter)
72 { }
73
74 void run() final
75 {
76 std::string fileName;
77 // write the actual data as vtu or vtp (plus the pieces file in the parallel case)
78 if (multiWriter_.commSize_ > 1)
79 fileName = multiWriter_.curWriter_->pwrite(/*name=*/multiWriter_.curOutFileName_,
80 /*path=*/multiWriter_.outputDir_,
81 /*extendPath=*/"",
82 static_cast<Dune::VTK::OutputType>(vtkFormat));
83 else
84 fileName = multiWriter_.curWriter_->write(/*name=*/multiWriter_.outputDir_ + "/" + multiWriter_.curOutFileName_,
85 static_cast<Dune::VTK::OutputType>(vtkFormat));
86
87 // determine name to write into the multi-file for the
88 // current time step
89 // The file names in the pvd file are relative, the path should therefore be stripped.
90 const std::filesystem::path fullPath{fileName};
91 const std::string localFileName = fullPath.filename();
92 multiWriter_.multiFile_.precision(16);
93 multiWriter_.multiFile_ << " <DataSet timestep=\"" << multiWriter_.curTime_ << "\" file=\""
94 << localFileName << "\"/>\n";
95 }
96
97 private:
98 VtkMultiWriter& multiWriter_;
99 };
100
101 enum { dim = GridView::dimension };
102
103 using VertexMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
104 using ElementMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
105
106public:
113
114 using VtkWriter = Dune::VTKWriter<GridView>;
115 using FunctionPtr = std::shared_ptr< Dune::VTKFunction< GridView > >;
116
117 VtkMultiWriter(bool asyncWriting,
118 const GridView& gridView,
119 const std::string& outputDir,
120 const std::string& simName = "",
121 std::string multiFileName = "")
122 : gridView_(gridView)
123 , elementMapper_(gridView, Dune::mcmgElementLayout())
124 , vertexMapper_(gridView, Dune::mcmgVertexLayout())
125 , curWriter_(nullptr)
126 , curWriterNum_(0)
127 , taskletRunner_(/*numThreads=*/asyncWriting?1:0)
128 {
129 outputDir_ = outputDir;
130 if (outputDir == "")
131 outputDir_ = ".";
132
133 simName_ = (simName.empty()) ? "sim" : simName;
134 multiFileName_ = multiFileName;
135 if (multiFileName_.empty())
136 multiFileName_ = outputDir_+"/"+simName_+".pvd";
137
138 commRank_ = gridView.comm().rank();
139 commSize_ = gridView.comm().size();
140 }
141
143 {
144 taskletRunner_.barrier();
145 releaseBuffers_();
146 finishMultiFile_();
147
148 if (commRank_ == 0)
149 multiFile_.close();
150 }
151
155 int curWriterNum() const
156 { return curWriterNum_; }
157
166 {
167#if DUNE_VERSION_NEWER(DUNE_GRID, 2, 8)
168 elementMapper_.update(gridView_);
169 vertexMapper_.update(gridView_);
170#else
171 elementMapper_.update();
172 vertexMapper_.update();
173#endif
174 }
175
179 void beginWrite(double t)
180 {
181 if (!multiFile_.is_open()) {
182 startMultiFile_(multiFileName_);
183 }
184
185 // make sure that all previous output has been written and no other thread
186 // accesses the memory used as the target for the extracted quantities
187 taskletRunner_.barrier();
188 releaseBuffers_();
189
190 curTime_ = t;
191 curOutFileName_ = fileName_();
192
193 curWriter_ = new VtkWriter(gridView_, Dune::VTK::conforming);
194 ++curWriterNum_;
195 }
196
204 {
205 ScalarBuffer *buf = new ScalarBuffer(numEntities);
206 managedScalarBuffers_.push_back(buf);
207 return buf;
208 }
209
216 VectorBuffer *allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
217 {
218 VectorBuffer *buf = new VectorBuffer(numOuter);
219 for (size_t i = 0; i < numOuter; ++ i)
220 (*buf)[i].resize(numInner);
221
222 managedVectorBuffers_.push_back(buf);
223 return buf;
224 }
225
242 void attachScalarVertexData(ScalarBuffer& buf, std::string name)
243 {
244 sanitizeScalarBuffer_(buf);
245
247 FunctionPtr fnPtr(new VtkFn(name,
248 gridView_,
249 vertexMapper_,
250 buf,
251 /*codim=*/dim));
252 curWriter_->addVertexData(fnPtr);
253 }
254
270 void attachScalarElementData(ScalarBuffer& buf, std::string name)
271 {
272 sanitizeScalarBuffer_(buf);
273
275 FunctionPtr fnPtr(new VtkFn(name,
276 gridView_,
277 elementMapper_,
278 buf,
279 /*codim=*/0));
280 curWriter_->addCellData(fnPtr);
281 }
282
299 void attachVectorVertexData(VectorBuffer& buf, std::string name)
300 {
301 sanitizeVectorBuffer_(buf);
302
304 FunctionPtr fnPtr(new VtkFn(name,
305 gridView_,
306 vertexMapper_,
307 buf,
308 /*codim=*/dim));
309 curWriter_->addVertexData(fnPtr);
310 }
311
315 void attachTensorVertexData(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 vertexMapper_,
326 buf,
327 /*codim=*/dim,
328 colIdx));
329 curWriter_->addVertexData(fnPtr);
330 }
331 }
332
348 void attachVectorElementData(VectorBuffer& buf, std::string name)
349 {
350 sanitizeVectorBuffer_(buf);
351
353 FunctionPtr fnPtr(new VtkFn(name,
354 gridView_,
355 elementMapper_,
356 buf,
357 /*codim=*/0));
358 curWriter_->addCellData(fnPtr);
359 }
360
364 void attachTensorElementData(TensorBuffer& buf, std::string name)
365 {
367
368 for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
369 std::ostringstream oss;
370 oss << name << "[" << colIdx << "]";
371
372 FunctionPtr fnPtr(new VtkFn(oss.str(),
373 gridView_,
374 elementMapper_,
375 buf,
376 /*codim=*/0,
377 colIdx));
378 curWriter_->addCellData(fnPtr);
379 }
380 }
381
389 void endWrite(bool onlyDiscard = false)
390 {
391 if (!onlyDiscard) {
392 auto tasklet = std::make_shared<WriteDataTasklet>(*this);
393 taskletRunner_.dispatch(tasklet);
394 }
395 else
396 --curWriterNum_;
397
398 // temporarily write the closing XML mumbo-jumbo to the mashup
399 // file so that the data set can be loaded even if the
400 // simulation is aborted (or not yet finished)
401 finishMultiFile_();
402 }
403
407 template <class Restarter>
408 void serialize(Restarter& res)
409 {
410 res.serializeSectionBegin("VTKMultiWriter");
411 res.serializeStream() << curWriterNum_ << "\n";
412
413 if (commRank_ == 0) {
414 std::streamsize fileLen = 0;
415 std::streamoff filePos = 0;
416 if (multiFile_.is_open()) {
417 // write the meta file into the restart file
418 filePos = multiFile_.tellp();
419 multiFile_.seekp(0, std::ios::end);
420 fileLen = multiFile_.tellp();
421 multiFile_.seekp(filePos);
422 }
423
424 res.serializeStream() << fileLen << " " << filePos << "\n";
425
426 if (fileLen > 0) {
427 std::ifstream multiFileIn(multiFileName_.c_str());
428 char *tmp = new char[fileLen];
429 multiFileIn.read(tmp, static_cast<long>(fileLen));
430 res.serializeStream().write(tmp, fileLen);
431 delete[] tmp;
432 }
433 }
434
435 res.serializeSectionEnd();
436 }
437
441 template <class Restarter>
442 void deserialize(Restarter& res)
443 {
444 res.deserializeSectionBegin("VTKMultiWriter");
445 res.deserializeStream() >> curWriterNum_;
446
447 if (commRank_ == 0) {
448 std::string dummy;
449 std::getline(res.deserializeStream(), dummy);
450
451 // recreate the meta file from the restart file
452 std::streamoff filePos;
453 std::streamsize fileLen;
454 res.deserializeStream() >> fileLen >> filePos;
455 std::getline(res.deserializeStream(), dummy);
456 if (multiFile_.is_open())
457 multiFile_.close();
458
459 if (fileLen > 0) {
460 multiFile_.open(multiFileName_.c_str());
461
462 char *tmp = new char[fileLen];
463 res.deserializeStream().read(tmp, fileLen);
464 multiFile_.write(tmp, fileLen);
465 delete[] tmp;
466 }
467
468 multiFile_.seekp(filePos);
469 }
470 else {
471 std::string tmp;
472 std::getline(res.deserializeStream(), tmp);
473 }
474 res.deserializeSectionEnd();
475 }
476
477private:
478 std::string fileName_()
479 {
480 // use a new file name for each time step
481 std::ostringstream oss;
482 oss << simName_ << "-"
483 << std::setw(5) << std::setfill('0') << curWriterNum_;
484 return oss.str();
485 }
486
487 std::string fileSuffix_()
488 { return (GridView::dimension == 1) ? "vtp" : "vtu"; }
489
490 void startMultiFile_(const std::string& multiFileName)
491 {
492 // only the first process writes to the multi-file
493 if (commRank_ == 0) {
494 // generate one meta vtk-file holding the individual time steps
495 multiFile_.open(multiFileName.c_str());
496 multiFile_ << "<?xml version=\"1.0\"?>\n"
497 "<VTKFile type=\"Collection\"\n"
498 " version=\"0.1\"\n"
499 " byte_order=\"LittleEndian\"\n"
500 " compressor=\"vtkZLibDataCompressor\">\n"
501 " <Collection>\n";
502 }
503 }
504
505 void finishMultiFile_()
506 {
507 // only the first process writes to the multi-file
508 if (commRank_ == 0) {
509 // make sure that we always have a working meta file
510 std::ofstream::pos_type pos = multiFile_.tellp();
511 multiFile_ << " </Collection>\n"
512 "</VTKFile>\n";
513 multiFile_.seekp(pos);
514 multiFile_.flush();
515 }
516 }
517
518 // make sure the field is well defined if running under valgrind
519 // and make sure that all values can be displayed by paraview
520 void sanitizeScalarBuffer_(ScalarBuffer&)
521 {
522 // nothing to do: this is done by VtkScalarFunction
523 }
524
525 void sanitizeVectorBuffer_(VectorBuffer&)
526 {
527 // nothing to do: this is done by VtkVectorFunction
528 }
529
530 // release the memory occupied by all buffer objects managed by the multi-writer
531 void releaseBuffers_()
532 {
533 // discard managed objects and the current VTK writer
534 delete curWriter_;
535 curWriter_ = nullptr;
536 while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
537 delete managedScalarBuffers_.front();
538 managedScalarBuffers_.pop_front();
539 }
540 while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
541 delete managedVectorBuffers_.front();
542 managedVectorBuffers_.pop_front();
543 }
544 }
545
546 const GridView gridView_;
547 ElementMapper elementMapper_;
548 VertexMapper vertexMapper_;
549
550 std::string outputDir_;
551 std::string simName_;
552 std::ofstream multiFile_;
553 std::string multiFileName_;
554
555 int commSize_; // number of processes in the communicator
556 int commRank_; // rank of the current process in the communicator
557
558 VtkWriter *curWriter_;
559 double curTime_;
560 std::string curOutFileName_;
561 int curWriterNum_;
562
563 std::list<ScalarBuffer *> managedScalarBuffers_;
564 std::list<VectorBuffer *> managedVectorBuffers_;
565
566 TaskletRunner taskletRunner_;
567};
568} // namespace Opm
569
570#endif
The base class for all output writers.
Definition: baseoutputwriter.hh:44
Dune::DynamicVector< double > Vector
Definition: baseoutputwriter.hh:47
std::vector< Tensor > TensorBuffer
Definition: baseoutputwriter.hh:51
Dune::DynamicMatrix< double > Tensor
Definition: baseoutputwriter.hh:48
std::vector< Vector > VectorBuffer
Definition: baseoutputwriter.hh:50
std::vector< Scalar > ScalarBuffer
Definition: baseoutputwriter.hh:49
double Scalar
Definition: baseoutputwriter.hh:46
The base class for tasklets.
Definition: tasklets.hh:47
void barrier()
Make sure that all tasklets have been completed after this method has been called.
Definition: tasklets.hh:278
void dispatch(std::shared_ptr< TaskletInterface > tasklet)
Add a new tasklet.
Definition: tasklets.hh:230
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:66
void endWrite(bool onlyDiscard=false)
Finalizes the current writer.
Definition: vtkmultiwriter.hh:389
void attachTensorVertexData(TensorBuffer &buf, std::string name)
Add a finished vertex-centered tensor field to the output.
Definition: vtkmultiwriter.hh:315
VtkMultiWriter(bool asyncWriting, const GridView &gridView, const std::string &outputDir, const std::string &simName="", std::string multiFileName="")
Definition: vtkmultiwriter.hh:117
void deserialize(Restarter &res)
Read the multi-writer's state from a restart file.
Definition: vtkmultiwriter.hh:442
void attachTensorElementData(TensorBuffer &buf, std::string name)
Add a finished element-centered tensor field to the output.
Definition: vtkmultiwriter.hh:364
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: vtkmultiwriter.hh:110
int curWriterNum() const
Returns the number of the current VTK file.
Definition: vtkmultiwriter.hh:155
VectorBuffer * allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
Allocate a managed buffer for a vector field.
Definition: vtkmultiwriter.hh:216
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition: vtkmultiwriter.hh:165
void beginWrite(double t)
Called whenever a new time step must be written.
Definition: vtkmultiwriter.hh:179
~VtkMultiWriter()
Definition: vtkmultiwriter.hh:142
ScalarBuffer * allocateManagedScalarBuffer(size_t numEntities)
Allocate a managed buffer for a scalar field.
Definition: vtkmultiwriter.hh:203
void attachScalarVertexData(ScalarBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:242
void attachVectorElementData(VectorBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:348
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: vtkmultiwriter.hh:111
void attachScalarElementData(ScalarBuffer &buf, std::string name)
Add a element centered quantity to the output.
Definition: vtkmultiwriter.hh:270
Dune::VTKWriter< GridView > VtkWriter
Definition: vtkmultiwriter.hh:114
std::shared_ptr< Dune::VTKFunction< GridView > > FunctionPtr
Definition: vtkmultiwriter.hh:115
void serialize(Restarter &res)
Write the multi-writer's state to a restart file.
Definition: vtkmultiwriter.hh:408
void attachVectorVertexData(VectorBuffer &buf, std::string name)
Add a finished vertex centered vector field to the output.
Definition: vtkmultiwriter.hh:299
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkscalarfunction.hh:50
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
Definition: vtktensorfunction.hh:47
Provides a vector-valued function using Dune::FieldVectors as elements.
Definition: vtkvectorfunction.hh:50
Definition: fvbaseprimaryvariables.hh:141
Definition: blackoilboundaryratevector.hh:37
Provides a mechanism to dispatch work to separate threads.