27#ifndef EWOMS_VTK_COMPOSITION_MODULE_HH
28#define EWOMS_VTK_COMPOSITION_MODULE_HH
36#include <opm/material/common/MathToolbox.hpp>
48template<
class TypeTag,
class MyTypeTag>
50template<
class TypeTag,
class MyTypeTag>
52template<
class TypeTag,
class MyTypeTag>
54template<
class TypeTag,
class MyTypeTag>
56template<
class TypeTag,
class MyTypeTag>
58template<
class TypeTag,
class MyTypeTag>
60template<
class TypeTag,
class MyTypeTag>
64template<
class TypeTag>
66template<
class TypeTag>
68template<
class TypeTag>
70template<
class TypeTag>
72template<
class TypeTag>
73struct VtkWriteMolarities<TypeTag, TTag::VtkComposition> {
static constexpr bool value =
false; };
74template<
class TypeTag>
75struct VtkWriteFugacities<TypeTag, TTag::VtkComposition> {
static constexpr bool value =
false; };
76template<
class TypeTag>
95template <
class TypeTag>
107 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
108 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
110 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
126 Parameters::registerParam<TypeTag, Properties::VtkWriteMassFractions>
127 (
"Include mass fractions in the VTK output files");
128 Parameters::registerParam<TypeTag, Properties::VtkWriteMoleFractions>
129 (
"Include mole fractions in the VTK output files");
130 Parameters::registerParam<TypeTag, Properties::VtkWriteTotalMassFractions>
131 (
"Include total mass fractions in the VTK output files");
132 Parameters::registerParam<TypeTag, Properties::VtkWriteTotalMoleFractions>
133 (
"Include total mole fractions in the VTK output files");
134 Parameters::registerParam<TypeTag, Properties::VtkWriteMolarities>
135 (
"Include component molarities in the VTK output files");
136 Parameters::registerParam<TypeTag, Properties::VtkWriteFugacities>
137 (
"Include component fugacities in the VTK output files");
138 Parameters::registerParam<TypeTag, Properties::VtkWriteFugacityCoeffs>
139 (
"Include component fugacity coefficients in the VTK output files");
148 if (moleFracOutput_())
150 if (massFracOutput_())
152 if (totalMassFracOutput_())
154 if (totalMoleFracOutput_())
156 if (molarityOutput_())
159 if (fugacityOutput_())
161 if (fugacityCoeffOutput_())
171 using Toolbox = MathToolbox<Evaluation>;
173 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
176 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
177 unsigned I = elemCtx.globalSpaceIndex(i, 0);
178 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
179 const auto& fs = intQuants.fluidState();
181 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
182 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
183 if (moleFracOutput_())
184 moleFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
185 if (massFracOutput_())
186 massFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.massFraction(phaseIdx, compIdx));
187 if (molarityOutput_())
188 molarity_[phaseIdx][compIdx][I] = Toolbox::value(fs.molarity(phaseIdx, compIdx));
190 if (fugacityCoeffOutput_())
191 fugacityCoeff_[phaseIdx][compIdx][I] =
192 Toolbox::value(fs.fugacityCoefficient(phaseIdx, compIdx));
196 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
197 if (totalMassFracOutput_()) {
199 Scalar totalMass = 0;
200 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
201 totalMass += Toolbox::value(fs.density(phaseIdx)) * Toolbox::value(fs.saturation(phaseIdx));
203 Toolbox::value(fs.density(phaseIdx))
204 *Toolbox::value(fs.saturation(phaseIdx))
205 *Toolbox::value(fs.massFraction(phaseIdx, compIdx));
207 totalMassFrac_[compIdx][I] = compMass / totalMass;
209 if (totalMoleFracOutput_()) {
210 Scalar compMoles = 0;
211 Scalar totalMoles = 0;
212 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
214 Toolbox::value(fs.molarDensity(phaseIdx))
215 *Toolbox::value(fs.saturation(phaseIdx));
217 Toolbox::value(fs.molarDensity(phaseIdx))
218 *Toolbox::value(fs.saturation(phaseIdx))
219 *Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
221 totalMoleFrac_[compIdx][I] = compMoles / totalMoles;
223 if (fugacityOutput_())
224 fugacity_[compIdx][I] = Toolbox::value(intQuants.fluidState().fugacity(0, compIdx));
239 if (moleFracOutput_())
241 if (massFracOutput_())
243 if (molarityOutput_())
245 if (totalMassFracOutput_())
247 if (totalMoleFracOutput_())
250 if (fugacityOutput_())
252 if (fugacityCoeffOutput_())
257 static bool massFracOutput_()
259 static bool val = Parameters::get<TypeTag, Properties::VtkWriteMassFractions>();
263 static bool moleFracOutput_()
265 static bool val = Parameters::get<TypeTag, Properties::VtkWriteMoleFractions>();
269 static bool totalMassFracOutput_()
271 static bool val = Parameters::get<TypeTag, Properties::VtkWriteTotalMassFractions>();
275 static bool totalMoleFracOutput_()
277 static bool val = Parameters::get<TypeTag, Properties::VtkWriteTotalMoleFractions>();
281 static bool molarityOutput_()
283 static bool val = Parameters::get<TypeTag, Properties::VtkWriteMolarities>();
287 static bool fugacityOutput_()
289 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFugacities>();
293 static bool fugacityCoeffOutput_()
295 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFugacityCoeffs>();
299 PhaseComponentBuffer moleFrac_;
300 PhaseComponentBuffer massFrac_;
301 PhaseComponentBuffer molarity_;
302 ComponentBuffer totalMassFrac_;
303 ComponentBuffer totalMoleFrac_;
305 ComponentBuffer fugacity_;
306 PhaseComponentBuffer fugacityCoeff_;
The base class for writer modules.
Definition: baseoutputmodule.hh:73
void commitComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a component-specific buffer to the result file.
Definition: baseoutputmodule.hh:442
std::array< ScalarBuffer, numComponents > ComponentBuffer
Definition: baseoutputmodule.hh:97
void resizeComponentBuffer_(ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a component specific quantity.
Definition: baseoutputmodule.hh:269
void resizePhaseComponentBuffer_(PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase and component specific buffer.
Definition: baseoutputmodule.hh:292
void commitPhaseComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase and component specific quantities to the output.
Definition: baseoutputmodule.hh:465
std::array< std::array< ScalarBuffer, numComponents >, numPhases > PhaseComponentBuffer
Definition: baseoutputmodule.hh:98
The base class for all output writers.
Definition: baseoutputwriter.hh:44
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:97
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkcompositionmodule.hh:169
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkcompositionmodule.hh:232
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:124
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkcompositionmodule.hh:146
VtkCompositionModule(const Simulator &simulator)
Definition: vtkcompositionmodule.hh:117
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:66
Definition: blackoilmodel.hh:72
Definition: blackoilboundaryratevector.hh:37
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:242
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Definition: vtkcompositionmodule.hh:43
a tag to mark properties as undefined
Definition: propertysystem.hh:40
Definition: vtkcompositionmodule.hh:59
Definition: vtkcompositionmodule.hh:61
Definition: vtkcompositionmodule.hh:49
Definition: vtkcompositionmodule.hh:57
Definition: vtkcompositionmodule.hh:51
Definition: vtkcompositionmodule.hh:53
Definition: vtkcompositionmodule.hh:55