27 #ifndef EWOMS_VTK_COMPOSITION_MODULE_HH 28 #define EWOMS_VTK_COMPOSITION_MODULE_HH 36 #include <opm/material/common/MathToolbox.hpp> 39 namespace Properties {
57 SET_BOOL_PROP(VtkComposition, VtkWriteTotalMassFractions,
false);
58 SET_BOOL_PROP(VtkComposition, VtkWriteTotalMoleFractions,
false);
76 template <
class TypeTag>
82 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
83 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
84 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
86 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
91 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
94 typedef typename ParentType::ComponentBuffer ComponentBuffer;
95 typedef typename ParentType::PhaseComponentBuffer PhaseComponentBuffer;
108 "Include mass fractions in the VTK output files");
110 "Include mole fractions in the VTK output files");
112 "Include total mass fractions in the VTK output files");
114 "Include total mole fractions in the VTK output files");
116 "Include component molarities in the VTK output files");
118 "Include component fugacities in the VTK output files");
120 "Include component fugacity coefficients in the VTK output files");
129 if (moleFracOutput_())
131 if (massFracOutput_())
133 if (totalMassFracOutput_())
135 if (totalMoleFracOutput_())
137 if (molarityOutput_())
140 if (fugacityOutput_())
142 if (fugacityCoeffOutput_())
152 typedef Opm::MathToolbox<Evaluation> Toolbox;
157 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
158 unsigned I = elemCtx.globalSpaceIndex(i, 0);
159 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
160 const auto& fs = intQuants.fluidState();
162 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
163 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
164 if (moleFracOutput_())
165 moleFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
166 if (massFracOutput_())
167 massFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.massFraction(phaseIdx, compIdx));
168 if (molarityOutput_())
169 molarity_[phaseIdx][compIdx][I] = Toolbox::value(fs.molarity(phaseIdx, compIdx));
171 if (fugacityCoeffOutput_())
172 fugacityCoeff_[phaseIdx][compIdx][I] =
173 Toolbox::value(fs.fugacityCoefficient(phaseIdx, compIdx));
177 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
178 if (totalMassFracOutput_()) {
180 Scalar totalMass = 0;
181 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
182 totalMass += Toolbox::value(fs.density(phaseIdx)) * Toolbox::value(fs.saturation(phaseIdx));
184 Toolbox::value(fs.density(phaseIdx))
185 *Toolbox::value(fs.saturation(phaseIdx))
186 *Toolbox::value(fs.massFraction(phaseIdx, compIdx));
188 totalMassFrac_[compIdx][I] = compMass / totalMass;
190 if (totalMoleFracOutput_()) {
191 Scalar compMoles = 0;
192 Scalar totalMoles = 0;
193 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
195 Toolbox::value(fs.molarDensity(phaseIdx))
196 *Toolbox::value(fs.saturation(phaseIdx));
198 Toolbox::value(fs.molarDensity(phaseIdx))
199 *Toolbox::value(fs.saturation(phaseIdx))
200 *Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
202 totalMoleFrac_[compIdx][I] = compMoles / totalMoles;
204 if (fugacityOutput_())
205 fugacity_[compIdx][I] = Toolbox::value(intQuants.fluidState().fugacity(0, compIdx));
220 if (moleFracOutput_())
222 if (massFracOutput_())
224 if (molarityOutput_())
226 if (totalMassFracOutput_())
228 if (totalMoleFracOutput_())
231 if (fugacityOutput_())
233 if (fugacityCoeffOutput_())
238 static bool massFracOutput_()
241 static bool moleFracOutput_()
244 static bool totalMassFracOutput_()
247 static bool totalMoleFracOutput_()
250 static bool molarityOutput_()
253 static bool fugacityOutput_()
256 static bool fugacityCoeffOutput_()
259 PhaseComponentBuffer moleFrac_;
260 PhaseComponentBuffer massFrac_;
261 PhaseComponentBuffer molarity_;
262 ComponentBuffer totalMassFrac_;
263 ComponentBuffer totalMoleFrac_;
265 ComponentBuffer fugacity_;
266 PhaseComponentBuffer fugacityCoeff_;
The base class for all output writers.
Definition: baseoutputwriter.hh:43
#define SET_BOOL_PROP(EffTypeTagName, PropTagName,...)
Set a property to a simple constant boolean value.
Definition: propertysystem.hh:361
Definition: baseauxiliarymodule.hh:37
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:105
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
The base class for writer modules.
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkcompositionmodule.hh:127
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkcompositionmodule.hh:213
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:68
This file provides the infrastructure to retrieve run-time parameters.
void resizeComponentBuffer_(ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a component specific quantity.
Definition: baseoutputmodule.hh:258
The base class for writer modules.
Definition: baseoutputmodule.hh:80
void resizePhaseComponentBuffer_(PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase and component specific buffer.
Definition: baseoutputmodule.hh:281
void commitComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a component-specific buffer to the result file.
Definition: baseoutputmodule.hh:431
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
Provides the magic behind the eWoms property system.
Simplifies writing multi-file VTK datasets.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:77
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element...
Definition: vtkcompositionmodule.hh:150
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:454