All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
vtkdiscretefracturemodule.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
28 #define EWOMS_VTK_DISCRETE_FRACTURE_MODULE_HH
29 
30 #include "vtkmultiwriter.hh"
31 #include "baseoutputmodule.hh"
32 
35 
36 #include <opm/common/Valgrind.hpp>
37 
38 #include <dune/common/fvector.hh>
39 
40 #include <cstdio>
41 
42 namespace Ewoms {
43 namespace Properties {
44 // create new type tag for the VTK multi-phase output
45 NEW_TYPE_TAG(VtkDiscreteFracture);
46 
47 // create the property tags needed for the multi phase module
48 NEW_PROP_TAG(GridManager);
49 NEW_PROP_TAG(VtkWriteFractureSaturations);
50 NEW_PROP_TAG(VtkWriteFractureMobilities);
51 NEW_PROP_TAG(VtkWriteFractureRelativePermeabilities);
52 NEW_PROP_TAG(VtkWriteFracturePorosity);
53 NEW_PROP_TAG(VtkWriteFractureIntrinsicPermeabilities);
54 NEW_PROP_TAG(VtkWriteFractureFilterVelocities);
55 NEW_PROP_TAG(VtkWriteFractureVolumeFraction);
56 NEW_PROP_TAG(VtkOutputFormat);
57 NEW_PROP_TAG(EnableVtkOutput);
58 NEW_PROP_TAG(DiscBaseOutputModule);
59 
60 // set default values for what quantities to output
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);
68 } // namespace Properties
69 
83 template <class TypeTag>
85 {
87 
88  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
89  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
90  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
91 
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;
95 
96  typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
97 
98  static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
100 
101  enum { dim = GridView::dimension };
102  enum { dimWorld = GridView::dimensionworld };
103  enum { numPhases = GET_PROP_VALUE(TypeTag, NumPhases) };
104 
105  typedef typename ParentType::ScalarBuffer ScalarBuffer;
106  typedef typename ParentType::PhaseBuffer PhaseBuffer;
107  typedef typename ParentType::PhaseVectorBuffer PhaseVectorBuffer;
108 
109 public:
110  VtkDiscreteFractureModule(const Simulator& simulator)
111  : ParentType(simulator)
112  { }
113 
117  static void registerParameters()
118  {
119  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureSaturations,
120  "Include the phase saturations in the VTK output files");
121  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureMobilities,
122  "Include the phase mobilities in the VTK output files");
123  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureRelativePermeabilities,
124  "Include the phase relative permeabilities in the "
125  "VTK output files");
126  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFracturePorosity,
127  "Include the porosity in the VTK output files");
128  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureIntrinsicPermeabilities,
129  "Include the intrinsic permeability in the VTK output files");
130  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureFilterVelocities,
131  "Include in the filter velocities of the phases "
132  "the VTK output files");
133  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteFractureVolumeFraction,
134  "Add the fraction of the total volume which is "
135  "occupied by fractures in the VTK output");
136  }
137 
143  {
144  if (saturationOutput_())
145  this->resizePhaseBuffer_(fractureSaturation_);
146  if (mobilityOutput_())
147  this->resizePhaseBuffer_(fractureMobility_);
148  if (relativePermeabilityOutput_())
149  this->resizePhaseBuffer_(fractureRelativePermeability_);
150 
151  if (porosityOutput_())
152  this->resizeScalarBuffer_(fracturePorosity_);
153  if (intrinsicPermeabilityOutput_())
154  this->resizeScalarBuffer_(fractureIntrinsicPermeability_);
155  if (volumeFractionOutput_())
156  this->resizeScalarBuffer_(fractureVolumeFraction_);
157 
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;
165  }
166  }
167  this->resizePhaseBuffer_(fractureVelocityWeight_);
168  }
169  }
170 
175  void processElement(const ElementContext& elemCtx)
176  {
177  if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
178  return;
179 
180  const auto& fractureMapper = elemCtx.simulator().gridManager().fractureMapper();
181 
182  for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
183  unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
184  if (!fractureMapper.isFractureVertex(I))
185  continue;
186 
187  const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
188  const auto& fs = intQuants.fractureFluidState();
189 
190  if (porosityOutput_()) {
191  Opm::Valgrind::CheckDefined(intQuants.fracturePorosity());
192  fracturePorosity_[I] = intQuants.fracturePorosity();
193  }
194  if (intrinsicPermeabilityOutput_()) {
195  const auto& K = intQuants.fractureIntrinsicPermeability();
196  fractureIntrinsicPermeability_[I] = K[0][0];
197  }
198 
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);
203  }
204  if (mobilityOutput_()) {
205  Opm::Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
206  fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
207  }
208  if (relativePermeabilityOutput_()) {
209  Opm::Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
210  fractureRelativePermeability_[phaseIdx][I] =
211  intQuants.fractureRelativePermeability(phaseIdx);
212  }
213  if (volumeFractionOutput_()) {
214  Opm::Valgrind::CheckDefined(intQuants.fractureVolume());
215  fractureVolumeFraction_[I] += intQuants.fractureVolume();
216  }
217  }
218  }
219 
220  if (velocityOutput_()) {
221  // calculate velocities if requested by the simulator
222  for (unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++ scvfIdx) {
223  const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, /*timeIdx=*/0);
224 
225  unsigned i = extQuants.interiorIndex();
226  unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
227 
228  unsigned j = extQuants.exteriorIndex();
229  unsigned J = elemCtx.globalSpaceIndex(j, /*timeIdx=*/0);
230 
231  if (!fractureMapper.isFractureEdge(I, J))
232  continue;
233 
234  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
235  Scalar weight =
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();
240 
241  Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
242  v *= weight;
243 
244  for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
245  fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
246  fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
247  }
248 
249  fractureVelocityWeight_[phaseIdx][I] += weight;
250  fractureVelocityWeight_[phaseIdx][J] += weight;
251  }
252  }
253  }
254  }
255 
259  void commitBuffers(BaseOutputWriter& baseWriter)
260  {
261  VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
262  if (!vtkWriter) {
263  return;
264  }
265 
266  if (saturationOutput_())
267  this->commitPhaseBuffer_(baseWriter, "fractureSaturation_%s", fractureSaturation_);
268  if (mobilityOutput_())
269  this->commitPhaseBuffer_(baseWriter, "fractureMobility_%s", fractureMobility_);
270  if (relativePermeabilityOutput_())
271  this->commitPhaseBuffer_(baseWriter, "fractureRelativePerm_%s", fractureRelativePermeability_);
272 
273  if (porosityOutput_())
274  this->commitScalarBuffer_(baseWriter, "fracturePorosity", fracturePorosity_);
275  if (intrinsicPermeabilityOutput_())
276  this->commitScalarBuffer_(baseWriter, "fractureIntrinsicPerm", fractureIntrinsicPermeability_);
277  if (volumeFractionOutput_()) {
278  // divide the fracture volume by the total volume of the finite volumes
279  for (unsigned I = 0; I < fractureVolumeFraction_.size(); ++I)
280  fractureVolumeFraction_[I] /= this->simulator_.model().dofTotalVolume(I);
281  this->commitScalarBuffer_(baseWriter, "fractureVolumeFraction", fractureVolumeFraction_);
282  }
283 
284  if (velocityOutput_()) {
285  size_t nDof = this->simulator_.model().numGridDof();
286 
287  for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
288  // first, divide the velocity field by the
289  // respective finite volume's surface area
290  for (unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx)
291  fractureVelocity_[phaseIdx][dofIdx] /=
292  std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
293  // commit the phase velocity
294  char name[512];
295  snprintf(name, 512, "fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx));
296 
297  DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
298  }
299  }
300  }
301 
302 private:
303  static bool saturationOutput_()
304  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureSaturations); }
305 
306  static bool mobilityOutput_()
307  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureMobilities); }
308 
309  static bool relativePermeabilityOutput_()
310  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureRelativePermeabilities); }
311 
312  static bool porosityOutput_()
313  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFracturePorosity); }
314 
315  static bool intrinsicPermeabilityOutput_()
316  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureIntrinsicPermeabilities); }
317 
318  static bool volumeFractionOutput_()
319  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureVolumeFraction); }
320 
321  static bool velocityOutput_()
322  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWriteFractureFilterVelocities); }
323 
324  PhaseBuffer fractureSaturation_;
325  PhaseBuffer fractureMobility_;
326  PhaseBuffer fractureRelativePermeability_;
327 
328  ScalarBuffer fracturePorosity_;
329  ScalarBuffer fractureVolumeFraction_;
330  ScalarBuffer fractureIntrinsicPermeability_;
331 
332  PhaseVectorBuffer fractureVelocity_;
333  PhaseBuffer fractureVelocityWeight_;
334 
335  PhaseVectorBuffer potentialGradient_;
336  PhaseBuffer potentialWeight_;
337 };
338 
339 } // namespace Ewoms
340 
341 #endif
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