27#ifndef OPM_VTK_DISCRETE_FRACTURE_MODULE_HPP
28#define OPM_VTK_DISCRETE_FRACTURE_MODULE_HPP
30#include <dune/common/fvector.hh>
32#include <opm/material/common/Valgrind.hpp>
60template <
class TypeTag>
78 enum { dim = GridView::dimension };
79 enum { dimWorld = GridView::dimensionworld };
82 using ScalarBuffer =
typename ParentType::ScalarBuffer;
83 using PhaseBuffer =
typename ParentType::PhaseBuffer;
84 using PhaseVectorBuffer =
typename ParentType::PhaseVectorBuffer;
107 if (params_.saturationOutput_) {
110 if (params_.mobilityOutput_) {
113 if (params_.relativePermeabilityOutput_) {
117 if (params_.porosityOutput_) {
120 if (params_.intrinsicPermeabilityOutput_) {
123 if (params_.volumeFractionOutput_) {
127 if (params_.velocityOutput_) {
128 size_t nDof = this->simulator_.model().numGridDof();
131 for (
unsigned dofIdx = 0; dofIdx <
nDof; ++dofIdx) {
132 fractureVelocity_[
phaseIdx][dofIdx].resize(dimWorld);
133 fractureVelocity_[
phaseIdx][dofIdx] = 0.0;
146 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
150 const auto& fractureMapper = elemCtx.simulator().vanguard().fractureMapper();
152 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
153 unsigned I = elemCtx.globalSpaceIndex(i, 0);
154 if (!fractureMapper.isFractureVertex(I)) {
158 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
159 const auto& fs = intQuants.fractureFluidState();
161 if (params_.porosityOutput_) {
162 Opm::Valgrind::CheckDefined(intQuants.fracturePorosity());
163 fracturePorosity_[I] = intQuants.fracturePorosity();
165 if (params_.intrinsicPermeabilityOutput_) {
166 const auto& K = intQuants.fractureIntrinsicPermeability();
167 fractureIntrinsicPermeability_[I] = K[0][0];
171 if (params_.saturationOutput_) {
172 Opm::Valgrind::CheckDefined(fs.saturation(
phaseIdx));
175 if (params_.mobilityOutput_) {
176 Opm::Valgrind::CheckDefined(intQuants.fractureMobility(
phaseIdx));
179 if (params_.relativePermeabilityOutput_) {
180 Opm::Valgrind::CheckDefined(intQuants.fractureRelativePermeability(
phaseIdx));
181 fractureRelativePermeability_[
phaseIdx][I] =
182 intQuants.fractureRelativePermeability(
phaseIdx);
184 if (params_.volumeFractionOutput_) {
185 Opm::Valgrind::CheckDefined(intQuants.fractureVolume());
186 fractureVolumeFraction_[I] += intQuants.fractureVolume();
191 if (params_.velocityOutput_) {
197 unsigned I = elemCtx.globalSpaceIndex(i, 0);
200 unsigned J = elemCtx.globalSpaceIndex(j, 0);
202 if (!fractureMapper.isFractureEdge(I, J)) {
209 Opm::Valgrind::CheckDefined(
extQuants.extrusionFactor());
221 fractureVelocityWeight_[
phaseIdx][I] += weight;
222 fractureVelocityWeight_[
phaseIdx][J] += weight;
238 if (params_.saturationOutput_) {
241 if (params_.mobilityOutput_) {
244 if (params_.relativePermeabilityOutput_) {
245 this->
commitPhaseBuffer_(baseWriter,
"fractureRelativePerm_%s", fractureRelativePermeability_);
248 if (params_.porosityOutput_) {
251 if (params_.intrinsicPermeabilityOutput_) {
252 this->
commitScalarBuffer_(baseWriter,
"fractureIntrinsicPerm", fractureIntrinsicPermeability_);
254 if (params_.volumeFractionOutput_) {
256 for (
unsigned I = 0; I < fractureVolumeFraction_.size(); ++I) {
257 fractureVolumeFraction_[I] /= this->simulator_.model().dofTotalVolume(I);
262 if (params_.velocityOutput_) {
263 size_t nDof = this->simulator_.model().numGridDof();
268 for (
unsigned dofIdx = 0; dofIdx <
nDof; ++dofIdx) {
269 fractureVelocity_[
phaseIdx][dofIdx] /=
270 std::max<Scalar>(1
e-20, fractureVelocityWeight_[
phaseIdx][dofIdx]);
274 snprintf(name, 512,
"fractureFilterVelocity_%s", FluidSystem::phaseName(
phaseIdx).data());
276 DiscBaseOutputModule::attachVectorDofData_(
baseWriter, fractureVelocity_[
phaseIdx], name);
283 PhaseBuffer fractureSaturation_{};
284 PhaseBuffer fractureMobility_{};
285 PhaseBuffer fractureRelativePermeability_{};
287 ScalarBuffer fracturePorosity_{};
288 ScalarBuffer fractureVolumeFraction_{};
289 ScalarBuffer fractureIntrinsicPermeability_{};
291 PhaseVectorBuffer fractureVelocity_{};
292 PhaseBuffer fractureVelocityWeight_{};
294 PhaseVectorBuffer potentialGradient_{};
295 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 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 discrete fractures i...
Definition vtkdiscretefracturemodule.hpp:62
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition vtkdiscretefracturemodule.hpp:144
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition vtkdiscretefracturemodule.hpp:96
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition vtkdiscretefracturemodule.hpp:231
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition vtkdiscretefracturemodule.hpp:105
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 VtkDiscreteFractureModule.
Definition vtkdiscretefractureparams.hpp:49
void read()
Reads the parameter values from the parameter system.
Definition vtkdiscretefractureparams.cpp:51
static void registerParameters()
Registers the parameters in parameter system.
Definition vtkdiscretefractureparams.cpp:31
VTK output module for quantities which make sense for all models which deal with discrete fractures i...
Simplifies writing multi-file VTK datasets.