25 #ifndef EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
26 #define EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
34 #include <dune/common/fvector.hh>
39 namespace Properties {
56 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureSaturations,
true);
57 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureMobilities,
false);
58 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureRelativePermeabilities,
true);
59 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFracturePorosity,
true);
60 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureIntrinsicPermeabilities,
false);
61 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureFilterVelocities,
false);
62 SET_BOOL_PROP(VtkDiscreteFracture, VtkWriteFractureVolumeFraction,
true);
78 template <
class TypeTag>
84 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
85 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
87 typedef typename GET_PROP_TYPE(TypeTag, GridManager) GridManager;
88 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
89 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
91 typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
93 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
96 enum { dim = GridView::dimension };
97 enum { dimWorld = GridView::dimensionworld };
106 : ParentType(simulator)
115 "Include the phase saturations in the VTK output files");
117 "Include the phase mobilities in the VTK output files");
119 "Include the phase relative permeabilities in the "
122 "Include the porosity in the VTK output files");
124 "Include the intrinsic permeability in the VTK output files");
126 "Include in the filter velocities of the phases "
127 "the VTK output files");
129 "Add the fraction of the total volume which is "
130 "occupied by fractures in the VTK output");
139 if (saturationOutput_())
141 if (mobilityOutput_())
143 if (relativePermeabilityOutput_())
146 if (porosityOutput_())
148 if (intrinsicPermeabilityOutput_())
150 if (volumeFractionOutput_())
153 if (velocityOutput_()) {
155 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
156 fractureVelocity_[phaseIdx].resize(nDof);
157 for (
int dofIdx = 0; dofIdx < nDof; ++dofIdx) {
158 fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
159 fractureVelocity_[phaseIdx][dofIdx] = 0.0;
172 const auto &fractureMapper = elemCtx.simulator().gridManager().fractureMapper();
174 for (
int i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
175 int I = elemCtx.globalSpaceIndex(i, 0);
176 if (!fractureMapper.isFractureVertex(I))
179 const auto &intQuants = elemCtx.intensiveQuantities(i, 0);
180 const auto &fs = intQuants.fractureFluidState();
182 if (porosityOutput_())
183 fracturePorosity_[I] = intQuants.fracturePorosity();
184 if (intrinsicPermeabilityOutput_()) {
185 const auto &K = intQuants.fractureIntrinsicPermeability();
186 fractureIntrinsicPermeability_[I] = K[0][0];
189 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
190 if (saturationOutput_())
191 fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
192 if (mobilityOutput_())
193 fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
194 if (relativePermeabilityOutput_())
195 fractureRelativePermeability_[phaseIdx][I] =
196 intQuants.fractureRelativePermeability(phaseIdx);
197 if (volumeFractionOutput_())
198 fractureVolumeFraction_[I] += intQuants.fractureVolume();
202 if (velocityOutput_()) {
204 for (
int scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(0); ++ scvfIdx) {
205 const auto &extQuants = elemCtx.extensiveQuantities(scvfIdx, 0);
207 int i = extQuants.interiorIndex();
208 int I = elemCtx.globalSpaceIndex(i, 0);
210 int j = extQuants.exteriorIndex();
211 int J = elemCtx.globalSpaceIndex(j, 0);
213 if (!fractureMapper.isFractureEdge(I, J))
216 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
217 Scalar weight = std::max(1e-16, std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
218 Valgrind::CheckDefined(extQuants.extrusionFactor());
219 assert(extQuants.extrusionFactor() > 0);
220 weight *= extQuants.extrusionFactor();
222 Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
225 for (
int dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
226 fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
227 fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
230 fractureVelocityWeight_[phaseIdx][I] += weight;
231 fractureVelocityWeight_[phaseIdx][J] += weight;
242 VtkMultiWriter *vtkWriter =
dynamic_cast<VtkMultiWriter*
>(&baseWriter);
247 if (saturationOutput_())
249 if (mobilityOutput_())
251 if (relativePermeabilityOutput_())
252 this->
commitPhaseBuffer_(baseWriter,
"fractureRelativePerm_%s", fractureRelativePermeability_);
254 if (porosityOutput_())
256 if (intrinsicPermeabilityOutput_())
257 this->
commitScalarBuffer_(baseWriter,
"fractureIntrinsicPerm", fractureIntrinsicPermeability_);
258 if (volumeFractionOutput_()) {
260 for (
unsigned I = 0; I < fractureVolumeFraction_.size(); ++I)
261 fractureVolumeFraction_[I] /= this->
simulator_.
model().dofTotalVolume(I);
265 if (velocityOutput_()) {
268 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
271 for (
int dofIdx = 0; dofIdx < nDof; ++dofIdx)
272 fractureVelocity_[phaseIdx][dofIdx] /=
273 std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
276 snprintf(name, 512,
"fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx));
278 DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
284 static bool saturationOutput_()
285 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureSaturations); }
287 static bool mobilityOutput_()
290 static bool relativePermeabilityOutput_()
291 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureRelativePermeabilities); }
293 static bool porosityOutput_()
296 static bool intrinsicPermeabilityOutput_()
297 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureIntrinsicPermeabilities); }
299 static bool volumeFractionOutput_()
300 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureVolumeFraction); }
302 static bool velocityOutput_()
303 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteFractureFilterVelocities); }
305 PhaseBuffer fractureSaturation_;
306 PhaseBuffer fractureMobility_;
307 PhaseBuffer fractureRelativePermeability_;
309 ScalarBuffer fracturePorosity_;
310 ScalarBuffer fractureVolumeFraction_;
311 ScalarBuffer fractureIntrinsicPermeability_;
313 PhaseVectorBuffer fractureVelocity_;
314 PhaseBuffer fractureVelocityWeight_;
316 PhaseVectorBuffer potentialGradient_;
317 PhaseBuffer potentialWeight_;
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:148
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:213
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiscretefracturemodule.hh:137
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkdiscretefracturemodule.hh:112
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkdiscretefracturemodule.hh:240
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:99
NEW_PROP_TAG(Grid)
The type of the DUNE grid.
This file provides the infrastructure to retrieve run-time parameters.
The base class for all output writers.
Definition: baseoutputwriter.hh:41
Simplifies writing multi-file VTK datasets.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:73
Definition: baseauxiliarymodule.hh:35
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:64
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:176
VtkDiscreteFractureModule(const Simulator &simulator)
Definition: vtkdiscretefracturemodule.hh:105
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:95
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:283
VTK output module for quantities which make sense for all models which deal with discrete fractures i...
Definition: vtkdiscretefracturemodule.hh:79
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element...
Definition: vtkdiscretefracturemodule.hh:170
Provides the magic behind the eWoms property system.
The base class for writer modules.
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:386
SET_BOOL_PROP(FvBaseDiscretization, EnableVtkOutput, true)
Enable the VTK output by default.
The base class for writer modules.
Definition: baseoutputmodule.hh:71
const Simulator & simulator_
Definition: baseoutputmodule.hh:487
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:61
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:95
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:90