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>
64template <
class TypeTag>
79 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
82 enum { dim = GridView::dimension };
83 enum { dimWorld = GridView::dimensionworld };
84 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
133 const std::size_t nDof = this->
simulator_.model().numGridDof();
134 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
135 fractureVelocity_[phaseIdx].resize(nDof);
136 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
137 fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
138 fractureVelocity_[phaseIdx][dofIdx] = 0.0;
151 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
155 const auto& fractureMapper = elemCtx.simulator().vanguard().fractureMapper();
157 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
158 const unsigned I = elemCtx.globalSpaceIndex(i, 0);
159 if (!fractureMapper.isFractureVertex(I)) {
163 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
164 const auto& fs = intQuants.fractureFluidState();
167 Valgrind::CheckDefined(intQuants.fracturePorosity());
168 fracturePorosity_[I] = intQuants.fracturePorosity();
171 const auto& K = intQuants.fractureIntrinsicPermeability();
172 fractureIntrinsicPermeability_[I] = K[0][0];
175 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
177 Valgrind::CheckDefined(fs.saturation(phaseIdx));
178 fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
181 Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
182 fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
185 Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
186 fractureRelativePermeability_[phaseIdx][I] =
187 intQuants.fractureRelativePermeability(phaseIdx);
190 Valgrind::CheckDefined(intQuants.fractureVolume());
191 fractureVolumeFraction_[I] += intQuants.fractureVolume();
198 for (
unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(0); ++scvfIdx) {
199 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, 0);
201 const unsigned i = extQuants.interiorIndex();
202 const unsigned I = elemCtx.globalSpaceIndex(i, 0);
204 const unsigned j = extQuants.exteriorIndex();
205 const unsigned J = elemCtx.globalSpaceIndex(j, 0);
207 if (!fractureMapper.isFractureEdge(I, J)) {
211 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
212 Scalar weight = std::max(Scalar{1e-16},
213 std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
214 Valgrind::CheckDefined(extQuants.extrusionFactor());
215 assert(extQuants.extrusionFactor() > 0);
216 weight *= extQuants.extrusionFactor();
218 Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
221 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
222 fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
223 fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
226 fractureVelocityWeight_[phaseIdx][I] += weight;
227 fractureVelocityWeight_[phaseIdx][J] += weight;
244 fractureSaturation_, BufferType::Dof);
248 fractureMobility_, BufferType::Dof);
252 fractureRelativePermeability_, BufferType::Dof);
257 fracturePorosity_, BufferType::Dof);
261 fractureIntrinsicPermeability_, BufferType::Dof);
265 for (
unsigned I = 0; I < fractureVolumeFraction_.size(); ++I) {
266 fractureVolumeFraction_[I] /= this->
simulator_.model().dofTotalVolume(I);
269 fractureVolumeFraction_, BufferType::Dof);
273 const std::size_t nDof = this->
simulator_.model().numGridDof();
274 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
277 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
278 fractureVelocity_[phaseIdx][dofIdx] /=
279 std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
283 snprintf(name, 512,
"fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
285 DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
292 PhaseBuffer fractureSaturation_{};
293 PhaseBuffer fractureMobility_{};
294 PhaseBuffer fractureRelativePermeability_{};
296 ScalarBuffer fracturePorosity_{};
297 ScalarBuffer fractureVolumeFraction_{};
298 ScalarBuffer fractureIntrinsicPermeability_{};
300 PhaseVectorBuffer fractureVelocity_{};
301 PhaseBuffer fractureVelocityWeight_{};
303 PhaseVectorBuffer potentialGradient_{};
304 PhaseBuffer potentialWeight_{};
The base class for writer modules.
Definition: baseoutputmodule.hh:68
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:86
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:157
const Simulator & simulator_
Definition: baseoutputmodule.hh:426
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:95
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:91
BufferType
Definition: baseoutputmodule.hh:143
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:337
void commitScalarBuffer_(BaseOutputWriter &baseWriter, const char *name, ScalarBuffer &buffer, BufferType bufferType)
Add a buffer containing scalar quantities to the result file.
Definition: baseoutputmodule.hh:238
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:198
The base class for all output writers.
Definition: baseoutputwriter.hh:46
VTK output module for quantities which make sense for all models which deal with discrete fractures i...
Definition: vtkdiscretefracturemodule.hpp:66
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkdiscretefracturemodule.hpp:149
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkdiscretefracturemodule.hpp:101
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkdiscretefracturemodule.hpp:236
VtkDiscreteFractureModule(const Simulator &simulator)
Definition: vtkdiscretefracturemodule.hpp:92
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiscretefracturemodule.hpp:110
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilboundaryratevector.hh:39
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:233
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.
bool velocityOutput_
Definition: vtkdiscretefractureparams.hpp:62
bool relativePermeabilityOutput_
Definition: vtkdiscretefractureparams.hpp:58
bool porosityOutput_
Definition: vtkdiscretefractureparams.hpp:59
static void registerParameters()
Registers the parameters in parameter system.
bool intrinsicPermeabilityOutput_
Definition: vtkdiscretefractureparams.hpp:60
bool volumeFractionOutput_
Definition: vtkdiscretefractureparams.hpp:61
bool saturationOutput_
Definition: vtkdiscretefractureparams.hpp:56
bool mobilityOutput_
Definition: vtkdiscretefractureparams.hpp:57