13#ifndef DUMUX_PNM_VTK_OUTPUT_MODULE_HH 
   14#define DUMUX_PNM_VTK_OUTPUT_MODULE_HH 
   25template<
class Gr
idVariables, 
class FluxVariables, 
class SolutionVector>
 
   29    using Scalar = 
typename GridVariables::Scalar;
 
   30    using FluxVarsCache = 
typename GridVariables::GridFluxVariablesCache::FluxVariablesCache;
 
   32    struct ThroatFluxDataInfo { std::function<Scalar(
const FluxVariables&, 
const FluxVarsCache&)> get; std::string 
name; };
 
   38                    const SolutionVector& 
sol,
 
   39                    const std::string& 
name,
 
   41                    Dune::VTK::DataMode dm = Dune::VTK::conforming,
 
   45        if constexpr (GridVariables::VolumeVariables::numFluidPhases() >= 1)
 
 
   56    void addFluxVariable(std::function<Scalar(
const FluxVariables&, 
const FluxVarsCache&)>&& f, 
const std::string& 
name)
 
   58        throatFluxDataInfo_.push_back(ThroatFluxDataInfo{f, 
name});
 
   59        const auto numElems = 
problem().gridGeometry().gridView().size(0);
 
   60        throatFluxData_.push_back(std::vector<Scalar>(numElems));
 
 
   65    void write(
double time, Dune::VTK::OutputType type = Dune::VTK::ascii)
 
   67        const auto gridView = 
problem().gridGeometry().gridView();
 
   68        const auto numElems = gridView.size(0);
 
   71        for (
auto& data : throatFluxData_)
 
   72            data.resize(numElems);
 
   78        for (
const auto& element : elements(gridView, Dune::Partitions::interior))
 
   80            const auto eIdx = 
problem().gridGeometry().elementMapper().index(element);
 
   83            fvElementGeometry.bindElement(element);
 
   84            elemVolVars.bind(element, fvElementGeometry, this->
sol());
 
   85            elemFluxVarsCache.bind(element, fvElementGeometry, elemVolVars);
 
   88            std::size_t dataIdx = 0;
 
   89            for (
auto&& scvf : scvfs(fvElementGeometry))
 
   93                    FluxVariables fluxVars;
 
   94                    fluxVars.init(
problem(), element, fvElementGeometry, elemVolVars, scvf, elemFluxVarsCache);
 
   97                    for(
auto& data : throatFluxData_)
 
   98                        data[eIdx] = throatFluxDataInfo_[dataIdx++].get(fluxVars, elemFluxVarsCache[scvf]);
 
  107        auto clearAndShrink = [] (
auto& data)
 
  110            data.shrink_to_fit();
 
  113        for (
auto& data : throatFluxData_)
 
  114            clearAndShrink(data);
 
 
  121    std::vector<ThroatFluxDataInfo> throatFluxDataInfo_;
 
  122    std::list<std::vector<Scalar>> throatFluxData_;
 
 
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Gather and process all required data and write them to a vtk file.
Definition pnmvtkoutputmodule.hh:65
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)
The constructor.
Definition pnmvtkoutputmodule.hh:37
void addFluxVariable(std::function< Scalar(const FluxVariables &, const FluxVarsCache &)> &&f, const std::string &name)
Definition pnmvtkoutputmodule.hh:56
const auto & problem() const
Return a reference to the problem.
Definition pnmvtkoutputmodule.hh:118
const std::string & paramGroup() const
Definition io/vtkoutputmodule.hh:79
void write(double time, Dune::VTK::OutputType type=Dune::VTK::ascii)
Definition io/vtkoutputmodule.hh:182
const std::string & name() const
Definition io/vtkoutputmodule.hh:204
bool verbose() const
Definition io/vtkoutputmodule.hh:203
void addField(const Vector &v, const std::string &name, Vtk::FieldType fieldType=Vtk::FieldType::automatic)
Definition io/vtkoutputmodule.hh:92
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition io/vtkoutputmodule.hh:343
const auto & problem() const
Definition io/vtkoutputmodule.hh:460
const GridVariables & gridVariables() const
Definition io/vtkoutputmodule.hh:461
void addVelocityOutput(std::shared_ptr< VelocityOutputType > velocityOutput)
Add a velocity output policy.
Definition io/vtkoutputmodule.hh:401
VelocityOutputType VelocityOutput
Definition io/vtkoutputmodule.hh:468
const SolutionVector & sol() const
Definition io/vtkoutputmodule.hh:463
const GridGeometry & gridGeometry() const
Definition io/vtkoutputmodule.hh:462
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition localview.hh:26
A VTK output module to simplify writing dumux simulation data to VTK format.
Definition discretization/porenetwork/fvelementgeometry.hh:24
Velocity output for pore-network models.