27 #ifndef EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
28 #define EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
36 #include <opm/common/Valgrind.hpp>
38 #include <dune/common/fvector.hh>
43 namespace Properties {
61 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureSaturations,
true);
62 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureMobilities,
false);
63 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureRelativePermeabilities,
true);
64 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFracturePorosity,
true);
65 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureIntrinsicPermeabilities,
false);
66 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureFilterVelocities,
false);
67 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureVolumeFraction,
true);
83 template <
class TypeTag>
89 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
90 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
92 typedef typename GET_PROP_TYPE(TypeTag, GridManager) GridManager;
93 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
94 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
96 typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
98 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
101 enum { dim = GridView::dimension };
102 enum { dimWorld = GridView::dimensionworld };
105 typedef typename ParentType::ScalarBuffer ScalarBuffer;
106 typedef typename ParentType::PhaseBuffer PhaseBuffer;
107 typedef typename ParentType::PhaseVectorBuffer PhaseVectorBuffer;
120 "Include the phase saturations in the VTK output files");
122 "Include the phase mobilities in the VTK output files");
124 "Include the phase relative permeabilities in the "
127 "Include the porosity in the VTK output files");
129 "Include the intrinsic permeability in the VTK output files");
131 "Include in the filter velocities of the phases "
132 "the VTK output files");
134 "Add the fraction of the total volume which is "
135 "occupied by fractures in the VTK output");
144 if (saturationOutput_())
146 if (mobilityOutput_())
148 if (relativePermeabilityOutput_())
151 if (porosityOutput_())
153 if (intrinsicPermeabilityOutput_())
155 if (volumeFractionOutput_())
158 if (velocityOutput_()) {
159 size_t nDof = this->simulator_.
model().numGridDof();
160 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
161 fractureVelocity_[phaseIdx].resize(nDof);
162 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
163 fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
164 fractureVelocity_[phaseIdx][dofIdx] = 0.0;
180 const auto& fractureMapper = elemCtx.simulator().gridManager().fractureMapper();
182 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
183 unsigned I = elemCtx.globalSpaceIndex(i, 0);
184 if (!fractureMapper.isFractureVertex(I))
187 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
188 const auto& fs = intQuants.fractureFluidState();
190 if (porosityOutput_()) {
191 Opm::Valgrind::CheckDefined(intQuants.fracturePorosity());
192 fracturePorosity_[I] = intQuants.fracturePorosity();
194 if (intrinsicPermeabilityOutput_()) {
195 const auto& K = intQuants.fractureIntrinsicPermeability();
196 fractureIntrinsicPermeability_[I] = K[0][0];
199 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
200 if (saturationOutput_()) {
201 Opm::Valgrind::CheckDefined(fs.saturation(phaseIdx));
202 fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
204 if (mobilityOutput_()) {
205 Opm::Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
206 fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
208 if (relativePermeabilityOutput_()) {
209 Opm::Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
210 fractureRelativePermeability_[phaseIdx][I] =
211 intQuants.fractureRelativePermeability(phaseIdx);
213 if (volumeFractionOutput_()) {
214 Opm::Valgrind::CheckDefined(intQuants.fractureVolume());
215 fractureVolumeFraction_[I] += intQuants.fractureVolume();
220 if (velocityOutput_()) {
222 for (
unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(0); ++ scvfIdx) {
223 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, 0);
225 unsigned i = extQuants.interiorIndex();
226 unsigned I = elemCtx.globalSpaceIndex(i, 0);
228 unsigned j = extQuants.exteriorIndex();
229 unsigned J = elemCtx.globalSpaceIndex(j, 0);
231 if (!fractureMapper.isFractureEdge(I, J))
234 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
236 std::max<Scalar>(1e-16, std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
237 Opm::Valgrind::CheckDefined(extQuants.extrusionFactor());
238 assert(extQuants.extrusionFactor() > 0);
239 weight *= extQuants.extrusionFactor();
241 Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
244 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
245 fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
246 fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
249 fractureVelocityWeight_[phaseIdx][I] += weight;
250 fractureVelocityWeight_[phaseIdx][J] += weight;
266 if (saturationOutput_())
268 if (mobilityOutput_())
270 if (relativePermeabilityOutput_())
271 this->
commitPhaseBuffer_(baseWriter,
"fractureRelativePerm_%s", fractureRelativePermeability_);
273 if (porosityOutput_())
275 if (intrinsicPermeabilityOutput_())
276 this->
commitScalarBuffer_(baseWriter,
"fractureIntrinsicPerm", fractureIntrinsicPermeability_);
277 if (volumeFractionOutput_()) {
279 for (
unsigned I = 0; I < fractureVolumeFraction_.size(); ++I)
280 fractureVolumeFraction_[I] /= this->simulator_.
model().dofTotalVolume(I);
284 if (velocityOutput_()) {
285 size_t nDof = this->simulator_.
model().numGridDof();
287 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
290 for (
unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx)
291 fractureVelocity_[phaseIdx][dofIdx] /=
292 std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
295 snprintf(name, 512,
"fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx));
297 DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
303 static bool saturationOutput_()
304 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureSaturations); }
306 static bool mobilityOutput_()
309 static bool relativePermeabilityOutput_()
310 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureRelativePermeabilities); }
312 static bool porosityOutput_()
315 static bool intrinsicPermeabilityOutput_()
316 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureIntrinsicPermeabilities); }
318 static bool volumeFractionOutput_()
319 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureVolumeFraction); }
321 static bool velocityOutput_()
322 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureFilterVelocities); }
324 PhaseBuffer fractureSaturation_;
325 PhaseBuffer fractureMobility_;
326 PhaseBuffer fractureRelativePermeability_;
328 ScalarBuffer fracturePorosity_;
329 ScalarBuffer fractureVolumeFraction_;
330 ScalarBuffer fractureIntrinsicPermeability_;
332 PhaseVectorBuffer fractureVelocity_;
333 PhaseBuffer fractureVelocityWeight_;
335 PhaseVectorBuffer potentialGradient_;
336 PhaseBuffer potentialWeight_;
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkdiscretefracturemodule.hh:117
The base class for all output writers.
Definition: baseoutputwriter.hh:43
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:408
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
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:305
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
The base class for writer modules.
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
This file provides the infrastructure to retrieve run-time parameters.
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiscretefracturemodule.hh:142
The base class for writer modules.
Definition: baseoutputmodule.hh:80
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element...
Definition: vtkdiscretefracturemodule.hh:175
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkdiscretefracturemodule.hh:259
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:235
Provides the magic behind the eWoms property system.
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:170
Simplifies writing multi-file VTK datasets.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
VTK output module for quantities which make sense for all models which deal with discrete fractures i...
Definition: vtkdiscretefracturemodule.hh:84