All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
vtkblackoilpolymermodule.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_BLACK_OIL_POLYMER_MODULE_HH
28 #define EWOMS_VTK_BLACK_OIL_POLYMER_MODULE_HH
29 
30 #include <opm/material/densead/Math.hpp>
31 
32 #include "vtkmultiwriter.hh"
33 #include "baseoutputmodule.hh"
34 
38 
39 #include <dune/common/fvector.hh>
40 
41 #include <cstdio>
42 
43 namespace Ewoms {
44 namespace Properties {
45 // create new type tag for the VTK multi-phase output
46 NEW_TYPE_TAG(VtkBlackOilPolymer);
47 
48 // create the property tags needed for the polymer output module
49 NEW_PROP_TAG(EnablePolymer);
50 NEW_PROP_TAG(EnableVtkOutput);
51 NEW_PROP_TAG(VtkWritePolymerConcentration);
52 NEW_PROP_TAG(VtkWritePolymerDeadPoreVolume);
53 NEW_PROP_TAG(VtkWritePolymerAdsorption);
54 NEW_PROP_TAG(VtkWritePolymerRockDensity);
55 NEW_PROP_TAG(VtkWritePolymerViscosityCorrection);
56 NEW_PROP_TAG(VtkWriteWaterViscosityCorrection);
57 
58 // set default values for what quantities to output
59 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerConcentration, true);
60 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerDeadPoreVolume, true);
61 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerViscosityCorrection, true);
62 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWriteWaterViscosityCorrection, true);
63 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerRockDensity, true);
64 SET_BOOL_PROP(VtkBlackOilPolymer, VtkWritePolymerAdsorption, true);
65 } // namespace Properties
66 } // namespace Ewoms
67 
68 namespace Ewoms {
74 template <class TypeTag>
76 {
78 
79  typedef typename GET_PROP_TYPE(TypeTag, Simulator) Simulator;
80  typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
81  typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
82  typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
83  typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
84 
85  static const int vtkFormat = GET_PROP_VALUE(TypeTag, VtkOutputFormat);
87 
88  enum { enablePolymer = GET_PROP_VALUE(TypeTag, EnablePolymer) };
89 
90  typedef typename ParentType::ScalarBuffer ScalarBuffer;
91 
92 public:
93  VtkBlackOilPolymerModule(const Simulator& simulator)
94  : ParentType(simulator)
95  { }
96 
101  static void registerParameters()
102  {
103  if (!enablePolymer)
104  return;
105 
106  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerConcentration,
107  "Include the concentration of the polymer component in the water phase "
108  "in the VTK output files");
109  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerDeadPoreVolume,
110  "Include the fraction of the \"dead\" pore volume "
111  "in the VTK output files");
112  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerRockDensity,
113  "Include the amount of already adsorbed polymer component"
114  "in the VTK output files");
115  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerAdsorption,
116  "Include the adsorption rate of the polymer component"
117  "in the VTK output files");
118  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection,
119  "Include the viscosity correction of the polymer component "
120  "in the VTK output files");
121  EWOMS_REGISTER_PARAM(TypeTag, bool, VtkWriteWaterViscosityCorrection,
122  "Include the viscosity correction of the water component "
123  "due to polymers in the VTK output files");
124  }
125 
131  {
132  if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
133  return;
134 
135  if (!enablePolymer)
136  return;
137 
138  if (polymerConcentrationOutput_())
139  this->resizeScalarBuffer_(polymerConcentration_);
140  if (polymerDeadPoreVolumeOutput_())
141  this->resizeScalarBuffer_(polymerDeadPoreVolume_);
142  if (polymerRockDensityOutput_())
143  this->resizeScalarBuffer_(polymerRockDensity_);
144  if (polymerAdsorptionOutput_())
145  this->resizeScalarBuffer_(polymerAdsorption_);
146  if (polymerViscosityCorrectionOutput_())
147  this->resizeScalarBuffer_(polymerViscosityCorrection_);
148  if (waterViscosityCorrectionOutput_())
149  this->resizeScalarBuffer_(waterViscosityCorrection_);
150  }
151 
156  void processElement(const ElementContext& elemCtx)
157  {
158  if (!EWOMS_GET_PARAM(TypeTag, bool, EnableVtkOutput))
159  return;
160 
161  if (!enablePolymer)
162  return;
163 
164  for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
165  const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
166  unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
167 
168  if (polymerConcentrationOutput_())
169  polymerConcentration_[globalDofIdx] =
170  Opm::scalarValue(intQuants.polymerConcentration());
171 
172  if (polymerDeadPoreVolumeOutput_())
173  polymerDeadPoreVolume_[globalDofIdx] =
174  Opm::scalarValue(intQuants.polymerDeadPoreVolume());
175 
176  if (polymerRockDensityOutput_())
177  polymerRockDensity_[globalDofIdx] =
178  Opm::scalarValue(intQuants.polymerRockDensity());
179 
180  if (polymerAdsorptionOutput_())
181  polymerAdsorption_[globalDofIdx] =
182  Opm::scalarValue(intQuants.polymerAdsorption());
183 
184  if (polymerViscosityCorrectionOutput_())
185  polymerViscosityCorrection_[globalDofIdx] =
186  Opm::scalarValue(intQuants.polymerViscosityCorrection());
187 
188  if (waterViscosityCorrectionOutput_())
189  waterViscosityCorrection_[globalDofIdx] =
190  Opm::scalarValue(intQuants.waterViscosityCorrection());
191  }
192  }
193 
197  void commitBuffers(BaseOutputWriter& baseWriter)
198  {
199  VtkMultiWriter *vtkWriter = dynamic_cast<VtkMultiWriter*>(&baseWriter);
200  if (!vtkWriter)
201  return;
202 
203  if (!enablePolymer)
204  return;
205 
206  if (polymerConcentrationOutput_())
207  this->commitScalarBuffer_(baseWriter, "polymer concentration", polymerConcentration_);
208 
209  if (polymerDeadPoreVolumeOutput_())
210  this->commitScalarBuffer_(baseWriter, "dead pore volume fraction", polymerDeadPoreVolume_);
211 
212  if (polymerRockDensityOutput_())
213  this->commitScalarBuffer_(baseWriter, "polymer rock density", polymerRockDensity_);
214 
215  if (polymerAdsorptionOutput_())
216  this->commitScalarBuffer_(baseWriter, "polymer adsorption", polymerAdsorption_);
217 
218  if (polymerViscosityCorrectionOutput_())
219  this->commitScalarBuffer_(baseWriter, "polymer viscosity correction", polymerViscosityCorrection_);
220 
221  if (waterViscosityCorrectionOutput_())
222  this->commitScalarBuffer_(baseWriter, "water viscosity correction", waterViscosityCorrection_);
223  }
224 
225 private:
226  static bool polymerConcentrationOutput_()
227  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerConcentration); }
228 
229  static bool polymerDeadPoreVolumeOutput_()
230  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerDeadPoreVolume); }
231 
232  static bool polymerRockDensityOutput_()
233  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerRockDensity); }
234 
235  static bool polymerAdsorptionOutput_()
236  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerAdsorption); }
237 
238  static bool polymerViscosityCorrectionOutput_()
239  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection); }
240 
241  static bool waterViscosityCorrectionOutput_()
242  { return EWOMS_GET_PARAM(TypeTag, bool, VtkWritePolymerViscosityCorrection); }
243 
244  ScalarBuffer polymerConcentration_;
245  ScalarBuffer polymerDeadPoreVolume_;
246  ScalarBuffer polymerRockDensity_;
247  ScalarBuffer polymerAdsorption_;
248  ScalarBuffer polymerViscosityCorrection_;
249  ScalarBuffer waterViscosityCorrection_;
250 };
251 } // namespace Ewoms
252 
253 #endif
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
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
Declares the properties required by the black oil model.
This file provides the infrastructure to retrieve run-time parameters.
The base class for writer modules.
Definition: baseoutputmodule.hh:80
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilpolymermodule.hh:130
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilpolymermodule.hh:101
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.
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element...
Definition: vtkblackoilpolymermodule.hh:156
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
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkblackoilpolymermodule.hh:197
VTK output module for the black oil model&#39;s polymer related quantities.
Definition: vtkblackoilpolymermodule.hh:75