vtkblackoilenergymodule.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_ENERGY_MODULE_HPP
28#define OPM_VTK_BLACK_OIL_ENERGY_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
63
64 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
66
67 enum { enableEnergy = getPropValue<TypeTag, Properties::EnableEnergy>() };
68 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
69
70 using BufferType = typename ParentType::BufferType;
71 using PhaseBuffer = typename ParentType::PhaseBuffer;
73
74public:
75 explicit VtkBlackOilEnergyModule(const Simulator& simulator)
76 : ParentType(simulator)
77 {
78 if constexpr (enableEnergy) {
79 params_.read();
80 }
81 }
82
87 static void registerParameters()
88 {
89 if constexpr (enableEnergy) {
91 }
92 }
93
98 void allocBuffers() override
99 {
100 if constexpr (enableEnergy) {
101 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
102 return;
103 }
104
105 if (params_.rockInternalEnergyOutput_) {
106 this->resizeScalarBuffer_(rockInternalEnergy_, BufferType::Dof);
107 }
109 this->resizeScalarBuffer_(totalThermalConductivity_, BufferType::Dof);
110 }
111 if (params_.fluidInternalEnergiesOutput_) {
112 this->resizePhaseBuffer_(fluidInternalEnergies_, BufferType::Dof);
113 }
114 if (params_.fluidEnthalpiesOutput_) {
115 this->resizePhaseBuffer_(fluidEnthalpies_, BufferType::Dof);
116 }
117 }
118 }
119
124 void processElement(const ElementContext& elemCtx) override
125 {
126 if constexpr (enableEnergy) {
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_.rockInternalEnergyOutput_) {
136 rockInternalEnergy_[globalDofIdx] =
137 scalarValue(intQuants.rockInternalEnergy());
138 }
139
141 totalThermalConductivity_[globalDofIdx] =
142 scalarValue(intQuants.totalThermalConductivity());
143 }
144
145 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
146 if (FluidSystem::phaseIsActive(phaseIdx)) {
147 if (params_.fluidInternalEnergiesOutput_) {
148 fluidInternalEnergies_[phaseIdx][globalDofIdx] =
149 scalarValue(intQuants.fluidState().internalEnergy(phaseIdx));
150 }
151
152 if (params_.fluidEnthalpiesOutput_) {
153 fluidEnthalpies_[phaseIdx][globalDofIdx] =
154 scalarValue(intQuants.fluidState().enthalpy(phaseIdx));
155 }
156 }
157 }
158 }
159 }
160 }
161
165 void commitBuffers(BaseOutputWriter& baseWriter) override
166 {
167 if constexpr (enableEnergy) {
168 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
169 return;
170 }
171
172 if (params_.rockInternalEnergyOutput_) {
173 this->commitScalarBuffer_(baseWriter, "volumetric internal energy rock",
174 rockInternalEnergy_, BufferType::Dof);
175 }
176
178 this->commitScalarBuffer_(baseWriter, "total thermal conductivity",
179 totalThermalConductivity_, BufferType::Dof);
180 }
181
182 if (params_.fluidInternalEnergiesOutput_) {
183 this->commitPhaseBuffer_(baseWriter, "internal energy_%s",
184 fluidInternalEnergies_, BufferType::Dof);
185 }
186
187 if (params_.fluidEnthalpiesOutput_) {
188 this->commitPhaseBuffer_(baseWriter, "enthalpy_%s",
189 fluidEnthalpies_, BufferType::Dof);
190 }
191 }
192 }
193
194private:
195 VtkBlackoilEnergyParams params_{};
196
197 ScalarBuffer rockInternalEnergy_{};
198 ScalarBuffer totalThermalConductivity_{};
199 PhaseBuffer fluidInternalEnergies_{};
200 PhaseBuffer fluidEnthalpies_{};
201};
202
203} // namespace Opm
204
205#endif // OPM_VTK_BLACKOIL_ENERGY_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
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:91
BufferType
Definition: baseoutputmodule.hh:143
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:337
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
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:198
The base class for all output writers.
Definition: baseoutputwriter.hh:46
VTK output module for the black oil model's energy related quantities.
Definition: vtkblackoilenergymodule.hpp:54
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkblackoilenergymodule.hpp:124
VtkBlackOilEnergyModule(const Simulator &simulator)
Definition: vtkblackoilenergymodule.hpp:75
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilenergymodule.hpp:98
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkblackoilenergymodule.hpp:165
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilenergymodule.hpp:87
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 VtkBlackoilEnergyOutputModule.
Definition: vtkblackoilenergyparams.hpp:46
bool fluidEnthalpiesOutput_
Definition: vtkblackoilenergyparams.hpp:56
bool rockInternalEnergyOutput_
Definition: vtkblackoilenergyparams.hpp:53
void read()
Reads the parameter values from the parameter system.
bool totalThermalConductivityOutput_
Definition: vtkblackoilenergyparams.hpp:54
static void registerParameters()
Registers the parameters in parameter system.
bool fluidInternalEnergiesOutput_
Definition: vtkblackoilenergyparams.hpp:55