My Project
Loading...
Searching...
No Matches
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 explicit 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:
107 using Scalar = BaseOutputWriter::Scalar;
108 using Vector = BaseOutputWriter::Vector;
109 using Tensor = BaseOutputWriter::Tensor;
110 using ScalarBuffer = BaseOutputWriter::ScalarBuffer;
111 using VectorBuffer = BaseOutputWriter::VectorBuffer;
112 using TensorBuffer = BaseOutputWriter::TensorBuffer;
113
114 using VtkWriter = Dune::VTKWriter<GridView>;
115 using FunctionPtr = std::shared_ptr< Dune::VTKFunction< GridView > >;
116
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
142 ~VtkMultiWriter() override
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 elementMapper_.update(gridView_);
168 vertexMapper_.update(gridView_);
169 }
170
174 void beginWrite(double t) override
175 {
176 if (!multiFile_.is_open()) {
177 startMultiFile_(multiFileName_);
178 }
179
180 // make sure that all previous output has been written and no other thread
181 // accesses the memory used as the target for the extracted quantities
182 taskletRunner_.barrier();
183 releaseBuffers_();
184
185 curTime_ = t;
186 curOutFileName_ = fileName_();
187
188 curWriter_ = new VtkWriter(gridView_, Dune::VTK::conforming);
189 ++curWriterNum_;
190 }
191
199 {
200 ScalarBuffer *buf = new ScalarBuffer(numEntities);
201 managedScalarBuffers_.push_back(buf);
202 return buf;
203 }
204
211 VectorBuffer *allocateManagedVectorBuffer(size_t numOuter, size_t numInner)
212 {
213 VectorBuffer *buf = new VectorBuffer(numOuter);
214 for (size_t i = 0; i < numOuter; ++ i)
215 (*buf)[i].resize(numInner);
216
217 managedVectorBuffers_.push_back(buf);
218 return buf;
219 }
220
237 void attachScalarVertexData(ScalarBuffer& buf, std::string name) override
238 {
239 sanitizeScalarBuffer_(buf);
240
242 FunctionPtr fnPtr(new VtkFn(name,
243 gridView_,
244 vertexMapper_,
245 buf,
246 /*codim=*/dim));
247 curWriter_->addVertexData(fnPtr);
248 }
249
265 void attachScalarElementData(ScalarBuffer& buf, std::string name) override
266 {
267 sanitizeScalarBuffer_(buf);
268
270 FunctionPtr fnPtr(new VtkFn(name,
271 gridView_,
272 elementMapper_,
273 buf,
274 /*codim=*/0));
275 curWriter_->addCellData(fnPtr);
276 }
277
294 void attachVectorVertexData(VectorBuffer& buf, std::string name) override
295 {
296 sanitizeVectorBuffer_(buf);
297
299 FunctionPtr fnPtr(new VtkFn(name,
300 gridView_,
301 vertexMapper_,
302 buf,
303 /*codim=*/dim));
304 curWriter_->addVertexData(fnPtr);
305 }
306
310 void attachTensorVertexData(TensorBuffer& buf, std::string name) override
311 {
313
314 for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
315 std::ostringstream oss;
316 oss << name << "[" << colIdx << "]";
317
318 FunctionPtr fnPtr(new VtkFn(oss.str(),
319 gridView_,
320 vertexMapper_,
321 buf,
322 /*codim=*/dim,
323 colIdx));
324 curWriter_->addVertexData(fnPtr);
325 }
326 }
327
343 void attachVectorElementData(VectorBuffer& buf, std::string name) override
344 {
345 sanitizeVectorBuffer_(buf);
346
348 FunctionPtr fnPtr(new VtkFn(name,
349 gridView_,
350 elementMapper_,
351 buf,
352 /*codim=*/0));
353 curWriter_->addCellData(fnPtr);
354 }
355
359 void attachTensorElementData(TensorBuffer& buf, std::string name) override
360 {
362
363 for (unsigned colIdx = 0; colIdx < buf[0].N(); ++colIdx) {
364 std::ostringstream oss;
365 oss << name << "[" << colIdx << "]";
366
367 FunctionPtr fnPtr(new VtkFn(oss.str(),
368 gridView_,
369 elementMapper_,
370 buf,
371 /*codim=*/0,
372 colIdx));
373 curWriter_->addCellData(fnPtr);
374 }
375 }
376
384 void endWrite(bool onlyDiscard = false) override
385 {
386 if (!onlyDiscard) {
387 auto tasklet = std::make_shared<WriteDataTasklet>(*this);
388 taskletRunner_.dispatch(tasklet);
389 }
390 else
391 --curWriterNum_;
392
393 // temporarily write the closing XML mumbo-jumbo to the mashup
394 // file so that the data set can be loaded even if the
395 // simulation is aborted (or not yet finished)
396 finishMultiFile_();
397 }
398
402 template <class Restarter>
404 {
405 res.serializeSectionBegin("VTKMultiWriter");
406 res.serializeStream() << curWriterNum_ << "\n";
407
408 if (commRank_ == 0) {
409 std::streamsize fileLen = 0;
410 std::streamoff filePos = 0;
411 if (multiFile_.is_open()) {
412 // write the meta file into the restart file
413 filePos = multiFile_.tellp();
414 multiFile_.seekp(0, std::ios::end);
415 fileLen = multiFile_.tellp();
416 multiFile_.seekp(filePos);
417 }
418
419 res.serializeStream() << fileLen << " " << filePos << "\n";
420
421 if (fileLen > 0) {
422 std::ifstream multiFileIn(multiFileName_.c_str());
423 char *tmp = new char[fileLen];
424 multiFileIn.read(tmp, static_cast<long>(fileLen));
425 res.serializeStream().write(tmp, fileLen);
426 delete[] tmp;
427 }
428 }
429
430 res.serializeSectionEnd();
431 }
432
436 template <class Restarter>
438 {
439 res.deserializeSectionBegin("VTKMultiWriter");
440 res.deserializeStream() >> curWriterNum_;
441
442 if (commRank_ == 0) {
443 std::string dummy;
444 std::getline(res.deserializeStream(), dummy);
445
446 // recreate the meta file from the restart file
447 std::streamoff filePos;
448 std::streamsize fileLen;
449 res.deserializeStream() >> fileLen >> filePos;
450 std::getline(res.deserializeStream(), dummy);
451 if (multiFile_.is_open())
452 multiFile_.close();
453
454 if (fileLen > 0) {
455 multiFile_.open(multiFileName_.c_str());
456
457 char *tmp = new char[fileLen];
458 res.deserializeStream().read(tmp, fileLen);
459 multiFile_.write(tmp, fileLen);
460 delete[] tmp;
461 }
462
463 multiFile_.seekp(filePos);
464 }
465 else {
466 std::string tmp;
467 std::getline(res.deserializeStream(), tmp);
468 }
469 res.deserializeSectionEnd();
470 }
471
472private:
473 std::string fileName_()
474 {
475 // use a new file name for each time step
476 std::ostringstream oss;
477 oss << simName_ << "-"
478 << std::setw(5) << std::setfill('0') << curWriterNum_;
479 return oss.str();
480 }
481
482 std::string fileSuffix_()
483 { return (GridView::dimension == 1) ? "vtp" : "vtu"; }
484
485 void startMultiFile_(const std::string& multiFileName)
486 {
487 // only the first process writes to the multi-file
488 if (commRank_ == 0) {
489 // generate one meta vtk-file holding the individual time steps
490 multiFile_.open(multiFileName.c_str());
491 multiFile_ << "<?xml version=\"1.0\"?>\n"
492 "<VTKFile type=\"Collection\"\n"
493 " version=\"0.1\"\n"
494 " byte_order=\"LittleEndian\"\n"
495 " compressor=\"vtkZLibDataCompressor\">\n"
496 " <Collection>\n";
497 }
498 }
499
500 void finishMultiFile_()
501 {
502 // only the first process writes to the multi-file
503 if (commRank_ == 0) {
504 // make sure that we always have a working meta file
505 std::ofstream::pos_type pos = multiFile_.tellp();
506 multiFile_ << " </Collection>\n"
507 "</VTKFile>\n";
508 multiFile_.seekp(pos);
509 multiFile_.flush();
510 }
511 }
512
513 // make sure the field is well defined if running under valgrind
514 // and make sure that all values can be displayed by paraview
515 void sanitizeScalarBuffer_(ScalarBuffer&)
516 {
517 // nothing to do: this is done by VtkScalarFunction
518 }
519
520 void sanitizeVectorBuffer_(VectorBuffer&)
521 {
522 // nothing to do: this is done by VtkVectorFunction
523 }
524
525 // release the memory occupied by all buffer objects managed by the multi-writer
526 void releaseBuffers_()
527 {
528 // discard managed objects and the current VTK writer
529 delete curWriter_;
530 curWriter_ = nullptr;
531 while (managedScalarBuffers_.begin() != managedScalarBuffers_.end()) {
532 delete managedScalarBuffers_.front();
533 managedScalarBuffers_.pop_front();
534 }
535 while (managedVectorBuffers_.begin() != managedVectorBuffers_.end()) {
536 delete managedVectorBuffers_.front();
537 managedVectorBuffers_.pop_front();
538 }
539 }
540
541 const GridView gridView_;
542 ElementMapper elementMapper_;
543 VertexMapper vertexMapper_;
544
545 std::string outputDir_;
546 std::string simName_;
547 std::ofstream multiFile_;
548 std::string multiFileName_;
549
550 int commSize_; // number of processes in the communicator
551 int commRank_; // rank of the current process in the communicator
552
553 VtkWriter *curWriter_;
554 double curTime_;
555 std::string curOutFileName_;
556 int curWriterNum_;
557
558 std::list<ScalarBuffer *> managedScalarBuffers_;
559 std::list<VectorBuffer *> managedVectorBuffers_;
560
561 TaskletRunner taskletRunner_;
562};
563} // namespace Opm
564
565#endif
The base class for all output writers.
The base class for all output writers.
Definition baseoutputwriter.hh:44
The base class for tasklets.
Definition tasklets.hpp:44
void barrier()
Make sure that all tasklets have been completed after this method has been called.
Definition tasklets.cpp:134
void dispatch(std::shared_ptr< TaskletInterface > tasklet)
Add a new tasklet.
Definition tasklets.cpp:101
Simplifies writing multi-file VTK datasets.
Definition vtkmultiwriter.hh:66
void attachVectorElementData(VectorBuffer &buf, std::string name) override
Add a element centered quantity to the output.
Definition vtkmultiwriter.hh:343
void endWrite(bool onlyDiscard=false) override
Finalizes the current writer.
Definition vtkmultiwriter.hh:384
void attachTensorVertexData(TensorBuffer &buf, std::string name) override
Add a finished vertex-centered tensor field to the output.
Definition vtkmultiwriter.hh:310
void attachTensorElementData(TensorBuffer &buf, std::string name) override
Add a finished element-centered tensor field to the output.
Definition vtkmultiwriter.hh:359
void deserialize(Restarter &res)
Read the multi-writer's state from a restart file.
Definition vtkmultiwriter.hh:437
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:211
void gridChanged()
Updates the internal data structures after mesh refinement.
Definition vtkmultiwriter.hh:165
void beginWrite(double t) override
Called whenever a new time step must be written.
Definition vtkmultiwriter.hh:174
ScalarBuffer * allocateManagedScalarBuffer(size_t numEntities)
Allocate a managed buffer for a scalar field.
Definition vtkmultiwriter.hh:198
void attachScalarElementData(ScalarBuffer &buf, std::string name) override
Add a element centered quantity to the output.
Definition vtkmultiwriter.hh:265
void attachScalarVertexData(ScalarBuffer &buf, std::string name) override
Add a finished vertex centered vector field to the output.
Definition vtkmultiwriter.hh:237
void serialize(Restarter &res)
Write the multi-writer's state to a restart file.
Definition vtkmultiwriter.hh:403
void attachVectorVertexData(VectorBuffer &buf, std::string name) override
Add a finished vertex centered vector field to the output.
Definition vtkmultiwriter.hh:294
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
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
Provides a mechanism to dispatch work to separate threads.
Provides a vector-valued function using Dune::FieldVectors as elements.
Provides a tensor-valued function using Dune::FieldMatrix objects as elements.
Provides a vector-valued function using Dune::FieldVectors as elements.