vtkblackoilpolymermodule.hpp
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 OPM_VTK_BLACK_OIL_POLYMER_MODULE_HPP
28#define OPM_VTK_BLACK_OIL_POLYMER_MODULE_HPP
29
30#include <dune/common/fvector.hh>
31
32#include <opm/material/densead/Math.hpp>
33
35
37
41
44
45namespace Opm::Properties::TTag {
46
47// create new type tag for the VTK multi-phase output
49
50} // namespace Opm::Properties::TTag
51
52namespace Opm {
53
59template <class TypeTag>
61{
63
69
70 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
72
73 enum { enablePolymer = getPropValue<TypeTag, Properties::EnablePolymer>() };
74
75 using BufferType = typename ParentType::BufferType;
77
78public:
79 explicit VtkBlackOilPolymerModule(const Simulator& simulator)
80 : ParentType(simulator)
81 {
82 if constexpr (enablePolymer) {
83 params_.read();
84 }
85 }
86
91 static void registerParameters()
92 {
93 if constexpr (enablePolymer) {
95 }
96 }
97
102 void allocBuffers() override
103 {
104 if constexpr (enablePolymer) {
105 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
106 return;
107 }
108
109 if (params_.polymerConcentrationOutput_) {
110 this->resizeScalarBuffer_(polymerConcentration_, BufferType::Dof);
111 }
112 if (params_.polymerDeadPoreVolumeOutput_) {
113 this->resizeScalarBuffer_(polymerDeadPoreVolume_, BufferType::Dof);
114 }
115 if (params_.polymerRockDensityOutput_) {
116 this->resizeScalarBuffer_(polymerRockDensity_, BufferType::Dof);
117 }
118 if (params_.polymerAdsorptionOutput_) {
119 this->resizeScalarBuffer_(polymerAdsorption_, BufferType::Dof);
120 }
122 this->resizeScalarBuffer_(polymerViscosityCorrection_, BufferType::Dof);
123 }
125 this->resizeScalarBuffer_(waterViscosityCorrection_, BufferType::Dof);
126 }
127 }
128 }
129
134 void processElement(const ElementContext& elemCtx) override
135 {
136 if constexpr (enablePolymer) {
137 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
138 return;
139 }
140
141 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
142 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
143 const unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
144
145 if (params_.polymerConcentrationOutput_) {
146 polymerConcentration_[globalDofIdx] =
147 scalarValue(intQuants.polymerConcentration());
148 }
149
150 if (params_.polymerDeadPoreVolumeOutput_) {
151 polymerDeadPoreVolume_[globalDofIdx] =
152 scalarValue(intQuants.polymerDeadPoreVolume());
153 }
154
155 if (params_.polymerRockDensityOutput_) {
156 polymerRockDensity_[globalDofIdx] =
157 scalarValue(intQuants.polymerRockDensity());
158 }
159
160 if (params_.polymerAdsorptionOutput_) {
161 polymerAdsorption_[globalDofIdx] =
162 scalarValue(intQuants.polymerAdsorption());
163 }
164
166 polymerViscosityCorrection_[globalDofIdx] =
167 scalarValue(intQuants.polymerViscosityCorrection());
168 }
169
171 waterViscosityCorrection_[globalDofIdx] =
172 scalarValue(intQuants.waterViscosityCorrection());
173 }
174 }
175 }
176 }
177
181 void commitBuffers(BaseOutputWriter& baseWriter) override
182 {
183 if constexpr (enablePolymer) {
184 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
185 return;
186 }
187
188 if (params_.polymerConcentrationOutput_) {
189 this->commitScalarBuffer_(baseWriter, "polymer concentration",
190 polymerConcentration_, BufferType::Dof);
191 }
192
193 if (params_.polymerDeadPoreVolumeOutput_) {
194 this->commitScalarBuffer_(baseWriter, "dead pore volume fraction",
195 polymerDeadPoreVolume_, BufferType::Dof);
196 }
197
198 if (params_.polymerRockDensityOutput_) {
199 this->commitScalarBuffer_(baseWriter, "polymer rock density",
200 polymerRockDensity_, BufferType::Dof);
201 }
202
203 if (params_.polymerAdsorptionOutput_) {
204 this->commitScalarBuffer_(baseWriter, "polymer adsorption",
205 polymerAdsorption_, BufferType::Dof);
206 }
207
209 this->commitScalarBuffer_(baseWriter, "polymer viscosity correction",
210 polymerViscosityCorrection_, BufferType::Dof);
211 }
212
214 this->commitScalarBuffer_(baseWriter, "water viscosity correction",
215 waterViscosityCorrection_, BufferType::Dof);
216 }
217 }
218 }
219
220private:
221 VtkBlackoilPolymerParams params_{};
222 ScalarBuffer polymerConcentration_{};
223 ScalarBuffer polymerDeadPoreVolume_{};
224 ScalarBuffer polymerRockDensity_{};
225 ScalarBuffer polymerAdsorption_{};
226 ScalarBuffer polymerViscosityCorrection_{};
227 ScalarBuffer waterViscosityCorrection_{};
228};
229
230} // namespace Opm
231
232#endif // OPM_VTK_BLACK_OIL_POLYMER_MODULE_HPP
Declares the properties required by the black oil model.
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
BufferType
Definition: baseoutputmodule.hh:143
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
The base class for all output writers.
Definition: baseoutputwriter.hh:46
VTK output module for the black oil model's polymer related quantities.
Definition: vtkblackoilpolymermodule.hpp:61
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilpolymermodule.hpp:91
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkblackoilpolymermodule.hpp:134
VtkBlackOilPolymerModule(const Simulator &simulator)
Definition: vtkblackoilpolymermodule.hpp:79
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilpolymermodule.hpp:102
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkblackoilpolymermodule.hpp:181
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
Declare the properties used by the infrastructure code of the finite volume discretizations.
The generic type tag for problems using the immiscible multi-phase model.
Definition: blackoilmodel.hh:81
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.
Definition: vtkblackoilpolymermodule.hpp:48
Struct holding the parameters for VtkBlackoilPolymerModule.
Definition: vtkblackoilpolymerparams.hpp:48
bool polymerRockDensityOutput_
Definition: vtkblackoilpolymerparams.hpp:57
bool polymerViscosityCorrectionOutput_
Definition: vtkblackoilpolymerparams.hpp:59
bool polymerDeadPoreVolumeOutput_
Definition: vtkblackoilpolymerparams.hpp:56
void read()
Reads the parameter values from the parameter system.
static void registerParameters()
Registers the parameters in parameter system.
bool waterViscosityCorrectionOutput_
Definition: vtkblackoilpolymerparams.hpp:60
bool polymerAdsorptionOutput_
Definition: vtkblackoilpolymerparams.hpp:58
bool polymerConcentrationOutput_
Definition: vtkblackoilpolymerparams.hpp:55