vtkdiffusionmodule.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*/
28#ifndef OPM_VTK_DIFFUSION_MODULE_HPP
29#define OPM_VTK_DIFFUSION_MODULE_HPP
30
31#include <opm/material/densead/Evaluation.hpp>
32#include <opm/material/densead/Math.hpp>
33
35
39
42
43namespace Opm {
44
56template <class TypeTag>
57class VtkDiffusionModule : public BaseOutputModule<TypeTag>
58{
60
65
66 using Toolbox = MathToolbox<Evaluation>;
67
68 using BufferType = typename ParentType::BufferType;
69 using PhaseBuffer = typename ParentType::PhaseBuffer;
71
72 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
74
75 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
76 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
77
78public:
79 explicit VtkDiffusionModule(const Simulator& simulator)
80 : ParentType(simulator)
81 {
82 params_.read();
83 }
84
88 static void registerParameters()
89 {
91 }
92
97 void allocBuffers() override
98 {
99 if (params_.tortuosityOutput_) {
100 this->resizePhaseBuffer_(tortuosity_, BufferType::Dof);
101 }
102 if (params_.diffusionCoefficientOutput_) {
103 this->resizePhaseComponentBuffer_(diffusionCoefficient_, BufferType::Dof);
104 }
106 this->resizePhaseComponentBuffer_(effectiveDiffusionCoefficient_, BufferType::Dof);
107 }
108 }
109
114 void processElement(const ElementContext& elemCtx) override
115 {
116 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
117 return;
118 }
119
120 for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
121 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
122 const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
123
124 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
125 if (params_.tortuosityOutput_) {
126 tortuosity_[phaseIdx][I] = Toolbox::value(intQuants.tortuosity(phaseIdx));
127 }
128 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
129 if (params_.diffusionCoefficientOutput_) {
130 diffusionCoefficient_[phaseIdx][compIdx][I] =
131 Toolbox::value(intQuants.diffusionCoefficient(phaseIdx, compIdx));
132 }
134 effectiveDiffusionCoefficient_[phaseIdx][compIdx][I] =
135 Toolbox::value(intQuants.effectiveDiffusionCoefficient(phaseIdx, compIdx));
136 }
137 }
138 }
139 }
140 }
141
145 void commitBuffers(BaseOutputWriter& baseWriter) override
146 {
147 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
148 return;
149 }
150
151 if (params_.tortuosityOutput_) {
152 this->commitPhaseBuffer_(baseWriter, "tortuosity", tortuosity_, BufferType::Dof);
153 }
154 if (params_.diffusionCoefficientOutput_) {
155 this->commitPhaseComponentBuffer_(baseWriter, "diffusionCoefficient",
156 diffusionCoefficient_, BufferType::Dof);
157 }
159 this->commitPhaseComponentBuffer_(baseWriter,
160 "effectiveDiffusionCoefficient",
161 effectiveDiffusionCoefficient_,
162 BufferType::Dof);
163 }
164 }
165
166private:
167 VtkDiffusionParams params_{};
168 PhaseBuffer tortuosity_{};
169 PhaseComponentBuffer diffusionCoefficient_{};
170 PhaseComponentBuffer effectiveDiffusionCoefficient_{};
171};
172
173} // namespace Opm
174
175#endif // OPM_VTK_DIFFUSION_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
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:91
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 commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:337
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 quantities which make sense for models which incorperate molecular diffusion.
Definition: vtkdiffusionmodule.hpp:58
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiffusionmodule.hpp:97
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkdiffusionmodule.hpp:145
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quanties relevant for an element.
Definition: vtkdiffusionmodule.hpp:114
VtkDiffusionModule(const Simulator &simulator)
Definition: vtkdiffusionmodule.hpp:79
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkdiffusionmodule.hpp:88
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 VtkDiffusionModule.
Definition: vtkdiffusionparams.hpp:46
static void registerParameters()
Registers the parameters in parameter system.
bool diffusionCoefficientOutput_
Definition: vtkdiffusionparams.hpp:54
bool tortuosityOutput_
Definition: vtkdiffusionparams.hpp:53
bool effectiveDiffusionCoefficientOutput_
Definition: vtkdiffusionparams.hpp:55
void read()
Reads the parameter values from the parameter system.