27#ifndef OPM_VTK_TRACER_MODULE_HPP
28#define OPM_VTK_TRACER_MODULE_HPP
30#include <dune/common/fvector.hh>
32#include <opm/models/blackoil/blackoilproperties.hh>
33#include <opm/models/io/baseoutputmodule.hh>
34#include <opm/models/io/vtkmultiwriter.hh>
35#include <opm/models/utils/parametersystem.hh>
36#include <opm/models/utils/propertysystem.hh>
49template<
class TypeTag,
class MyTypeTag>
51 using type = UndefinedProperty;
55template<
class TypeTag>
57 static constexpr bool value =
false;
68 template <
class TypeTag>
71 using ParentType = BaseOutputModule<TypeTag>;
73 using Simulator = GetPropType<TypeTag, Properties::Simulator>;
74 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
75 using Evaluation = GetPropType<TypeTag, Properties::Evaluation>;
76 using ElementContext = GetPropType<TypeTag, Properties::ElementContext>;
78 using GridView = GetPropType<TypeTag, Properties::GridView>;
79 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
81 static constexpr int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
82 using VtkMultiWriter = ::Opm::VtkMultiWriter<GridView, vtkFormat>;
85 using ScalarBuffer =
typename ParentType::ScalarBuffer;
89 : ParentType(simulator)
98 Parameters::registerParam<TypeTag, Properties::VtkWriteTracerConcentration>
99 (
"Include the tracer concentration in the VTK output files");
108 if (eclTracerConcentrationOutput_()){
109 const auto& tracerModel = this->simulator_.problem().tracerModel();
110 eclTracerConcentration_.resize(tracerModel.numTracers());
111 for (std::size_t tracerIdx = 0; tracerIdx < eclTracerConcentration_.size(); ++tracerIdx) {
113 this->resizeScalarBuffer_(eclTracerConcentration_[tracerIdx]);
125 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
128 const auto& tracerModel = elemCtx.problem().tracerModel();
130 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
131 unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
133 if (eclTracerConcentrationOutput_()){
134 for (std::size_t tracerIdx = 0; tracerIdx < eclTracerConcentration_.size(); ++tracerIdx) {
135 eclTracerConcentration_[tracerIdx][globalDofIdx] = tracerModel.tracerConcentration(tracerIdx, globalDofIdx);
146 VtkMultiWriter *vtkWriter =
dynamic_cast<VtkMultiWriter*
>(&baseWriter);
150 if (eclTracerConcentrationOutput_()){
151 const auto& tracerModel = this->simulator_.problem().tracerModel();
152 for (std::size_t tracerIdx = 0; tracerIdx < eclTracerConcentration_.size(); ++tracerIdx) {
153 const std::string tmp =
"tracerConcentration_" + tracerModel.name(tracerIdx);
154 this->commitScalarBuffer_(baseWriter,tmp.c_str(), eclTracerConcentration_[tracerIdx]);
163 static bool eclTracerConcentrationOutput_()
165 static bool val = Parameters::get<TypeTag, Properties::VtkWriteTracerConcentration>();
170 std::vector<ScalarBuffer> eclTracerConcentration_;
VTK output module for the tracer model's parameters.
Definition: VtkTracerModule.hpp:70
static void registerParameters()
Register all run-time parameters for the tracer VTK output module.
Definition: VtkTracerModule.hpp:96
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: VtkTracerModule.hpp:144
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: VtkTracerModule.hpp:123
VtkTracerModule(const Simulator &simulator)
Definition: VtkTracerModule.hpp:88
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: VtkTracerModule.hpp:106
Definition: AluGridVanguard.hpp:57
Definition: BlackoilPhases.hpp:27
Definition: VtkTracerModule.hpp:45
Definition: VtkTracerModule.hpp:50
UndefinedProperty type
Definition: VtkTracerModule.hpp:51