27#ifndef EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
28#define EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
36#include <opm/material/common/Valgrind.hpp>
38#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>
65template<
class TypeTag,
class MyTypeTag>
69template<
class TypeTag>
71template<
class TypeTag>
73template<
class TypeTag>
75template<
class TypeTag>
77template<
class TypeTag>
79template<
class TypeTag>
81template<
class TypeTag>
100template <
class TypeTag>
115 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
118 enum { dim = GridView::dimension };
119 enum { dimWorld = GridView::dimensionworld };
120 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
136 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureSaturations>
137 (
"Include the phase saturations in the VTK output files");
138 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureMobilities>
139 (
"Include the phase mobilities in the VTK output files");
140 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureRelativePermeabilities>
141 (
"Include the phase relative permeabilities in the "
143 Parameters::registerParam<TypeTag, Properties::VtkWriteFracturePorosity>
144 (
"Include the porosity in the VTK output files");
145 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureIntrinsicPermeabilities>
146 (
"Include the intrinsic permeability in the VTK output files");
147 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureFilterVelocities>
148 (
"Include in the filter velocities of the phases in the VTK output files");
149 Parameters::registerParam<TypeTag, Properties::VtkWriteFractureVolumeFraction>
150 (
"Add the fraction of the total volume which is "
151 "occupied by fractures in the VTK output");
160 if (saturationOutput_())
162 if (mobilityOutput_())
164 if (relativePermeabilityOutput_())
167 if (porosityOutput_())
169 if (intrinsicPermeabilityOutput_())
171 if (volumeFractionOutput_())
174 if (velocityOutput_()) {
175 size_t nDof = this->
simulator_.model().numGridDof();
176 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
177 fractureVelocity_[phaseIdx].resize(nDof);
178 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
179 fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
180 fractureVelocity_[phaseIdx][dofIdx] = 0.0;
193 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
196 const auto& fractureMapper = elemCtx.simulator().vanguard().fractureMapper();
198 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
199 unsigned I = elemCtx.globalSpaceIndex(i, 0);
200 if (!fractureMapper.isFractureVertex(I))
203 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
204 const auto& fs = intQuants.fractureFluidState();
206 if (porosityOutput_()) {
207 Opm::Valgrind::CheckDefined(intQuants.fracturePorosity());
208 fracturePorosity_[I] = intQuants.fracturePorosity();
210 if (intrinsicPermeabilityOutput_()) {
211 const auto& K = intQuants.fractureIntrinsicPermeability();
212 fractureIntrinsicPermeability_[I] = K[0][0];
215 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
216 if (saturationOutput_()) {
217 Opm::Valgrind::CheckDefined(fs.saturation(phaseIdx));
218 fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
220 if (mobilityOutput_()) {
221 Opm::Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
222 fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
224 if (relativePermeabilityOutput_()) {
225 Opm::Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
226 fractureRelativePermeability_[phaseIdx][I] =
227 intQuants.fractureRelativePermeability(phaseIdx);
229 if (volumeFractionOutput_()) {
230 Opm::Valgrind::CheckDefined(intQuants.fractureVolume());
231 fractureVolumeFraction_[I] += intQuants.fractureVolume();
236 if (velocityOutput_()) {
238 for (
unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(0); ++ scvfIdx) {
239 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, 0);
241 unsigned i = extQuants.interiorIndex();
242 unsigned I = elemCtx.globalSpaceIndex(i, 0);
244 unsigned j = extQuants.exteriorIndex();
245 unsigned J = elemCtx.globalSpaceIndex(j, 0);
247 if (!fractureMapper.isFractureEdge(I, J))
250 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
252 std::max<Scalar>(1e-16, std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
253 Opm::Valgrind::CheckDefined(extQuants.extrusionFactor());
254 assert(extQuants.extrusionFactor() > 0);
255 weight *= extQuants.extrusionFactor();
257 Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
260 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
261 fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
262 fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
265 fractureVelocityWeight_[phaseIdx][I] += weight;
266 fractureVelocityWeight_[phaseIdx][J] += weight;
282 if (saturationOutput_())
284 if (mobilityOutput_())
286 if (relativePermeabilityOutput_())
287 this->
commitPhaseBuffer_(baseWriter,
"fractureRelativePerm_%s", fractureRelativePermeability_);
289 if (porosityOutput_())
291 if (intrinsicPermeabilityOutput_())
292 this->
commitScalarBuffer_(baseWriter,
"fractureIntrinsicPerm", fractureIntrinsicPermeability_);
293 if (volumeFractionOutput_()) {
295 for (
unsigned I = 0; I < fractureVolumeFraction_.size(); ++I)
296 fractureVolumeFraction_[I] /= this->
simulator_.model().dofTotalVolume(I);
300 if (velocityOutput_()) {
301 size_t nDof = this->
simulator_.model().numGridDof();
303 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
306 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx)
307 fractureVelocity_[phaseIdx][dofIdx] /=
308 std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
311 snprintf(name, 512,
"fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
313 DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
319 static bool saturationOutput_()
321 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureSaturations>();
325 static bool mobilityOutput_()
327 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureMobilities>();
331 static bool relativePermeabilityOutput_()
333 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureRelativePermeabilities>();
337 static bool porosityOutput_()
339 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFracturePorosity>();
343 static bool intrinsicPermeabilityOutput_()
345 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureIntrinsicPermeabilities>();
349 static bool volumeFractionOutput_()
351 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureVolumeFraction>();
355 static bool velocityOutput_()
357 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFractureFilterVelocities>();
361 PhaseBuffer fractureSaturation_;
362 PhaseBuffer fractureMobility_;
363 PhaseBuffer fractureRelativePermeability_;
365 ScalarBuffer fracturePorosity_;
366 ScalarBuffer fractureVolumeFraction_;
367 ScalarBuffer fractureIntrinsicPermeability_;
369 PhaseVectorBuffer fractureVelocity_;
370 PhaseBuffer fractureVelocityWeight_;
372 PhaseVectorBuffer potentialGradient_;
373 PhaseBuffer potentialWeight_;
The base class for writer modules.
Definition: baseoutputmodule.hh:73
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:91
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:419
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:246
const Simulator & simulator_
Definition: baseoutputmodule.hh:519
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:100
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:96
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 quantities which make sense for all models which deal with discrete fractures i...
Definition: vtkdiscretefracturemodule.hh:102
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiscretefracturemodule.hh:158
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkdiscretefracturemodule.hh:134
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkdiscretefracturemodule.hh:275
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkdiscretefracturemodule.hh:191
VtkDiscreteFractureModule(const Simulator &simulator)
Definition: vtkdiscretefracturemodule.hh:127
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: vtkdiscretefracturemodule.hh:48
a tag to mark properties as undefined
Definition: propertysystem.hh:40
Definition: vtkdiscretefracturemodule.hh:64
Definition: vtkdiscretefracturemodule.hh:62
Definition: vtkdiscretefracturemodule.hh:56
Definition: vtkdiscretefracturemodule.hh:60
Definition: vtkdiscretefracturemodule.hh:58
Definition: vtkdiscretefracturemodule.hh:54
Definition: vtkdiscretefracturemodule.hh:66