25 #ifndef EWOMS_VTK_COMPOSITION_MODULE_HH
26 #define EWOMS_VTK_COMPOSITION_MODULE_HH
34 #include <opm/material/common/MathToolbox.hpp>
37 namespace Properties {
54 SET_BOOL_PROP(VtkComposition, VtkWriteTotalMassFractions,
false);
55 SET_BOOL_PROP(VtkComposition, VtkWriteTotalMoleFractions,
false);
73 template <
class TypeTag>
79 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
80 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
81 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
83 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
88 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
96 : ParentType(simulator)
105 "Include mass fractions in the VTK output files");
107 "Include mole fractions in the VTK output files");
109 "Include total mass fractions in the VTK output files");
111 "Include total mole fractions in the VTK output files");
113 "Include component molarities in the VTK output files");
115 "Include component fugacities in the VTK output files");
117 "Include component fugacity coefficients in the VTK output files");
126 if (moleFracOutput_())
128 if (massFracOutput_())
130 if (totalMassFracOutput_())
132 if (totalMoleFracOutput_())
134 if (molarityOutput_())
137 if (fugacityOutput_())
139 if (fugacityCoeffOutput_())
149 typedef Opm::MathToolbox<Evaluation> Toolbox;
151 for (
int i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
152 int I = elemCtx.globalSpaceIndex(i, 0);
153 const auto &intQuants = elemCtx.intensiveQuantities(i, 0);
154 const auto &fs = intQuants.fluidState();
156 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
157 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
158 if (moleFracOutput_())
159 moleFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
160 if (massFracOutput_())
161 massFrac_[phaseIdx][compIdx][I] = Toolbox::value(fs.massFraction(phaseIdx, compIdx));
162 if (molarityOutput_())
163 molarity_[phaseIdx][compIdx][I] = Toolbox::value(fs.molarity(phaseIdx, compIdx));
165 if (fugacityCoeffOutput_())
166 fugacityCoeff_[phaseIdx][compIdx][I] =
167 Toolbox::value(fs.fugacityCoefficient(phaseIdx, compIdx));
171 for (
int compIdx = 0; compIdx < numComponents; ++compIdx) {
172 if (totalMassFracOutput_()) {
174 Scalar totalMass = 0;
175 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
176 totalMass += Toolbox::value(fs.density(phaseIdx)) * Toolbox::value(fs.saturation(phaseIdx));
178 Toolbox::value(fs.density(phaseIdx))
179 *Toolbox::value(fs.saturation(phaseIdx))
180 *Toolbox::value(fs.massFraction(phaseIdx, compIdx));
182 totalMassFrac_[compIdx][I] = compMass / totalMass;
184 if (totalMoleFracOutput_()) {
185 Scalar compMoles = 0;
186 Scalar totalMoles = 0;
187 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
189 Toolbox::value(fs.molarDensity(phaseIdx))
190 *Toolbox::value(fs.saturation(phaseIdx));
192 Toolbox::value(fs.molarDensity(phaseIdx))
193 *Toolbox::value(fs.saturation(phaseIdx))
194 *Toolbox::value(fs.moleFraction(phaseIdx, compIdx));
196 totalMoleFrac_[compIdx][I] = compMoles / totalMoles;
198 if (fugacityOutput_())
199 fugacity_[compIdx][I] = Toolbox::value(intQuants.fluidState().fugacity(0, compIdx));
209 VtkMultiWriter *vtkWriter =
dynamic_cast<VtkMultiWriter*
>(&baseWriter);
214 if (moleFracOutput_())
216 if (massFracOutput_())
218 if (molarityOutput_())
220 if (totalMassFracOutput_())
222 if (totalMoleFracOutput_())
225 if (fugacityOutput_())
227 if (fugacityCoeffOutput_())
232 static bool massFracOutput_()
235 static bool moleFracOutput_()
238 static bool totalMassFracOutput_()
241 static bool totalMoleFracOutput_()
244 static bool molarityOutput_()
247 static bool fugacityOutput_()
250 static bool fugacityCoeffOutput_()
253 PhaseComponentBuffer moleFrac_;
254 PhaseComponentBuffer massFrac_;
255 PhaseComponentBuffer molarity_;
256 ComponentBuffer totalMassFrac_;
257 ComponentBuffer totalMoleFrac_;
259 ComponentBuffer fugacity_;
260 PhaseComponentBuffer fugacityCoeff_;
std::array< ScalarBuffer, numComponents > ComponentBuffer
Definition: baseoutputmodule.hh:96
void resizePhaseComponentBuffer_(PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase and component specific buffer.
Definition: baseoutputmodule.hh:259
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element...
Definition: vtkcompositionmodule.hh:147
void resizeComponentBuffer_(ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a component specific quantity.
Definition: baseoutputmodule.hh:236
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
std::array< std::array< ScalarBuffer, numComponents >, numPhases > PhaseComponentBuffer
Definition: baseoutputmodule.hh:97
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:432
NEW_PROP_TAG(Grid)
The type of the DUNE grid.
This file provides the infrastructure to retrieve run-time parameters.
void commitComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a component-specific buffer to the result file.
Definition: baseoutputmodule.hh:409
The base class for all output writers.
Definition: baseoutputwriter.hh:41
Simplifies writing multi-file VTK datasets.
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:73
Definition: baseauxiliarymodule.hh:35
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkcompositionmodule.hh:207
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:64
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkcompositionmodule.hh:124
static void registerParameters()
Register all run-time parameters for the Vtk output module.
Definition: vtkcompositionmodule.hh:102
VtkCompositionModule(const Simulator &simulator)
Definition: vtkcompositionmodule.hh:95
Provides the magic behind the eWoms property system.
The base class for writer modules.
SET_BOOL_PROP(FvBaseDiscretization, EnableVtkOutput, true)
Enable the VTK output by default.
The base class for writer modules.
Definition: baseoutputmodule.hh:71
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:61
VTK output module for the fluid composition.
Definition: vtkcompositionmodule.hh:74
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:95