27#ifndef EWOMS_VTK_BLACK_OIL_POLYMER_MODULE_HH
28#define EWOMS_VTK_BLACK_OIL_POLYMER_MODULE_HH
30#include <opm/material/densead/Math.hpp>
39#include <dune/common/fvector.hh>
53template<
class TypeTag,
class MyTypeTag>
55template<
class TypeTag,
class MyTypeTag>
57template<
class TypeTag,
class MyTypeTag>
59template<
class TypeTag,
class MyTypeTag>
61template<
class TypeTag,
class MyTypeTag>
63template<
class TypeTag,
class MyTypeTag>
67template<
class TypeTag>
69template<
class TypeTag>
71template<
class TypeTag>
73template<
class TypeTag>
75template<
class TypeTag>
77template<
class TypeTag>
88template <
class TypeTag>
99 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
102 enum { enablePolymer = getPropValue<TypeTag, Properties::EnablePolymer>() };
120 Parameters::registerParam<TypeTag, Properties::VtkWritePolymerConcentration>
121 (
"Include the concentration of the polymer component in the water phase "
122 "in the VTK output files");
123 Parameters::registerParam<TypeTag, Properties::VtkWritePolymerDeadPoreVolume>
124 (
"Include the fraction of the \"dead\" pore volume "
125 "in the VTK output files");
126 Parameters::registerParam<TypeTag, Properties::VtkWritePolymerRockDensity>
127 (
"Include the amount of already adsorbed polymer component"
128 "in the VTK output files");
129 Parameters::registerParam<TypeTag, Properties::VtkWritePolymerAdsorption>
130 (
"Include the adsorption rate of the polymer component"
131 "in the VTK output files");
132 Parameters::registerParam<TypeTag, Properties::VtkWritePolymerViscosityCorrection>
133 (
"Include the viscosity correction of the polymer component "
134 "in the VTK output files");
135 Parameters::registerParam<TypeTag, Properties::VtkWriteWaterViscosityCorrection>
136 (
"Include the viscosity correction of the water component "
137 "due to polymers in the VTK output files");
146 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
152 if (polymerConcentrationOutput_())
154 if (polymerDeadPoreVolumeOutput_())
156 if (polymerRockDensityOutput_())
158 if (polymerAdsorptionOutput_())
160 if (polymerViscosityCorrectionOutput_())
162 if (waterViscosityCorrectionOutput_())
172 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
178 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
179 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, 0);
180 unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
182 if (polymerConcentrationOutput_())
183 polymerConcentration_[globalDofIdx] =
184 scalarValue(intQuants.polymerConcentration());
186 if (polymerDeadPoreVolumeOutput_())
187 polymerDeadPoreVolume_[globalDofIdx] =
188 scalarValue(intQuants.polymerDeadPoreVolume());
190 if (polymerRockDensityOutput_())
191 polymerRockDensity_[globalDofIdx] =
192 scalarValue(intQuants.polymerRockDensity());
194 if (polymerAdsorptionOutput_())
195 polymerAdsorption_[globalDofIdx] =
196 scalarValue(intQuants.polymerAdsorption());
198 if (polymerViscosityCorrectionOutput_())
199 polymerViscosityCorrection_[globalDofIdx] =
200 scalarValue(intQuants.polymerViscosityCorrection());
202 if (waterViscosityCorrectionOutput_())
203 waterViscosityCorrection_[globalDofIdx] =
204 scalarValue(intQuants.waterViscosityCorrection());
220 if (polymerConcentrationOutput_())
223 if (polymerDeadPoreVolumeOutput_())
224 this->
commitScalarBuffer_(baseWriter,
"dead pore volume fraction", polymerDeadPoreVolume_);
226 if (polymerRockDensityOutput_())
229 if (polymerAdsorptionOutput_())
232 if (polymerViscosityCorrectionOutput_())
233 this->
commitScalarBuffer_(baseWriter,
"polymer viscosity correction", polymerViscosityCorrection_);
235 if (waterViscosityCorrectionOutput_())
236 this->
commitScalarBuffer_(baseWriter,
"water viscosity correction", waterViscosityCorrection_);
240 static bool polymerConcentrationOutput_()
242 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerConcentration>();
246 static bool polymerDeadPoreVolumeOutput_()
248 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerDeadPoreVolume>();
252 static bool polymerRockDensityOutput_()
254 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerRockDensity>();
258 static bool polymerAdsorptionOutput_()
260 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerAdsorption>();
264 static bool polymerViscosityCorrectionOutput_()
266 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerViscosityCorrection>();
270 static bool waterViscosityCorrectionOutput_()
272 static bool val = Parameters::get<TypeTag, Properties::VtkWritePolymerViscosityCorrection>();
276 ScalarBuffer polymerConcentration_;
277 ScalarBuffer polymerDeadPoreVolume_;
278 ScalarBuffer polymerRockDensity_;
279 ScalarBuffer polymerAdsorption_;
280 ScalarBuffer polymerViscosityCorrection_;
281 ScalarBuffer waterViscosityCorrection_;
Declares the properties required by the black oil model.
The base class for writer modules.
Definition: baseoutputmodule.hh:73
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:91
void commitScalarBuffer_(BaseOutputWriter &baseWriter, const char *name, ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer containing scalar quantities to the result file.
Definition: baseoutputmodule.hh:316
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:162
The base class for all output writers.
Definition: baseoutputwriter.hh:44
VTK output module for the black oil model's polymer related quantities.
Definition: vtkblackoilpolymermodule.hh:90
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilpolymermodule.hh:144
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilpolymermodule.hh:115
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkblackoilpolymermodule.hh:170
VtkBlackOilPolymerModule(const Simulator &simulator)
Definition: vtkblackoilpolymermodule.hh:107
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkblackoilpolymermodule.hh:211
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:66
Definition: blackoilmodel.hh:72
Definition: blackoilboundaryratevector.hh:37
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition: propertysystem.hh:242
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Definition: vtkblackoilpolymermodule.hh:48
a tag to mark properties as undefined
Definition: propertysystem.hh:40
Definition: vtkblackoilpolymermodule.hh:58
Definition: vtkblackoilpolymermodule.hh:54
Definition: vtkblackoilpolymermodule.hh:56
Definition: vtkblackoilpolymermodule.hh:60
Definition: vtkblackoilpolymermodule.hh:62
Definition: vtkblackoilpolymermodule.hh:64