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