12#ifndef DUMUX_IO_VTK_OUTPUT_MODULE_HH 
   13#define DUMUX_IO_VTK_OUTPUT_MODULE_HH 
   19#include <dune/common/timer.hh> 
   20#include <dune/common/fvector.hh> 
   21#include <dune/common/typetraits.hh> 
   23#include <dune/geometry/type.hh> 
   24#include <dune/geometry/multilineargeometry.hh> 
   26#include <dune/grid/common/mcmgmapper.hh> 
   27#include <dune/grid/common/partitionset.hh> 
   28#include <dune/grid/io/file/vtk/vtkwriter.hh> 
   29#include <dune/grid/io/file/vtk/vtksequencewriter.hh> 
   30#include <dune/grid/common/partitionset.hh> 
   47template<
class Gr
idGeometry>
 
   50    using GridView = 
typename GridGeometry::GridView;
 
   51    static constexpr int dim = GridView::dimension;
 
   58                        const std::string& 
name,
 
   60                        Dune::VTK::DataMode dm = Dune::VTK::conforming,
 
   71        writer_ = std::make_shared<Dune::VTKWriter<GridView>>(
gridGeometry.gridView(), dm, coordPrecision);
 
   72        sequenceWriter_ = std::make_unique<Dune::VTKSequenceWriter<GridView>>(writer_, 
name);
 
 
   80    { 
return paramGroup_; }
 
 
   91    template<
typename Vector>
 
   93                  const std::string& 
name,
 
 
  107    template<
typename Vector>
 
  109                  const std::string& 
name,
 
  114        const auto nComp = getNumberOfComponents_(v);
 
  116        const auto numElemDofs = 
gridGeometry().elementMapper().size();
 
  117        const auto numVertexDofs = 
gridGeometry().vertexMapper().size();
 
  122            if(numElemDofs == numVertexDofs)
 
  123                DUNE_THROW(Dune::InvalidStateException, 
"Automatic deduction of FieldType failed. Please explicitly specify FieldType::element or FieldType::vertex.");
 
  125            if(v.size() == numElemDofs)
 
  127            else if(v.size() == numVertexDofs)
 
  130                DUNE_THROW(Dune::RangeError, 
"Size mismatch of added field!");
 
  136                if(v.size() != numElemDofs)
 
  137                    DUNE_THROW(Dune::RangeError, 
"Size mismatch of added field!");
 
  140                if(v.size() != numVertexDofs)
 
  141                    DUNE_THROW(Dune::RangeError, 
"Size mismatch of added field!");
 
 
  159        for (
auto i = 0UL; i < fields_.size(); ++i)
 
  161            if (fields_[i].
name() == field.name() && fields_[i].codim() == field.codim())
 
  163                fields_[i] = std::move(field);
 
  164                std::cout << Fmt::format(
 
  165                    "VtkOutputModule: Replaced field \"{}\" (codim {}). " 
  166                    "A field by the same name & codim had already been registered previously.\n",
 
  167                    field.name(), field.codim()
 
  174        fields_.push_back(std::move(field));
 
 
  182    void write(
double time, Dune::VTK::OutputType type = Dune::VTK::ascii)
 
  187        if (dm_ == Dune::VTK::conforming)
 
  189        else if (dm_ == Dune::VTK::nonconforming)
 
  192            DUNE_THROW(Dune::NotImplemented, 
"Output for provided VTK data mode");
 
  197            std::cout << Fmt::format(
"Writing output for problem \"{}\". Took {:.2g} seconds.\n", name_, timer.elapsed());
 
 
  204    const std::string& 
name()
 const { 
return name_; }
 
  205    Dune::VTK::DataMode 
dataMode()
 const { 
return dm_; }
 
  206    Dumux::Vtk::Precision 
precision()
 const { 
return precision_; }
 
  208    Dune::VTKWriter<GridView>& 
writer() { 
return *writer_; }
 
  211    const std::vector<Field>& 
fields()
 const { 
return fields_; }
 
  245        std::vector<int> rank;
 
  248        if (!fields_.empty() || addProcessRank_)
 
  250            const auto numCells = gridGeometry_.gridView().size(0);
 
  255                rank.resize(numCells);
 
  257                for (
const auto& element : elements(gridGeometry_.gridView(), Dune::Partitions::interior))
 
  259                    const auto eIdxGlobal = gridGeometry_.elementMapper().index(element);
 
  260                    rank[eIdxGlobal] = gridGeometry_.gridView().comm().rank();
 
  270                this->
addCellData(
Field(gridGeometry_.gridView(), gridGeometry_.elementMapper(), rank, 
"process rank", 1, 0));
 
  273            for (
auto&& field : fields_)
 
  275                if (field.codim() == 0)
 
  277                else if (field.codim() == dim)
 
  280                    DUNE_THROW(Dune::RangeError, 
"Cannot add wrongly sized vtk scalar field!");
 
  294        this->addedCellData_.clear();
 
  295        this->addedVertexData_.clear();
 
 
  301        DUNE_THROW(Dune::NotImplemented, 
"Non-conforming VTK output");
 
 
  305    template<
class Vector>
 
  306    std::size_t getNumberOfComponents_(
const Vector& v)
 
  308        if constexpr (Dune::IsIndexable<decltype(std::declval<Vector>()[0])>::value)
 
  314    const GridGeometry& gridGeometry_;
 
  316    const std::string paramGroup_;
 
  317    Dune::VTK::DataMode dm_;
 
  319    Dumux::Vtk::Precision precision_;
 
  321    std::shared_ptr<Dune::VTKWriter<GridView>> writer_;
 
  322    std::unique_ptr<Dune::VTKSequenceWriter<GridView>> sequenceWriter_;
 
  324    std::vector<Field> fields_; 
 
  326    bool addProcessRank_ = 
true;
 
 
  341template<
class Gr
idVariables, 
class SolutionVector>
 
  345    using GridGeometry = 
typename GridVariables::GridGeometry;
 
  347    using VV = 
typename GridVariables::VolumeVariables;
 
  348    using Scalar = 
typename GridVariables::Scalar;
 
  350    using GridView = 
typename GridGeometry::GridView;
 
  353        dim = GridView::dimension,
 
  354        dimWorld = GridView::dimensionworld
 
  357    using Element = 
typename GridView::template Codim<0>::Entity;
 
  358    using VolVarsVector = Dune::FieldVector<Scalar, dimWorld>;
 
  364    struct VolVarScalarDataInfo { std::function<Scalar(
const VV&)> get; std::string 
name; Dumux::Vtk::Precision precision_; };
 
  365    struct VolVarVectorDataInfo { std::function<VolVarsVector(
const VV&)> get; std::string 
name; Dumux::Vtk::Precision precision_; };
 
  376                    const SolutionVector& 
sol,
 
  377                    const std::string& 
name,
 
  379                    Dune::VTK::DataMode dm = Dune::VTK::conforming,
 
  384    , velocityOutput_(std::make_shared<VelocityOutputType>())
 
 
  408                                                const std::string& 
name)
 
  412        for (
auto i = 0UL; i < volVarScalarDataInfo_.size(); ++i)
 
  414            if (volVarScalarDataInfo_[i].
name == 
name)
 
  416                volVarScalarDataInfo_[i] = VolVarScalarDataInfo{f, 
name, this->
precision()};
 
  417                std::cout << Fmt::format(
 
  418                    "VtkOutputModule: Replaced volume variable output \"{}\". " 
  419                    "A field by the same name had already been registered previously.\n",
 
  427        volVarScalarDataInfo_.push_back(VolVarScalarDataInfo{f, 
name, this->
precision()});
 
 
  434    template<
class VVV = VolVarsVector, 
typename std::enable_if_t<(VVV::dimension > 1), 
int> = 0>
 
  436                                                       const std::string& 
name)
 
  440        for (
auto i = 0UL; i < volVarVectorDataInfo_.size(); ++i)
 
  442            if (volVarVectorDataInfo_[i].
name == 
name)
 
  444                volVarVectorDataInfo_[i] = VolVarVectorDataInfo{f, 
name, this->
precision()};
 
  445                std::cout << Fmt::format(
 
  446                    "VtkOutputModule: Replaced volume variable output \"{}\". " 
  447                    "A field by the same name had already been registered previously.\n",
 
  455        volVarVectorDataInfo_.push_back(VolVarVectorDataInfo{f, 
name, this->
precision()});
 
 
  460    const auto& 
problem()
 const { 
return gridVariables_.curGridVolVars().problem(); }
 
  462    const GridGeometry& 
gridGeometry()
 const { 
return gridVariables_.gridGeometry(); }
 
  463    const SolutionVector& 
sol()
 const { 
return sol_; }
 
  476        const Dune::VTK::DataMode dm = Dune::VTK::conforming;
 
  483        std::vector<VelocityVector> velocity(velocityOutput_->numFluidPhases());
 
  486        std::vector<double> rank;
 
  489        std::vector<std::vector<Scalar>> volVarScalarData;
 
  490        std::vector<std::vector<VolVarsVector>> volVarVectorData;
 
  493        if (!volVarScalarDataInfo_.empty()
 
  494            || !volVarVectorDataInfo_.empty()
 
  495            || !this->fields().empty()
 
  496            || velocityOutput_->enableOutput()
 
  499            const auto numCells = 
gridGeometry().gridView().size(0);
 
  500            const auto numDofs = numDofs_();
 
  503            if (!volVarScalarDataInfo_.empty())
 
  504                volVarScalarData.resize(volVarScalarDataInfo_.size(), std::vector<Scalar>(numDofs));
 
  505            if (!volVarVectorDataInfo_.empty())
 
  506                volVarVectorData.resize(volVarVectorDataInfo_.size(), std::vector<VolVarsVector>(numDofs));
 
  508            if (velocityOutput_->enableOutput())
 
  510                for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  513                        velocity[phaseIdx].resize(numCells);
 
  515                        velocity[phaseIdx].resize(numDofs);
 
  518                        if(isBox && dim == 1)
 
  519                            velocity[phaseIdx].resize(numCells);
 
  521                            velocity[phaseIdx].resize(numDofs);
 
  527            if (addProcessRank_) rank.resize(numCells);
 
  530            auto elemVolVars = 
localView(gridVariables_.curGridVolVars());
 
  531            for (
const auto& element : elements(
gridGeometry().gridView(), Dune::Partitions::interior))
 
  533                const auto eIdxGlobal = 
gridGeometry().elementMapper().index(element);
 
  536                if (velocityOutput_->enableOutput())
 
  538                    fvGeometry.bind(element);
 
  539                    elemVolVars.bind(element, fvGeometry, sol_);
 
  543                    fvGeometry.bindElement(element);
 
  544                    elemVolVars.bindElement(element, fvGeometry, sol_);
 
  547                if (!volVarScalarDataInfo_.empty() || !volVarVectorDataInfo_.empty())
 
  549                    for (
const auto& scv : scvs(fvGeometry))
 
  551                        const auto dofIdxGlobal = scv.dofIndex();
 
  552                        const auto& volVars = elemVolVars[scv];
 
  555                        for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  556                            volVarScalarData[i][dofIdxGlobal] = volVarScalarDataInfo_[i].get(volVars);
 
  559                        for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  560                            volVarVectorData[i][dofIdxGlobal] = volVarVectorDataInfo_[i].get(volVars);
 
  565                if (velocityOutput_->enableOutput())
 
  567                    const auto elemFluxVarsCache = 
localView(gridVariables_.gridFluxVarsCache()).bind(element, fvGeometry, elemVolVars);
 
  569                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  570                        velocityOutput_->calculateVelocity(velocity[phaseIdx], element, fvGeometry, elemVolVars, elemFluxVarsCache, phaseIdx);
 
  575                    rank[eIdxGlobal] = 
static_cast<double>(
gridGeometry().gridView().comm().rank());
 
  583            if constexpr (isBox || isPQ1Bubble)
 
  585                for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  587                                          volVarScalarDataInfo_[i].
name, 1, dim, dm, this->
precision()) );
 
  588                for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  590                                          volVarVectorDataInfo_[i].
name, dimWorld, dim, dm, this->
precision()) );
 
  592                if constexpr (isPQ1Bubble)
 
  594                    for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  596                                            volVarScalarDataInfo_[i].
name, 1, 0,dm, this->
precision()) );
 
  597                    for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  599                                            volVarVectorDataInfo_[i].
name, dimWorld, 0,dm, this->
precision()) );
 
  605                for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  607                                        volVarScalarDataInfo_[i].
name, 1, 0,dm, this->
precision()) );
 
  608                for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  610                                        volVarVectorDataInfo_[i].
name, dimWorld, 0,dm, this->
precision()) );
 
  614            if (velocityOutput_->enableOutput())
 
  619                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  621                                              "velocity_" + velocityOutput_->phaseName(phaseIdx) + 
" (m/s)",
 
  622                                              dimWorld, dim, dm, this->precision()) );
 
  627                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  629                                            "velocity_" + velocityOutput_->phaseName(phaseIdx) + 
" (m/s)",
 
  630                                            dimWorld, 0, dm, this->precision()) );
 
  639            for (
auto&& field : this->
fields())
 
  641                if (field.codim() == 0)
 
  643                else if (field.codim() == dim)
 
  646                    DUNE_THROW(Dune::RangeError, 
"Cannot add wrongly sized vtk scalar field!");
 
  660        this->addedCellData_.clear();
 
  661        this->addedVertexData_.clear();
 
 
  667        const Dune::VTK::DataMode dm = Dune::VTK::nonconforming;
 
  670        if(!isBox && !isDiamond)
 
  671            DUNE_THROW(Dune::NotImplemented,
 
  672                "Non-conforming output for discretization scheme " << GridGeometry::discMethod
 
  680        if (enableVelocityOutput_ && !velocityOutput_->enableOutput())
 
  681            std::cerr << 
"Warning! Velocity output was enabled in the input file" 
  682                    << 
" but no velocity output policy was set for the VTK output module:" 
  683                    << 
" There will be no velocity output." 
  684                    << 
" Use the addVelocityOutput member function of the VTK output module." << std::endl;
 
  686        std::vector<VelocityVector> velocity(velocityOutput_->numFluidPhases());
 
  689        std::vector<double> rank;
 
  692        using ScalarDataContainer = std::vector< std::vector<Scalar> >;
 
  693        using VectorDataContainer = std::vector< std::vector<VolVarsVector> >;
 
  694        std::vector< ScalarDataContainer > volVarScalarData;
 
  695        std::vector< VectorDataContainer > volVarVectorData;
 
  698        if (!volVarScalarDataInfo_.empty()
 
  699            || !volVarVectorDataInfo_.empty()
 
  700            || !this->fields().empty()
 
  701            || velocityOutput_->enableOutput()
 
  704            const auto numCells = 
gridGeometry().gridView().size(0);
 
  705            const auto outputSize = numDofs_();
 
  708            if (!volVarScalarDataInfo_.empty())
 
  709                volVarScalarData.resize(volVarScalarDataInfo_.size(), ScalarDataContainer(numCells));
 
  710            if (!volVarVectorDataInfo_.empty())
 
  711                volVarVectorData.resize(volVarVectorDataInfo_.size(), VectorDataContainer(numCells));
 
  713            if (velocityOutput_->enableOutput())
 
  715                for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  717                    if((isBox && dim == 1) || isDiamond)
 
  718                        velocity[phaseIdx].resize(numCells);
 
  720                        velocity[phaseIdx].resize(outputSize);
 
  725            if (addProcessRank_) rank.resize(numCells);
 
  729            auto elemVolVars = 
localView(gridVariables_.curGridVolVars());
 
  730            for (
const auto& element : elements(
gridGeometry().gridView(), Dune::Partitions::interior))
 
  732                const auto eIdxGlobal = 
gridGeometry().elementMapper().index(element);
 
  735                if (velocityOutput_->enableOutput())
 
  737                    fvGeometry.bind(element);
 
  738                    elemVolVars.bind(element, fvGeometry, sol_);
 
  742                    fvGeometry.bindElement(element);
 
  743                    elemVolVars.bindElement(element, fvGeometry, sol_);
 
  746                const auto numLocalDofs = fvGeometry.numScv();
 
  748                for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  749                    volVarScalarData[i][eIdxGlobal].resize(numLocalDofs);
 
  750                for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  751                    volVarVectorData[i][eIdxGlobal].resize(numLocalDofs);
 
  753                if (!volVarScalarDataInfo_.empty() || !volVarVectorDataInfo_.empty())
 
  755                    for (
const auto& scv : scvs(fvGeometry))
 
  757                        const auto& volVars = elemVolVars[scv];
 
  760                        for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  761                            volVarScalarData[i][eIdxGlobal][scv.localDofIndex()] = volVarScalarDataInfo_[i].get(volVars);
 
  764                        for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  765                            volVarVectorData[i][eIdxGlobal][scv.localDofIndex()] = volVarVectorDataInfo_[i].get(volVars);
 
  770                if (velocityOutput_->enableOutput())
 
  772                    const auto elemFluxVarsCache = 
localView(gridVariables_.gridFluxVarsCache()).bind(element, fvGeometry, elemVolVars);
 
  773                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  774                        velocityOutput_->calculateVelocity(velocity[phaseIdx], element, fvGeometry, elemVolVars, elemFluxVarsCache, phaseIdx);
 
  779                    rank[eIdxGlobal] = 
static_cast<double>(
gridGeometry().gridView().comm().rank());
 
  787            static constexpr int dofLocCodim = isDiamond ? 1 : dim;
 
  788            for (std::size_t i = 0; i < volVarScalarDataInfo_.size(); ++i)
 
  791                    volVarScalarData[i], volVarScalarDataInfo_[i].
name,
 
  795            for (std::size_t i = 0; i < volVarVectorDataInfo_.size(); ++i)
 
  798                    volVarVectorData[i], volVarVectorDataInfo_[i].
name,
 
  799                    dimWorld, dofLocCodim, dm, this->
precision()
 
  803            if (velocityOutput_->enableOutput())
 
  806                if (dim > 1 && !isDiamond)
 
  807                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  810                            "velocity_" + velocityOutput_->phaseName(phaseIdx) + 
" (m/s)",
 
  811                            dimWorld, dofLocCodim, dm, this->precision()
 
  816                    for (
int phaseIdx = 0; phaseIdx < velocityOutput_->numFluidPhases(); ++phaseIdx)
 
  819                            "velocity_" + velocityOutput_->phaseName(phaseIdx) + 
" (m/s)",
 
  820                            dimWorld, 0, dm, this->precision()
 
  829                rank, 
"process rank", 1, 0
 
  833        for (
const auto& field : this->
fields())
 
  835            if (field.codim() == 0)
 
  837            else if (field.codim() == dim || field.codim() == 1)
 
  840                DUNE_THROW(Dune::RangeError, 
"Cannot add wrongly sized vtk scalar field!");
 
  853        this->addedCellData_.clear();
 
  854        this->addedVertexData_.clear();
 
 
  858    std::size_t numDofs_()
 const 
  862        if constexpr (isBox || isDiamond || isPQ1Bubble)
 
  868    const GridVariables& gridVariables_;
 
  869    const SolutionVector& sol_;
 
  871    std::vector<VolVarScalarDataInfo> volVarScalarDataInfo_; 
 
  872    std::vector<VolVarVectorDataInfo> volVarVectorDataInfo_; 
 
  874    std::shared_ptr<VelocityOutput> velocityOutput_; 
 
  875    bool enableVelocityOutput_ = 
false;
 
  876    bool addProcessRank_ = 
true;
 
 
Velocity output for implicit (porous media) models.
Definition io/velocityoutput.hh:29
@ vertex
Definition io/velocityoutput.hh:45
@ automatic
Definition io/velocityoutput.hh:45
@ element
Definition io/velocityoutput.hh:45
std::vector< Dune::FieldVector< Scalar, dimWorld > > VelocityVector
Definition io/velocityoutput.hh:38
Dumux::Vtk::Precision precision() const
Definition io/vtkoutputmodule.hh:206
const std::string & paramGroup() const
Definition io/vtkoutputmodule.hh:79
void addCellData(const Field &field)
Definition io/vtkoutputmodule.hh:213
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Definition io/vtkoutputmodule.hh:182
VtkOutputModuleBase(const GridGeometry &gridGeometry, const std::string &name, const std::string ¶mGroup="", Dune::VTK::DataMode dm=Dune::VTK::conforming, bool verbose=true)
Definition io/vtkoutputmodule.hh:57
Dune::VTK::DataMode dataMode() const
Definition io/vtkoutputmodule.hh:205
Dune::VTKWriter< GridView > & writer()
Definition io/vtkoutputmodule.hh:208
Dune::VTKSequenceWriter< GridView > & sequenceWriter()
Definition io/vtkoutputmodule.hh:209
const std::string & name() const
Definition io/vtkoutputmodule.hh:204
const std::vector< Field > & fields() const
Definition io/vtkoutputmodule.hh:211
virtual ~VtkOutputModuleBase()=default
virtual void writeNonConforming_(double time, Dune::VTK::OutputType type)
Assembles the fields and adds them to the writer (nonconforming output)
Definition io/vtkoutputmodule.hh:299
std::vector< std::string > addedCellData_
Definition io/vtkoutputmodule.hh:233
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:108
bool verbose() const
Definition io/vtkoutputmodule.hh:203
Vtk::template Field< GridView > Field
the type of Field that can be added to this writer
Definition io/vtkoutputmodule.hh:55
std::vector< std::string > addedVertexData_
Definition io/vtkoutputmodule.hh:234
virtual void writeConforming_(double time, Dune::VTK::OutputType type)
Assembles the fields and adds them to the writer (conforming output)
Definition io/vtkoutputmodule.hh:238
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:92
void addVertexData(const Field &field)
Definition io/vtkoutputmodule.hh:222
const GridGeometry & gridGeometry() const
Definition io/vtkoutputmodule.hh:201
void addField(Field &&field)
Add a scalar or vector valued vtk field.
Definition io/vtkoutputmodule.hh:155
const auto & problem() const
Definition io/vtkoutputmodule.hh:460
const GridVariables & gridVariables() const
Definition io/vtkoutputmodule.hh:461
void addVolumeVariable(std::function< VolVarsVector(const VolumeVariables &)> &&f, const std::string &name)
Definition io/vtkoutputmodule.hh:435
VV VolumeVariables
export type of the volume variables for the outputfields
Definition io/vtkoutputmodule.hh:373
void addVolumeVariable(std::function< Scalar(const VolumeVariables &)> &&f, const std::string &name)
Definition io/vtkoutputmodule.hh:407
void addVelocityOutput(std::shared_ptr< VelocityOutputType > velocityOutput)
Add a velocity output policy.
Definition io/vtkoutputmodule.hh:401
const VelocityOutput & velocityOutput() const
Definition io/vtkoutputmodule.hh:469
VelocityOutputType VelocityOutput
Definition io/vtkoutputmodule.hh:468
Vtk::template Field< GridView > Field
the type of Field that can be added to this writer
Definition io/vtkoutputmodule.hh:371
void writeNonConforming_(double time, Dune::VTK::OutputType type) override
Assembles the fields and adds them to the writer (nonconforming output)
Definition io/vtkoutputmodule.hh:665
const SolutionVector & sol() const
Definition io/vtkoutputmodule.hh:463
const GridGeometry & gridGeometry() const
Definition io/vtkoutputmodule.hh:462
void writeConforming_(double time, Dune::VTK::OutputType type) override
Assembles the fields and adds them to the writer (conforming output)
Definition io/vtkoutputmodule.hh:474
VtkOutputModule(const GridVariables &gridVariables, const SolutionVector &sol, const std::string &name, const std::string ¶mGroup="", Dune::VTK::DataMode dm=Dune::VTK::conforming, bool verbose=true)
Definition io/vtkoutputmodule.hh:375
const std::vector< VolVarScalarDataInfo > & volVarScalarDataInfo() const
Definition io/vtkoutputmodule.hh:465
const std::vector< VolVarVectorDataInfo > & volVarVectorDataInfo() const
Definition io/vtkoutputmodule.hh:466
Vtk field types available in Dumux.
Dune style VTK functions.
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
T getParamFromGroup(Args &&... args)
A free function to get a parameter from the parameter tree singleton with a model group.
Definition parameters.hh:149
Default velocity output policy for porous media models.
The available discretization methods in Dumux.
constexpr FCDiamond fcdiamond
Definition method.hh:152
constexpr Box box
Definition method.hh:147
constexpr PQ1Bubble pq1bubble
Definition method.hh:148
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.