VtkTracerModule.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_TRACER_MODULE_HPP
28#define OPM_VTK_TRACER_MODULE_HPP
29
30#include <dune/common/fvector.hh>
31
38
39#include <string>
40#include <vector>
41
42namespace Opm::Parameters {
43
44// set default values for what quantities to output
45struct VtkWriteTracerConcentration { static constexpr bool value = false; };
46
47} // namespace Opm::Parameters
48
49namespace Opm {
50
56 template <class TypeTag>
57 class VtkTracerModule : public BaseOutputModule<TypeTag>
58 {
60
65
68
69 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
71
72 using BufferType = typename ParentType::BufferType;
73 using ScalarBuffer = typename ParentType::ScalarBuffer;
74
75 public:
76 explicit VtkTracerModule(const Simulator& simulator)
77 : ParentType(simulator)
78 {}
79
84 static void registerParameters()
85 {
86 Parameters::Register<Parameters::VtkWriteTracerConcentration>
87 ("Include the tracer concentration in the VTK output files");
88 }
89
94 void allocBuffers() override
95 {
96 if (eclTracerConcentrationOutput_()) {
97 const auto& tracerModel = this->simulator_.problem().tracerModel();
98 eclFreeTracerConcentration_.resize(tracerModel.numTracers());
99 eclSolTracerConcentration_.resize(tracerModel.numTracers());
100 const auto& enableSolTracers = tracerModel.enableSolTracers();
101
102 for (std::size_t tracerIdx = 0; tracerIdx < eclFreeTracerConcentration_.size(); ++tracerIdx) {
103 this->resizeScalarBuffer_(eclFreeTracerConcentration_[tracerIdx], BufferType::Dof);
104 if (enableSolTracers[tracerIdx]) {
105 this->resizeScalarBuffer_(eclSolTracerConcentration_[tracerIdx], BufferType::Dof);
106 }
107 }
108 }
109 }
110
115 void processElement(const ElementContext& elemCtx) override
116 {
117 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
118 return;
119 }
120
121 if (eclTracerConcentrationOutput_()) {
122 const auto& tracerModel = elemCtx.problem().tracerModel();
123 const auto& enableSolTracers = tracerModel.enableSolTracers();
124
125 for (std::size_t tracerIdx = 0; tracerIdx < eclFreeTracerConcentration_.size(); ++tracerIdx) {
126 // free tracer
127 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
128 const unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
129 eclFreeTracerConcentration_[tracerIdx][globalDofIdx] =
130 tracerModel.freeTracerConcentration(tracerIdx, globalDofIdx);
131 }
132 // solution tracer (only if it exist)
133 if (enableSolTracers[tracerIdx]) {
134 for (unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++dofIdx) {
135 const unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, /*timeIdx=*/0);
136 eclSolTracerConcentration_[tracerIdx][globalDofIdx] =
137 tracerModel.solTracerConcentration(tracerIdx, globalDofIdx);
138 }
139 }
140 }
141 }
142 }
143
147 void commitBuffers(BaseOutputWriter& baseWriter) override
148 {
149 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter))
150 return;
151
152 if (eclTracerConcentrationOutput_()){
153 const auto& tracerModel = this->simulator_.problem().tracerModel();
154 const auto& enableSolTracers = tracerModel.enableSolTracers();
155
156 for (std::size_t tracerIdx = 0; tracerIdx < eclFreeTracerConcentration_.size(); ++tracerIdx) {
157 const std::string tmp = "freeTracerConcentration_" + tracerModel.name(tracerIdx);
158 this->commitScalarBuffer_(baseWriter, tmp.c_str(),
159 eclFreeTracerConcentration_[tracerIdx], BufferType::Dof);
160 if (enableSolTracers[tracerIdx]) {
161 const std::string tmp2 = "solTracerConcentration_" + tracerModel.name(tracerIdx);
162 this->commitScalarBuffer_(baseWriter, tmp2.c_str(),
163 eclSolTracerConcentration_[tracerIdx], BufferType::Dof);
164 }
165 }
166 }
167 }
168
169 private:
170 static bool eclTracerConcentrationOutput_()
171 {
172 static bool val = Parameters::Get<Parameters::VtkWriteTracerConcentration>();
173 return val;
174 }
175
176 std::vector<ScalarBuffer> eclFreeTracerConcentration_;
177 std::vector<ScalarBuffer> eclSolTracerConcentration_;
178 };
179
180} // namespace Opm
181
182#endif // OPM_VTK_TRACER_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
const Simulator & simulator_
Definition: baseoutputmodule.hh:426
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
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
VTK output module for the tracer model's parameters.
Definition: VtkTracerModule.hpp:58
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: VtkTracerModule.hpp:147
static void registerParameters()
Register all run-time parameters for the tracer VTK output module.
Definition: VtkTracerModule.hpp:84
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: VtkTracerModule.hpp:115
VtkTracerModule(const Simulator &simulator)
Definition: VtkTracerModule.hpp:76
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: VtkTracerModule.hpp:94
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilnewtonmethodparams.hpp:31
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: VtkTracerModule.hpp:45
static constexpr bool value
Definition: VtkTracerModule.hpp:45