3.4
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
io/vtkoutputmodule.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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_IO_VTK_OUTPUT_MODULE_HH
25#define DUMUX_IO_VTK_OUTPUT_MODULE_HH
26
27#include <functional>
28#include <memory>
29#include <string>
30
31#include <dune/common/timer.hh>
32#include <dune/common/fvector.hh>
33#include <dune/common/typetraits.hh>
34
35#include <dune/geometry/type.hh>
36#include <dune/geometry/multilineargeometry.hh>
37
38#include <dune/grid/common/mcmgmapper.hh>
39#include <dune/grid/common/partitionset.hh>
40#include <dune/grid/io/file/vtk/vtkwriter.hh>
41#include <dune/grid/io/file/vtk/vtksequencewriter.hh>
42#include <dune/grid/common/partitionset.hh>
43
44#include <dumux/common/parameters.hh>
45#include <dumux/io/format.hh>
46#include <dumux/discretization/method.hh>
47
48#include <dumux/io/vtk/function.hh>
49#include <dumux/io/vtk/fieldtype.hh>
50#include "velocityoutput.hh"
51
52namespace Dumux {
53
59template<class GridGeometry>
60class VtkOutputModuleBase
61{
62 using GridView = typename GridGeometry::GridView;
63 using Field = Vtk::template Field<GridView>;
64 static constexpr int dim = GridView::dimension;
65
66public:
68 struct [[deprecated("use Vtk::FieldType instead")]] FieldType
69 {
70 static constexpr auto element = Vtk::FieldType::element;
71 static constexpr auto vertex = Vtk::FieldType::vertex;
72 static constexpr auto automatic = Vtk::FieldType::automatic;
73 };
74
75 VtkOutputModuleBase(const GridGeometry& gridGeometry,
76 const std::string& name,
77 const std::string& paramGroup = "",
78 Dune::VTK::DataMode dm = Dune::VTK::conforming,
79 bool verbose = true)
80 : gridGeometry_(gridGeometry)
81 , name_(name)
82 , paramGroup_(paramGroup)
83 , dm_(dm)
84 , verbose_(gridGeometry.gridView().comm().rank() == 0 && verbose)
85 {
86 const auto precisionString = getParamFromGroup<std::string>(paramGroup, "Vtk.Precision", "Float32");
87 precision_ = Dumux::Vtk::stringToPrecision(precisionString);
88 const auto coordPrecision = Dumux::Vtk::stringToPrecision(getParamFromGroup<std::string>(paramGroup, "Vtk.CoordPrecision", precisionString));
89 writer_ = std::make_shared<Dune::VTKWriter<GridView>>(gridGeometry.gridView(), dm, coordPrecision);
90 sequenceWriter_ = std::make_unique<Dune::VTKSequenceWriter<GridView>>(writer_, name);
91 }
92
93 virtual ~VtkOutputModuleBase() = default;
94
96 const std::string& paramGroup() const
97 { return paramGroup_; }
98
108 template<typename Vector>
109 void addField(const Vector& v,
110 const std::string& name,
111 Vtk::FieldType fieldType = Vtk::FieldType::automatic)
112 { addField(v, name, this->precision(), fieldType); }
113
124 template<typename Vector>
125 void addField(const Vector& v,
126 const std::string& name,
127 Dumux::Vtk::Precision precision,
128 Vtk::FieldType fieldType = Vtk::FieldType::automatic)
129 {
130 // Deduce the number of components from the given vector type
131 const auto nComp = getNumberOfComponents_(v);
132
133 const auto numElemDofs = gridGeometry().elementMapper().size();
134 const auto numVertexDofs = gridGeometry().vertexMapper().size();
135
136 // Automatically deduce the field type ...
137 if(fieldType == Vtk::FieldType::automatic)
138 {
139 if(numElemDofs == numVertexDofs)
140 DUNE_THROW(Dune::InvalidStateException, "Automatic deduction of FieldType failed. Please explicitly specify FieldType::element or FieldType::vertex.");
141
142 if(v.size() == numElemDofs)
143 fieldType = Vtk::FieldType::element;
144 else if(v.size() == numVertexDofs)
145 fieldType = Vtk::FieldType::vertex;
146 else
147 DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
148 }
149 // ... or check if the user-specified type matches the size of v
150 else
151 {
152 if(fieldType == Vtk::FieldType::element)
153 if(v.size() != numElemDofs)
154 DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
155
156 if(fieldType == Vtk::FieldType::vertex)
157 if(v.size() != numVertexDofs)
158 DUNE_THROW(Dune::RangeError, "Size mismatch of added field!");
159 }
160
161 // add the appropriate field
162 if (fieldType == Vtk::FieldType::element)
163 fields_.emplace_back(gridGeometry_.gridView(), gridGeometry_.elementMapper(), v, name, nComp, 0, dm_, precision);
164 else
165 fields_.emplace_back(gridGeometry_.gridView(), gridGeometry_.vertexMapper(), v, name, nComp, dim, dm_, precision);
166 }
167
173 void write(double time, Dune::VTK::OutputType type = Dune::VTK::ascii)
174 {
175 Dune::Timer timer;
176
177 // write to file depending on data mode
178 if (dm_ == Dune::VTK::conforming)
179 writeConforming_(time, type);
180 else if (dm_ == Dune::VTK::nonconforming)
181 writeNonConforming_(time, type);
182 else
183 DUNE_THROW(Dune::NotImplemented, "Output for provided VTK data mode");
184
186 timer.stop();
187 if (verbose_)
188 std::cout << Fmt::format("Writing output for problem \"{}\". Took {:.2g} seconds.\n", name_, timer.elapsed());
189 }
190
191protected:
192 const GridGeometry& gridGeometry() const { return gridGeometry_; }
193
194 bool verbose() const { return verbose_; }
195 const std::string& name() const { return name_; }
196 Dune::VTK::DataMode dataMode() const { return dm_; }
197 Dumux::Vtk::Precision precision() const { return precision_; }
198
199 Dune::VTKWriter<GridView>& writer() { return *writer_; }
200 Dune::VTKSequenceWriter<GridView>& sequenceWriter() { return *sequenceWriter_; }
201
202 const std::vector<Field>& fields() const { return fields_; }
203
204private:
206 virtual void writeConforming_(double time, Dune::VTK::OutputType type)
207 {
211
212 // process rank
213 static bool addProcessRank = getParamFromGroup<bool>(this->paramGroup(), "Vtk.AddProcessRank");
214 std::vector<int> rank;
215
217 if (!fields_.empty() || addProcessRank)
218 {
219 const auto numCells = gridGeometry_.gridView().size(0);
220
221 // maybe allocate space for the process rank
222 if (addProcessRank)
223 {
224 rank.resize(numCells);
225
226 for (const auto& element : elements(gridGeometry_.gridView(), Dune::Partitions::interior))
227 {
228 const auto eIdxGlobal = gridGeometry_.elementMapper().index(element);
229 rank[eIdxGlobal] = gridGeometry_.gridView().comm().rank();
230 }
231 }
232
236
237 // the process rank
238 if (addProcessRank)
239 this->sequenceWriter().addCellData(Field(gridGeometry_.gridView(), gridGeometry_.elementMapper(), rank, "process rank", 1, 0).get());
240
241 // also register additional (non-standardized) user fields if any
242 for (auto&& field : fields_)
243 {
244 if (field.codim() == 0)
245 this->sequenceWriter().addCellData(field.get());
246 else if (field.codim() == dim)
247 this->sequenceWriter().addVertexData(field.get());
248 else
249 DUNE_THROW(Dune::RangeError, "Cannot add wrongly sized vtk scalar field!");
250 }
251 }
252
256 this->sequenceWriter().write(time, type);
257
261 this->writer().clear();
262 }
263
265 virtual void writeNonConforming_(double time, Dune::VTK::OutputType type)
266 {
267 DUNE_THROW(Dune::NotImplemented, "Non-conforming VTK output");
268 }
269
271 template<class Vector>
272 std::size_t getNumberOfComponents_(const Vector& v)
273 {
274 if constexpr (Dune::IsIndexable<decltype(std::declval<Vector>()[0])>::value)
275 return v[0].size();
276 else
277 return 1;
278 }
279
280 const GridGeometry& gridGeometry_;
281 std::string name_;
282 const std::string paramGroup_;
283 Dune::VTK::DataMode dm_;
284 bool verbose_;
285 Dumux::Vtk::Precision precision_;
286
287 std::shared_ptr<Dune::VTKWriter<GridView>> writer_;
288 std::unique_ptr<Dune::VTKSequenceWriter<GridView>> sequenceWriter_;
289
290 std::vector<Field> fields_;
291};
292
305template<class GridVariables, class SolutionVector>
306class VtkOutputModule : public VtkOutputModuleBase<typename GridVariables::GridGeometry>
307{
308 using ParentType = VtkOutputModuleBase<typename GridVariables::GridGeometry>;
309 using GridGeometry = typename GridVariables::GridGeometry;
310
311 using VV = typename GridVariables::VolumeVariables;
312 using Scalar = typename GridVariables::Scalar;
313
314 using GridView = typename GridGeometry::GridView;
315
316 enum {
317 dim = GridView::dimension,
318 dimWorld = GridView::dimensionworld
319 };
320
321 using Element = typename GridView::template Codim<0>::Entity;
322 using VolVarsVector = Dune::FieldVector<Scalar, dimWorld>;
323
324 static constexpr bool isBox = GridGeometry::discMethod == DiscretizationMethod::box;
325 static constexpr int dofCodim = isBox ? dim : 0;
326
327 struct VolVarScalarDataInfo { std::function<Scalar(const VV&)> get; std::string name; Dumux::Vtk::Precision precision_; };
328 struct VolVarVectorDataInfo { std::function<VolVarsVector(const VV&)> get; std::string name; Dumux::Vtk::Precision precision_; };
329 using Field = Vtk::template Field<GridView>;
330
331 using VelocityOutputType = Dumux::VelocityOutput<GridVariables>;
332
333public:
335 using VolumeVariables = VV;
336
337 VtkOutputModule(const GridVariables& gridVariables,
338 const SolutionVector& sol,
339 const std::string& name,
340 const std::string& paramGroup = "",
341 Dune::VTK::DataMode dm = Dune::VTK::conforming,
342 bool verbose = true)
343 : ParentType(gridVariables.gridGeometry(), name, paramGroup, dm, verbose)
344 , gridVariables_(gridVariables)
345 , sol_(sol)
346 , velocityOutput_(std::make_shared<VelocityOutputType>())
347 {}
348
353
360 void addVelocityOutput(std::shared_ptr<VelocityOutputType> velocityOutput)
361 { velocityOutput_ = velocityOutput; }
362
366 void addVolumeVariable(std::function<Scalar(const VolumeVariables&)>&& f,
367 const std::string& name)
368 {
369 volVarScalarDataInfo_.push_back(VolVarScalarDataInfo{f, name, this->precision()});
370 }
371
376 template<class VVV = VolVarsVector, typename std::enable_if_t<(VVV::dimension > 1), int> = 0>
377 void addVolumeVariable(std::function<VolVarsVector(const VolumeVariables&)>&& f,
378 const std::string& name)
379 {
380 volVarVectorDataInfo_.push_back(VolVarVectorDataInfo{f, name, this->precision()});
381 }
382
383protected:
384 // some return functions for differing implementations to use
385 const auto& problem() const { return gridVariables_.curGridVolVars().problem(); }
386 const GridVariables& gridVariables() const { return gridVariables_; }
387 const GridGeometry& gridGeometry() const { return gridVariables_.gridGeometry(); }
388 const SolutionVector& sol() const { return sol_; }
389
390 const std::vector<VolVarScalarDataInfo>& volVarScalarDataInfo() const { return volVarScalarDataInfo_; }
391 const std::vector<VolVarVectorDataInfo>& volVarVectorDataInfo() const { return volVarVectorDataInfo_; }
392
393 using VelocityOutput = VelocityOutputType;
394 const VelocityOutput& velocityOutput() const { return *velocityOutput_; }
395
396private:
397
399 void writeConforming_(double time, Dune::VTK::OutputType type) override
400 {
401 const Dune::VTK::DataMode dm = Dune::VTK::conforming;
405
406 // instantiate the velocity output
407 using VelocityVector = typename VelocityOutput::VelocityVector;
408 std::vector<VelocityVector> velocity(velocityOutput_->numFluidPhases());
409
410 // process rank
411 static bool addProcessRank = getParamFromGroup<bool>(this->paramGroup(), "Vtk.AddProcessRank");
412 std::vector<double> rank;
413
414 // volume variable data
415 std::vector<std::vector<Scalar>> volVarScalarData;
416 std::vector<std::vector<VolVarsVector>> volVarVectorData;
417
419 if (!volVarScalarDataInfo_.empty()
420 || !volVarVectorDataInfo_.empty()
421 || !this->fields().empty()
422 || velocityOutput_->enableOutput()
423 || addProcessRank)
424 {
425 const auto numCells = gridGeometry().gridView().size(0);
426 const auto numDofs = numDofs_();
427
428 // get fields for all volume variables
429 if (!volVarScalarDataInfo_.empty())
430 volVarScalarData.resize(volVarScalarDataInfo_.size(), std::vector<Scalar>(numDofs));
431 if (!volVarVectorDataInfo_.empty())
432 volVarVectorData.resize(volVarVectorDataInfo_.size(), std::vector<VolVarsVector>(numDofs));
433
434 if (velocityOutput_->enableOutput())
435 {
436 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
437 {
438 if(isBox && dim == 1)
439 velocity[phaseIdx].resize(numCells);
440 else
441 velocity[phaseIdx].resize(numDofs);
442 }
443 }
444
445 // maybe allocate space for the process rank
446 if (addProcessRank) rank.resize(numCells);
447
448 for (const auto& element : elements(gridGeometry().gridView(), Dune::Partitions::interior))
449 {
450 const auto eIdxGlobal = gridGeometry().elementMapper().index(element);
451
452 auto fvGeometry = localView(gridGeometry());
453 auto elemVolVars = localView(gridVariables_.curGridVolVars());
454
455 // If velocity output is enabled we need to bind to the whole stencil
456 // otherwise element-local data is sufficient
457 if (velocityOutput_->enableOutput())
458 {
459 fvGeometry.bind(element);
460 elemVolVars.bind(element, fvGeometry, sol_);
461 }
462 else
463 {
464 fvGeometry.bindElement(element);
465 elemVolVars.bindElement(element, fvGeometry, sol_);
466 }
467
468 if (!volVarScalarDataInfo_.empty()
469 || !volVarVectorDataInfo_.empty())
470 {
471 for (auto&& scv : scvs(fvGeometry))
472 {
473 const auto dofIdxGlobal = scv.dofIndex();
474 const auto& volVars = elemVolVars[scv];
475
476 // get the scalar-valued data
477 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
478 volVarScalarData[i][dofIdxGlobal] = volVarScalarDataInfo_[i].get(volVars);
479
480 // get the vector-valued data
481 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
482 volVarVectorData[i][dofIdxGlobal] = volVarVectorDataInfo_[i].get(volVars);
483 }
484 }
485
486 // velocity output
487 if (velocityOutput_->enableOutput())
488 {
489 auto elemFluxVarsCache = localView(gridVariables_.gridFluxVarsCache());
490 elemFluxVarsCache.bind(element, fvGeometry, elemVolVars);
491
492 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
493 velocityOutput_->calculateVelocity(velocity[phaseIdx], element, fvGeometry, elemVolVars, elemFluxVarsCache, phaseIdx);
494 }
495
497 if (addProcessRank)
498 rank[eIdxGlobal] = static_cast<double>(gridGeometry().gridView().comm().rank());
499 }
500
504
505 // volume variables if any
506 if (isBox)
507 {
508 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
509 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().vertexMapper(), volVarScalarData[i],
510 volVarScalarDataInfo_[i].name, /*numComp*/1, /*codim*/dim, dm, this->precision()).get() );
511 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
512 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().vertexMapper(), volVarVectorData[i],
513 volVarVectorDataInfo_[i].name, /*numComp*/dimWorld, /*codim*/dim, dm, this->precision()).get() );
514 }
515 else
516 {
517 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
518 this->sequenceWriter().addCellData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), volVarScalarData[i],
519 volVarScalarDataInfo_[i].name, /*numComp*/1, /*codim*/0,dm, this->precision()).get() );
520 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
521 this->sequenceWriter().addCellData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), volVarVectorData[i],
522 volVarVectorDataInfo_[i].name, /*numComp*/dimWorld, /*codim*/0,dm, this->precision()).get() );
523 }
524
525 // the velocity field
526 if (velocityOutput_->enableOutput())
527 {
528 if (isBox && dim > 1)
529 {
530 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
531 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().vertexMapper(), velocity[phaseIdx],
532 "velocity_" + velocityOutput_->phaseName(phaseIdx) + " (m/s)",
533 /*numComp*/dimWorld, /*codim*/dim, dm, this->precision()).get() );
534 }
535 // cell-centered models
536 else
537 {
538 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
539 this->sequenceWriter().addCellData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), velocity[phaseIdx],
540 "velocity_" + velocityOutput_->phaseName(phaseIdx) + " (m/s)",
541 /*numComp*/dimWorld, /*codim*/0, dm, this->precision()).get() );
542 }
543 }
544
545 // the process rank
546 if (addProcessRank)
547 this->sequenceWriter().addCellData(Field(gridGeometry().gridView(), gridGeometry().elementMapper(), rank, "process rank", 1, 0).get());
548
549 // also register additional (non-standardized) user fields if any
550 for (auto&& field : this->fields())
551 {
552 if (field.codim() == 0)
553 this->sequenceWriter().addCellData(field.get());
554 else if (field.codim() == dim)
555 this->sequenceWriter().addVertexData(field.get());
556 else
557 DUNE_THROW(Dune::RangeError, "Cannot add wrongly sized vtk scalar field!");
558 }
559 }
560
564 this->sequenceWriter().write(time, type);
565
569 this->writer().clear();
570 }
571
573 void writeNonConforming_(double time, Dune::VTK::OutputType type) override
574 {
575 const Dune::VTK::DataMode dm = Dune::VTK::nonconforming;
576
577 if(!isBox)
578 DUNE_THROW(Dune::InvalidStateException, "Non-conforming output makes no sense for cell-centered schemes!");
579
583
584 // check the velocity output
585 bool enableVelocityOutput = getParamFromGroup<bool>(this->paramGroup(), "Vtk.AddVelocity");
586 if (enableVelocityOutput == true && !velocityOutput_->enableOutput())
587 std::cerr << "Warning! Velocity output was enabled in the input file"
588 << " but no velocity output policy was set for the VTK output module:"
589 << " There will be no velocity output."
590 << " Use the addVelocityOutput member function of the VTK output module." << std::endl;
591 using VelocityVector = typename VelocityOutput::VelocityVector;
592 std::vector<VelocityVector> velocity(velocityOutput_->numFluidPhases());
593
594 // process rank
595 static bool addProcessRank = getParamFromGroup<bool>(this->paramGroup(), "Vtk.AddProcessRank");
596 std::vector<double> rank;
597
598 // volume variable data (indexing: volvardata/element/localcorner)
599 using ScalarDataContainer = std::vector< std::vector<Scalar> >;
600 using VectorDataContainer = std::vector< std::vector<VolVarsVector> >;
601 std::vector< ScalarDataContainer > volVarScalarData;
602 std::vector< VectorDataContainer > volVarVectorData;
603
605 if (!volVarScalarDataInfo_.empty()
606 || !volVarVectorDataInfo_.empty()
607 || !this->fields().empty()
608 || velocityOutput_->enableOutput()
609 || addProcessRank)
610 {
611 const auto numCells = gridGeometry().gridView().size(0);
612 const auto numDofs = numDofs_();
613
614 // get fields for all volume variables
615 if (!volVarScalarDataInfo_.empty())
616 volVarScalarData.resize(volVarScalarDataInfo_.size(), ScalarDataContainer(numCells));
617 if (!volVarVectorDataInfo_.empty())
618 volVarVectorData.resize(volVarVectorDataInfo_.size(), VectorDataContainer(numCells));
619
620 if (velocityOutput_->enableOutput())
621 {
622 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
623 {
624 if(isBox && dim == 1)
625 velocity[phaseIdx].resize(numCells);
626 else
627 velocity[phaseIdx].resize(numDofs);
628 }
629 }
630
631 // maybe allocate space for the process rank
632 if (addProcessRank) rank.resize(numCells);
633
634 for (const auto& element : elements(gridGeometry().gridView(), Dune::Partitions::interior))
635 {
636 const auto eIdxGlobal = gridGeometry().elementMapper().index(element);
637 const auto numCorners = element.subEntities(dim);
638
639 auto fvGeometry = localView(gridGeometry());
640 auto elemVolVars = localView(gridVariables_.curGridVolVars());
641
642 // resize element-local data containers
643 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
644 volVarScalarData[i][eIdxGlobal].resize(numCorners);
645 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
646 volVarVectorData[i][eIdxGlobal].resize(numCorners);
647
648 // If velocity output is enabled we need to bind to the whole stencil
649 // otherwise element-local data is sufficient
650 if (velocityOutput_->enableOutput())
651 {
652 fvGeometry.bind(element);
653 elemVolVars.bind(element, fvGeometry, sol_);
654 }
655 else
656 {
657 fvGeometry.bindElement(element);
658 elemVolVars.bindElement(element, fvGeometry, sol_);
659 }
660
661 if (!volVarScalarDataInfo_.empty()
662 || !volVarVectorDataInfo_.empty())
663 {
664 for (auto&& scv : scvs(fvGeometry))
665 {
666 const auto& volVars = elemVolVars[scv];
667
668 // get the scalar-valued data
669 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
670 volVarScalarData[i][scv.elementIndex()][scv.localDofIndex()] = volVarScalarDataInfo_[i].get(volVars);
671
672 // get the vector-valued data
673 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
674 volVarVectorData[i][scv.elementIndex()][scv.localDofIndex()] = volVarVectorDataInfo_[i].get(volVars);
675 }
676 }
677
678 // velocity output
679 if (velocityOutput_->enableOutput())
680 {
681 auto elemFluxVarsCache = localView(gridVariables_.gridFluxVarsCache());
682 elemFluxVarsCache.bind(element, fvGeometry, elemVolVars);
683
684 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
685 velocityOutput_->calculateVelocity(velocity[phaseIdx], element, fvGeometry, elemVolVars, elemFluxVarsCache, phaseIdx);
686 }
687
689 if (addProcessRank)
690 rank[eIdxGlobal] = static_cast<double>(gridGeometry().gridView().comm().rank());
691 }
692
696
697 // volume variables if any
698 for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
699 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), volVarScalarData[i],
700 volVarScalarDataInfo_[i].name, /*numComp*/1, /*codim*/dim, /*nonconforming*/dm, this->precision()).get() );
701
702 for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
703 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), volVarVectorData[i],
704 volVarVectorDataInfo_[i].name, /*numComp*/dimWorld, /*codim*/dim, /*nonconforming*/dm, this->precision()).get() );
705
706 // the velocity field
707 if (velocityOutput_->enableOutput())
708 {
709 // node-wise velocities
710 if (dim > 1)
711 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
712 this->sequenceWriter().addVertexData( Field(gridGeometry().gridView(), gridGeometry().vertexMapper(), velocity[phaseIdx],
713 "velocity_" + velocityOutput_->phaseName(phaseIdx) + " (m/s)",
714 /*numComp*/dimWorld, /*codim*/dim, dm, this->precision()).get() );
715
716 // cell-wise velocities
717 else
718 for (int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
719 this->sequenceWriter().addCellData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), velocity[phaseIdx],
720 "velocity_" + velocityOutput_->phaseName(phaseIdx) + " (m/s)",
721 /*numComp*/dimWorld, /*codim*/0,dm, this->precision()).get());
722 }
723
724 // the process rank
725 if (addProcessRank)
726 this->sequenceWriter().addCellData( Field(gridGeometry().gridView(), gridGeometry().elementMapper(), rank, "process rank", 1, 0).get() );
727
728 // also register additional (non-standardized) user fields if any
729 for (auto&& field : this->fields())
730 {
731 if (field.codim() == 0)
732 this->sequenceWriter().addCellData(field.get());
733 else if (field.codim() == dim)
734 this->sequenceWriter().addVertexData(field.get());
735 else
736 DUNE_THROW(Dune::RangeError, "Cannot add wrongly sized vtk scalar field!");
737 }
738 }
739
743 this->sequenceWriter().write(time, type);
744
748 this->writer().clear();
749 }
750
752 std::size_t numDofs_() const { return dofCodim == dim ? gridGeometry().vertexMapper().size() : gridGeometry().elementMapper().size(); }
753
754 const GridVariables& gridVariables_;
755 const SolutionVector& sol_;
756
757 std::vector<VolVarScalarDataInfo> volVarScalarDataInfo_;
758 std::vector<VolVarVectorDataInfo> volVarVectorDataInfo_;
759
760 std::shared_ptr<VelocityOutput> velocityOutput_;
761};
762
763} // end namespace Dumux
764
765#endif
parameters.hh
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
method.hh
The available discretization methods in Dumux.
fieldtype.hh
Vtk field types available in Dumux.
function.hh
Dune style VTK functions.
format.hh
Formatting based on the fmt-library which implements std::format of C++20.
Dumux::localView
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:38
Dumux::Vtk::stringToPrecision
Precision stringToPrecision(std::string_view precisionName)
Maps a string (e.g. from input) to a Dune precision type.
Definition: precision.hh:39
Dumux::Vtk::FieldType
FieldType
Identifier for vtk field types.
Definition: fieldtype.hh:34
Dumux
Definition: adapt.hh:29
Dumux::PoreNetwork::Throat::numCorners
std::size_t numCorners(Shape shape)
Returns the number of corners of a given geometry.
Definition: throatproperties.hh:215
Dumux::VelocityOutput
Velocity output for implicit (porous media) models.
Definition: io/velocityoutput.hh:41
Dumux::VelocityOutput::VelocityVector
std::vector< Dune::FieldVector< Scalar, dimWorld > > VelocityVector
Definition: io/velocityoutput.hh:50
Dumux::VtkOutputModuleBase
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition: io/vtkoutputmodule.hh:61
Dumux::VtkOutputModuleBase::precision
Dumux::Vtk::Precision precision() const
Definition: io/vtkoutputmodule.hh:197
Dumux::VtkOutputModuleBase::paramGroup
const std::string & paramGroup() const
the parameter group for getting parameter from the parameter tree
Definition: io/vtkoutputmodule.hh:96
Dumux::VtkOutputModuleBase::write
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Definition: io/vtkoutputmodule.hh:173
Dumux::VtkOutputModuleBase::VtkOutputModuleBase
VtkOutputModuleBase(const GridGeometry &gridGeometry, const std::string &name, const std::string &paramGroup="", Dune::VTK::DataMode dm=Dune::VTK::conforming, bool verbose=true)
Definition: io/vtkoutputmodule.hh:75
Dumux::VtkOutputModuleBase::dataMode
Dune::VTK::DataMode dataMode() const
Definition: io/vtkoutputmodule.hh:196
Dumux::VtkOutputModuleBase::writer
Dune::VTKWriter< GridView > & writer()
Definition: io/vtkoutputmodule.hh:199
Dumux::VtkOutputModuleBase::sequenceWriter
Dune::VTKSequenceWriter< GridView > & sequenceWriter()
Definition: io/vtkoutputmodule.hh:200
Dumux::VtkOutputModuleBase::name
const std::string & name() const
Definition: io/vtkoutputmodule.hh:195
Dumux::VtkOutputModuleBase::fields
const std::vector< Field > & fields() const
Definition: io/vtkoutputmodule.hh:202
Dumux::VtkOutputModuleBase::~VtkOutputModuleBase
virtual ~VtkOutputModuleBase()=default
Dumux::VtkOutputModuleBase::addField
void addField(const Vector &v, const std::string &name, Dumux::Vtk::Precision precision, Vtk::FieldType fieldType=Vtk::FieldType::automatic)
Add a scalar or vector valued vtk field.
Definition: io/vtkoutputmodule.hh:125
Dumux::VtkOutputModuleBase::verbose
bool verbose() const
Definition: io/vtkoutputmodule.hh:194
Dumux::VtkOutputModuleBase::addField
void addField(const Vector &v, const std::string &name, Vtk::FieldType fieldType=Vtk::FieldType::automatic)
Add a scalar or vector valued vtk field.
Definition: io/vtkoutputmodule.hh:109
Dumux::VtkOutputModuleBase::gridGeometry
const GridGeometry & gridGeometry() const
Definition: io/vtkoutputmodule.hh:192
Dumux::VtkOutputModuleBase::FieldType
export field type
Definition: io/vtkoutputmodule.hh:69
Dumux::VtkOutputModule
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition: io/vtkoutputmodule.hh:307
Dumux::VtkOutputModule::problem
const auto & problem() const
Definition: io/vtkoutputmodule.hh:385
Dumux::VtkOutputModule::gridVariables
const GridVariables & gridVariables() const
Definition: io/vtkoutputmodule.hh:386
Dumux::VtkOutputModule::addVolumeVariable
void addVolumeVariable(std::function< VolVarsVector(const VolumeVariables &)> &&f, const std::string &name)
Definition: io/vtkoutputmodule.hh:377
Dumux::VtkOutputModule::VolumeVariables
VV VolumeVariables
export type of the volume variables for the outputfields
Definition: io/vtkoutputmodule.hh:335
Dumux::VtkOutputModule::addVolumeVariable
void addVolumeVariable(std::function< Scalar(const VolumeVariables &)> &&f, const std::string &name)
Definition: io/vtkoutputmodule.hh:366
Dumux::VtkOutputModule::addVelocityOutput
void addVelocityOutput(std::shared_ptr< VelocityOutputType > velocityOutput)
Add a velocity output policy.
Definition: io/vtkoutputmodule.hh:360
Dumux::VtkOutputModule::velocityOutput
const VelocityOutput & velocityOutput() const
Definition: io/vtkoutputmodule.hh:394
Dumux::VtkOutputModule::sol
const SolutionVector & sol() const
Definition: io/vtkoutputmodule.hh:388
Dumux::VtkOutputModule::gridGeometry
const GridGeometry & gridGeometry() const
Definition: io/vtkoutputmodule.hh:387
Dumux::VtkOutputModule::VtkOutputModule
VtkOutputModule(const GridVariables &gridVariables, const SolutionVector &sol, const std::string &name, const std::string &paramGroup="", Dune::VTK::DataMode dm=Dune::VTK::conforming, bool verbose=true)
Definition: io/vtkoutputmodule.hh:337
Dumux::VtkOutputModule::volVarScalarDataInfo
const std::vector< VolVarScalarDataInfo > & volVarScalarDataInfo() const
Definition: io/vtkoutputmodule.hh:390
Dumux::VtkOutputModule::volVarVectorDataInfo
const std::vector< VolVarVectorDataInfo > & volVarVectorDataInfo() const
Definition: io/vtkoutputmodule.hh:391
velocityoutput.hh
Velocity output for porous media models.
Generated by doxygen 1.9.3