vtkcompositionmodule.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_COMPOSITION_MODULE_HPP
28#define OPM_VTK_COMPOSITION_MODULE_HPP
29
30#include <opm/material/common/MathToolbox.hpp>
31
33
37
40
41namespace Opm {
42
55template <class TypeTag>
57{
59
64
66
67 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
68 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
69
70 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
72
73 using BufferType = typename ParentType::BufferType;
76
77public:
78 explicit VtkCompositionModule(const Simulator& simulator)
79 : ParentType(simulator)
80 {
81 params_.read();
82 }
83
87 static void registerParameters()
88 {
90 }
91
96 void allocBuffers() override
97 {
98 if (params_.moleFracOutput_) {
99 this->resizePhaseComponentBuffer_(moleFrac_, BufferType::Dof);
100 }
101 if (params_.massFracOutput_) {
102 this->resizePhaseComponentBuffer_(massFrac_, BufferType::Dof);
103 }
104 if (params_.totalMassFracOutput_) {
105 this->resizeComponentBuffer_(totalMassFrac_, BufferType::Dof);
106 }
107 if (params_.totalMoleFracOutput_) {
108 this->resizeComponentBuffer_(totalMoleFrac_, BufferType::Dof);
109 }
110 if (params_.molarityOutput_) {
111 this->resizePhaseComponentBuffer_(molarity_, BufferType::Dof);
112 }
113
114 if (params_.fugacityOutput_) {
115 this->resizeComponentBuffer_(fugacity_, BufferType::Dof);
116 }
117 if (params_.fugacityCoeffOutput_) {
118 this->resizePhaseComponentBuffer_(fugacityCoeff_, BufferType::Dof);
119 }
120 }
121
126 void processElement(const ElementContext& elemCtx) override
127 {
128 using Toolbox = MathToolbox<Evaluation>;
129
130 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
131 return;
132 }
133
134 for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
135 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
136 const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
137 const auto& fs = intQuants.fluidState();
138
139 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
140 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
141 if (params_.moleFracOutput_) {
142 moleFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
143 }
144 if (params_.massFracOutput_) {
145 massFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.massFraction(phaseIdx, compIdx));
146 }
147 if (params_.molarityOutput_) {
148 molarity_[phaseIdx][compIdx][I] = Toolbox::value(fs.molarity(phaseIdx, compIdx));
149 }
150
151 if (params_.fugacityCoeffOutput_) {
152 fugacityCoeff_[phaseIdx][compIdx][I] =
153 Toolbox::value(fs.fugacityCoefficient(phaseIdx, compIdx));
154 }
155 }
156 }
157
158 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
159 if (params_.totalMassFracOutput_) {
160 Scalar compMass = 0;
161 Scalar totalMass = 0;
162 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
163 totalMass += Toolbox::value(fs.density(phaseIdx)) * Toolbox::value(fs.saturation(phaseIdx));
164 compMass +=
165 Toolbox::value(fs.density(phaseIdx)) *
166 Toolbox::value(fs.saturation(phaseIdx)) *
167 Toolbox::value(fs.massFraction(phaseIdx, compIdx));
168 }
169 totalMassFrac_[compIdx][I] = compMass / totalMass;
170 }
171 if (params_.totalMoleFracOutput_) {
172 Scalar compMoles = 0;
173 Scalar totalMoles = 0;
174 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
175 totalMoles +=
176 Toolbox::value(fs.molarDensity(phaseIdx)) *
177 Toolbox::value(fs.saturation(phaseIdx));
178 compMoles +=
179 Toolbox::value(fs.molarDensity(phaseIdx)) *
180 Toolbox::value(fs.saturation(phaseIdx)) *
181 Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
182 }
183 totalMoleFrac_[compIdx][I] = compMoles / totalMoles;
184 }
185 if (params_.fugacityOutput_) {
186 fugacity_[compIdx][I] = Toolbox::value(intQuants.fluidState().fugacity(/*phaseIdx=*/0, compIdx));
187 }
188 }
189 }
190 }
191
195 void commitBuffers(BaseOutputWriter& baseWriter) override
196 {
197 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
198 return;
199 }
200
201 if (params_.moleFracOutput_) {
202 this->commitPhaseComponentBuffer_(baseWriter, "moleFrac_%s^%s",
203 moleFrac_, BufferType::Dof);
204 }
205 if (params_.massFracOutput_) {
206 this->commitPhaseComponentBuffer_(baseWriter, "massFrac_%s^%s",
207 massFrac_, BufferType::Dof);
208 }
209 if (params_.molarityOutput_) {
210 this->commitPhaseComponentBuffer_(baseWriter, "molarity_%s^%s",
211 molarity_, BufferType::Dof);
212 }
213 if (params_.totalMassFracOutput_) {
214 this->commitComponentBuffer_(baseWriter, "totalMassFrac^%s",
215 totalMassFrac_, BufferType::Dof);
216 }
217 if (params_.totalMoleFracOutput_) {
218 this->commitComponentBuffer_(baseWriter, "totalMoleFrac^%s",
219 totalMoleFrac_, BufferType::Dof);
220 }
221
222 if (params_.fugacityOutput_) {
223 this->commitComponentBuffer_(baseWriter, "fugacity^%s",
224 fugacity_, BufferType::Dof);
225 }
226 if (params_.fugacityCoeffOutput_) {
227 this->commitPhaseComponentBuffer_(baseWriter, "fugacityCoeff_%s^%s",
228 fugacityCoeff_, BufferType::Dof);
229 }
230 }
231
232private:
233 VtkCompositionParams params_{};
234 PhaseComponentBuffer moleFrac_{};
235 PhaseComponentBuffer massFrac_{};
236 PhaseComponentBuffer molarity_{};
237 ComponentBuffer totalMassFrac_{};
238 ComponentBuffer totalMoleFrac_{};
239
240 ComponentBuffer fugacity_{};
241 PhaseComponentBuffer fugacityCoeff_{};
242};
243
244} // namespace Opm
245
246#endif // OPM_VTK_COMPOSITION_MODULE_HPP
The base class for writer modules.
Definition: baseoutputmodule.hh:68
void commitPhaseComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseComponentBuffer &buffer, BufferType bufferType)
Add a phase and component specific quantities to the output.
Definition: baseoutputmodule.hh:369
void commitComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, ComponentBuffer &buffer, BufferType bufferType)
Add a component-specific buffer to the result file.
Definition: baseoutputmodule.hh:353
std::array< ScalarBuffer, numComponents > ComponentBuffer
Definition: baseoutputmodule.hh:92
BufferType
Definition: baseoutputmodule.hh:143
std::array< std::array< ScalarBuffer, numComponents >, numPhases > PhaseComponentBuffer
Definition: baseoutputmodule.hh:93
void resizePhaseComponentBuffer_(PhaseComponentBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a phase and component specific buffer.
Definition: baseoutputmodule.hh:224
void resizeComponentBuffer_(ComponentBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a component specific quantity.
Definition: baseoutputmodule.hh:211
The base class for all output writers.
Definition: baseoutputwriter.hh:46
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hpp:57
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkcompositionmodule.hpp:126
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkcompositionmodule.hpp:195
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkcompositionmodule.hpp:96
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hpp:87
VtkCompositionModule(const Simulator &simulator)
Definition: vtkcompositionmodule.hpp:78
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 VtkCompositionModule.
Definition: vtkcompositionparams.hpp:49
bool fugacityCoeffOutput_
Definition: vtkcompositionparams.hpp:62
bool totalMoleFracOutput_
Definition: vtkcompositionparams.hpp:59
bool moleFracOutput_
Definition: vtkcompositionparams.hpp:57
bool molarityOutput_
Definition: vtkcompositionparams.hpp:60
bool massFracOutput_
Definition: vtkcompositionparams.hpp:56
bool totalMassFracOutput_
Definition: vtkcompositionparams.hpp:58
bool fugacityOutput_
Definition: vtkcompositionparams.hpp:61
static void registerParameters()
Registers the parameters in parameter system.
void read()
Reads the parameter values from the parameter system.