27#ifndef OPM_VTK_MULTI_PHASE_MODULE_HPP
28#define OPM_VTK_MULTI_PHASE_MODULE_HPP
30#include <dune/common/fvector.hh>
32#include <opm/material/common/MathToolbox.hpp>
33#include <opm/material/common/Valgrind.hpp>
66template<
class TypeTag>
82 enum { dimWorld = GridView::dimensionworld };
85 using ScalarBuffer =
typename ParentType::ScalarBuffer;
86 using VectorBuffer =
typename ParentType::VectorBuffer;
87 using TensorBuffer =
typename ParentType::TensorBuffer;
88 using PhaseBuffer =
typename ParentType::PhaseBuffer;
90 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
92 using PhaseVectorBuffer = std::array<VectorBuffer, numPhases>;
115 if (params_.extrusionFactorOutput_) {
118 if (params_.pressureOutput_) {
121 if (params_.densityOutput_) {
124 if (params_.saturationOutput_) {
127 if (params_.mobilityOutput_) {
130 if (params_.relativePermeabilityOutput_) {
133 if (params_.viscosityOutput_) {
136 if (params_.averageMolarMassOutput_) {
140 if (params_.porosityOutput_) {
143 if (params_.intrinsicPermeabilityOutput_) {
147 if (params_.velocityOutput_) {
148 size_t nDof = this->simulator_.model().numGridDof();
159 if (params_.potentialGradientOutput_) {
160 size_t nDof = this->simulator_.model().numGridDof();
179 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
183 const auto& problem =
elemCtx.problem();
184 for (
unsigned i = 0; i <
elemCtx.numPrimaryDof(0); ++i) {
185 unsigned I =
elemCtx.globalSpaceIndex(i, 0);
189 if (params_.extrusionFactorOutput_) {
190 extrusionFactor_[
I] =
intQuants.extrusionFactor();
192 if (params_.porosityOutput_) {
196 if (params_.intrinsicPermeabilityOutput_) {
197 const auto&
K = problem.intrinsicPermeability(
elemCtx, i, 0);
206 if (!FluidSystem::phaseIsActive(
phaseIdx)) {
209 if (params_.pressureOutput_) {
212 if (params_.densityOutput_) {
215 if (params_.saturationOutput_) {
218 if (params_.mobilityOutput_) {
221 if (params_.relativePermeabilityOutput_) {
224 if (params_.viscosityOutput_) {
227 if (params_.averageMolarMassOutput_) {
233 if (params_.potentialGradientOutput_) {
239 unsigned I =
elemCtx.globalSpaceIndex(i, 0);
242 Scalar weight =
extQuants.extrusionFactor();
256 if (params_.velocityOutput_) {
262 unsigned I =
elemCtx.globalSpaceIndex(i, 0);
265 unsigned J =
elemCtx.globalSpaceIndex(j, 0);
268 Scalar weight = std::max<Scalar>(1
e-16,
270 Valgrind::CheckDefined(
extQuants.extrusionFactor());
276 for (
unsigned k = 0;
k < dimWorld; ++
k) {
279 if (
v.two_norm() > 1
e-20) {
280 weight /=
v.two_norm();
304 if (params_.extrusionFactorOutput_) {
307 if (params_.pressureOutput_) {
310 if (params_.densityOutput_) {
313 if (params_.saturationOutput_) {
316 if (params_.mobilityOutput_) {
319 if (params_.relativePermeabilityOutput_) {
322 if (params_.viscosityOutput_) {
325 if (params_.averageMolarMassOutput_) {
329 if (params_.porosityOutput_) {
332 if (params_.intrinsicPermeabilityOutput_) {
336 if (params_.velocityOutput_) {
337 size_t numDof = this->simulator_.model().numGridDof();
342 for (
unsigned i = 0; i < numDof; ++i) {
347 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(
phaseIdx).data());
353 if (params_.potentialGradientOutput_) {
354 size_t numDof = this->simulator_.model().numGridDof();
359 for (
unsigned i = 0; i < numDof; ++i) {
366 DiscBaseOutputModule::attachVectorDofData_(
baseWriter,
383 return params_.velocityOutput_ || params_.potentialGradientOutput_;
388 ScalarBuffer extrusionFactor_{};
389 PhaseBuffer pressure_{};
390 PhaseBuffer density_{};
391 PhaseBuffer saturation_{};
392 PhaseBuffer mobility_{};
393 PhaseBuffer relativePermeability_{};
394 PhaseBuffer viscosity_{};
395 PhaseBuffer averageMolarMass_{};
397 ScalarBuffer porosity_{};
398 TensorBuffer intrinsicPermeability_{};
400 PhaseVectorBuffer velocity_{};
401 PhaseBuffer velocityWeight_{};
403 PhaseVectorBuffer potentialGradient_{};
404 PhaseBuffer potentialWeight_{};
The base class for writer modules.
The base class for writer modules.
Definition baseoutputmodule.hh:67
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition baseoutputmodule.hh:345
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition baseoutputmodule.hh:201
void commitTensorBuffer_(BaseOutputWriter &baseWriter, const char *name, TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer containing tensorial quantities to the result file.
Definition baseoutputmodule.hh:288
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition baseoutputmodule.hh:166
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:244
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition baseoutputmodule.hh:156
The base class for all output writers.
Definition baseoutputwriter.hh:44
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition vtkmultiphasemodule.hpp:68
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition vtkmultiphasemodule.hpp:177
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition vtkmultiphasemodule.hpp:104
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition vtkmultiphasemodule.hpp:113
bool needExtensiveQuantities() const final
Returns true iff the module needs to access the extensive quantities of a context to do its job.
Definition vtkmultiphasemodule.hpp:381
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition vtkmultiphasemodule.hpp:297
Simplifies writing multi-file VTK datasets.
Definition vtkmultiwriter.hh:66
Declare the properties used by the infrastructure code of the finite volume discretizations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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:235
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Struct holding the parameters for VtkMultiPhaseModule.
Definition vtkmultiphaseparams.hpp:54
static void registerParameters()
Registers the parameters in parameter system.
Definition vtkmultiphaseparams.cpp:31
void read()
Reads the parameter values from the parameter system.
Definition vtkmultiphaseparams.cpp:59
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Simplifies writing multi-file VTK datasets.