vtkblackoilmicpmodule.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_MICP_MODULE_HPP
28#define OPM_VTK_BLACK_OIL_MICP_MODULE_HPP
29
30#include <dune/common/fvector.hh>
31
32#include <opm/material/densead/Math.hpp>
33
35
37
41
44
45namespace Opm {
46
52template <class TypeTag>
54{
56
62
63 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
65
66 enum { enableMICP = getPropValue<TypeTag, Properties::EnableMICP>() };
67
68 using BufferType = typename ParentType::BufferType;
70
71public:
72 explicit VtkBlackOilMICPModule(const Simulator& simulator)
73 : ParentType(simulator)
74 {
75 if constexpr (enableMICP) {
76 params_.read();
77 }
78 }
79
84 static void registerParameters()
85 {
86 if constexpr (enableMICP) {
88 }
89 }
90
95 void allocBuffers() override
96 {
97 if constexpr (enableMICP) {
98 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
99 return;
100 }
101
102 if (params_.microbialConcentrationOutput_) {
103 this->resizeScalarBuffer_(microbialConcentration_, BufferType::Dof);
104 }
105 if (params_.oxygenConcentrationOutput_) {
106 this->resizeScalarBuffer_(oxygenConcentration_, BufferType::Dof);
107 }
108 if (params_.ureaConcentrationOutput_) {
109 this->resizeScalarBuffer_(ureaConcentration_, BufferType::Dof);
110 }
111 if (params_.biofilmConcentrationOutput_) {
112 this->resizeScalarBuffer_(biofilmConcentration_, BufferType::Dof);
113 }
114 if (params_.calciteConcentrationOutput_) {
115 this->resizeScalarBuffer_(calciteConcentration_, BufferType::Dof);
116 }
117 }
118 }
119
124 void processElement(const ElementContext& elemCtx) override
125 {
126 if constexpr (enableMICP) {
127 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
128 return;
129 }
130
131 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
132 const auto& intQuants = elemCtx.intensiveQuantities(dofIdx, /*timeIdx=*/0);
133 const unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
134
135 if (params_.microbialConcentrationOutput_) {
136 microbialConcentration_[globalDofIdx] =
137 scalarValue(intQuants.microbialConcentration());
138 }
139
140 if (params_.oxygenConcentrationOutput_) {
141 oxygenConcentration_[globalDofIdx] =
142 scalarValue(intQuants.oxygenConcentration());
143 }
144
145 if (params_.ureaConcentrationOutput_) {
146 ureaConcentration_[globalDofIdx] =
147 scalarValue(intQuants.ureaConcentration());
148 }
149
150 if (params_.biofilmConcentrationOutput_) {
151 biofilmConcentration_[globalDofIdx] =
152 scalarValue(intQuants.biofilmConcentration());
153 }
154
155 if (params_.calciteConcentrationOutput_) {
156 calciteConcentration_[globalDofIdx] =
157 scalarValue(intQuants.calciteConcentration());
158 }
159 }
160 }
161 }
162
166 void commitBuffers(BaseOutputWriter& baseWriter) override
167 {
168 if constexpr (enableMICP) {
169 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
170 return;
171 }
172
173 if (params_.microbialConcentrationOutput_) {
174 this->commitScalarBuffer_(baseWriter, "microbial concentration",
175 microbialConcentration_, BufferType::Dof);
176 }
177
178 if (params_.oxygenConcentrationOutput_) {
179 this->commitScalarBuffer_(baseWriter, "oxygen concentration",
180 oxygenConcentration_, BufferType::Dof);
181 }
182
183 if (params_.ureaConcentrationOutput_) {
184 this->commitScalarBuffer_(baseWriter, "urea concentration",
185 ureaConcentration_, BufferType::Dof);
186 }
187
188 if (params_.biofilmConcentrationOutput_) {
189 this->commitScalarBuffer_(baseWriter, "biofilm fraction",
190 biofilmConcentration_, BufferType::Dof);
191 }
192
193 if (params_.calciteConcentrationOutput_) {
194 this->commitScalarBuffer_(baseWriter, "calcite fraction",
195 calciteConcentration_, BufferType::Dof);
196 }
197 }
198 }
199
200private:
201 VtkBlackoilMICPParams params_{};
202 ScalarBuffer microbialConcentration_{};
203 ScalarBuffer oxygenConcentration_{};
204 ScalarBuffer ureaConcentration_{};
205 ScalarBuffer biofilmConcentration_{};
206 ScalarBuffer calciteConcentration_{};
207};
208
209} // namespace Opm
210
211#endif // OPM_VTK_BLACKOIL_MICP_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 MICP model's related quantities.
Definition: vtkblackoilmicpmodule.hpp:54
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkblackoilmicpmodule.hpp:124
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkblackoilmicpmodule.hpp:166
VtkBlackOilMICPModule(const Simulator &simulator)
Definition: vtkblackoilmicpmodule.hpp:72
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilmicpmodule.hpp:95
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilmicpmodule.hpp:84
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
Declare the properties used by the infrastructure code of the finite volume discretizations.
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.
Struct holding the parameters for VtkBlackoilMICPModule.
Definition: vtkblackoilmicpparams.hpp:47
bool oxygenConcentrationOutput_
Definition: vtkblackoilmicpparams.hpp:55
static void registerParameters()
Registers the parameters in parameter system.
void read()
Reads the parameter values from the parameter system.
bool ureaConcentrationOutput_
Definition: vtkblackoilmicpparams.hpp:56
bool calciteConcentrationOutput_
Definition: vtkblackoilmicpparams.hpp:58
bool microbialConcentrationOutput_
Definition: vtkblackoilmicpparams.hpp:54
bool biofilmConcentrationOutput_
Definition: vtkblackoilmicpparams.hpp:57