27#ifndef OPM_VTK_BLACK_OIL_MODULE_HPP
28#define OPM_VTK_BLACK_OIL_MODULE_HPP
30#include <dune/common/fvector.hh>
32#include <opm/material/densead/Math.hpp>
51template <
class TypeTag>
64 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
67 enum { oilPhaseIdx = FluidSystem::oilPhaseIdx };
68 enum { gasPhaseIdx = FluidSystem::gasPhaseIdx };
69 enum { waterPhaseIdx = FluidSystem::waterPhaseIdx };
71 enum { gasCompIdx = FluidSystem::gasCompIdx };
72 enum { oilCompIdx = FluidSystem::oilCompIdx };
73 enum { waterCompIdx = FluidSystem::waterCompIdx };
143 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
147 for (
unsigned dofIdx = 0; dofIdx < elemCtx.numPrimaryDof(0); ++dofIdx) {
148 const auto& fs = elemCtx.intensiveQuantities(dofIdx, 0).fluidState();
149 using FluidState =
typename std::remove_const<
typename std::remove_reference<
decltype(fs)>::type>::type;
150 unsigned globalDofIdx = elemCtx.globalSpaceIndex(dofIdx, 0);
152 const auto& primaryVars = elemCtx.primaryVars(dofIdx, 0);
154 unsigned pvtRegionIdx = elemCtx.primaryVars(dofIdx, 0).pvtRegionIndex();
156 if (FluidSystem::phaseIsActive(oilPhaseIdx))
157 SoMax = std::max(getValue(fs.saturation(oilPhaseIdx)),
158 elemCtx.problem().maxOilSaturation(globalDofIdx));
160 if (FluidSystem::phaseIsActive(gasPhaseIdx) && FluidSystem::phaseIsActive(oilPhaseIdx)) {
161 Scalar x_oG = getValue(fs.moleFraction(oilPhaseIdx, gasCompIdx));
162 Scalar x_gO = getValue(fs.moleFraction(gasPhaseIdx, oilCompIdx));
163 Scalar X_oG = getValue(fs.massFraction(oilPhaseIdx, gasCompIdx));
164 Scalar X_gO = getValue(fs.massFraction(gasPhaseIdx, oilCompIdx));
165 Scalar Rs = FluidSystem::convertXoGToRs(X_oG, pvtRegionIdx);
166 Scalar Rv = FluidSystem::convertXgOToRv(X_gO, pvtRegionIdx);
169 FluidSystem::template saturatedDissolutionFactor<FluidState, Scalar>(fs,
173 Scalar X_oG_sat = FluidSystem::convertRsToXoG(RsSat, pvtRegionIdx);
174 Scalar x_oG_sat = FluidSystem::convertXoGToxoG(X_oG_sat, pvtRegionIdx);
177 FluidSystem::template saturatedDissolutionFactor<FluidState, Scalar>(fs,
181 Scalar X_gO_sat = FluidSystem::convertRvToXgO(RvSat, pvtRegionIdx);
182 Scalar x_gO_sat = FluidSystem::convertXgOToxgO(X_gO_sat, pvtRegionIdx);
184 gasDissolutionFactor_[globalDofIdx] = Rs;
187 oilVaporizationFactor_[globalDofIdx] = Rv;
190 oilSaturationPressure_[globalDofIdx] =
191 FluidSystem::template saturationPressure<FluidState, Scalar>(fs, oilPhaseIdx, pvtRegionIdx);
194 gasSaturationPressure_[globalDofIdx] =
195 FluidSystem::template saturationPressure<FluidState, Scalar>(fs, gasPhaseIdx, pvtRegionIdx);
198 saturatedOilGasDissolutionFactor_[globalDofIdx] = RsSat;
201 saturatedGasOilVaporizationFactor_[globalDofIdx] = RvSat;
204 if (x_oG_sat <= 0.0) {
205 oilSaturationRatio_[globalDofIdx] = 1.0;
208 oilSaturationRatio_[globalDofIdx] = x_oG / x_oG_sat;
211 if (x_gO_sat <= 0.0) {
212 gasSaturationRatio_[globalDofIdx] = 1.0;
215 gasSaturationRatio_[globalDofIdx] = x_gO / x_gO_sat;
220 oilFormationVolumeFactor_[globalDofIdx] =
221 1.0 / FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, oilPhaseIdx, pvtRegionIdx);
224 gasFormationVolumeFactor_[globalDofIdx] =
225 1.0 / FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, gasPhaseIdx, pvtRegionIdx);
228 waterFormationVolumeFactor_[globalDofIdx] =
229 1.0 / FluidSystem::template inverseFormationVolumeFactor<FluidState, Scalar>(fs, waterPhaseIdx, pvtRegionIdx);
233 primaryVarsMeaningWater_[globalDofIdx] =
234 static_cast<int>(primaryVars.primaryVarsMeaningWater());
235 primaryVarsMeaningGas_[globalDofIdx] =
236 static_cast<int>(primaryVars.primaryVarsMeaningGas());
237 primaryVarsMeaningPressure_[globalDofIdx] =
238 static_cast<int>(primaryVars.primaryVarsMeaningPressure());
286 this->
commitScalarBuffer_(baseWriter,
"primary vars meaning water", primaryVarsMeaningWater_);
288 this->
commitScalarBuffer_(baseWriter,
"primary vars meaning pressure", primaryVarsMeaningPressure_);
294 ScalarBuffer gasDissolutionFactor_{};
295 ScalarBuffer oilVaporizationFactor_{};
296 ScalarBuffer oilFormationVolumeFactor_{};
297 ScalarBuffer gasFormationVolumeFactor_{};
298 ScalarBuffer waterFormationVolumeFactor_{};
299 ScalarBuffer oilSaturationPressure_{};
300 ScalarBuffer gasSaturationPressure_{};
302 ScalarBuffer saturatedOilGasDissolutionFactor_{};
303 ScalarBuffer saturatedGasOilVaporizationFactor_{};
304 ScalarBuffer oilSaturationRatio_{};
305 ScalarBuffer gasSaturationRatio_{};
307 ScalarBuffer primaryVarsMeaningPressure_{};
308 ScalarBuffer primaryVarsMeaningWater_{};
309 ScalarBuffer primaryVarsMeaningGas_{};
Declares the properties required by the black oil model.
The base class for writer modules.
Definition: baseoutputmodule.hh:67
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:85
void commitScalarBuffer_(BaseOutputWriter &baseWriter, const char *name, ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer containing scalar quantities to the result file.
Definition: baseoutputmodule.hh:244
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:156
The base class for all output writers.
Definition: baseoutputwriter.hh:44
VTK output module for the black oil model's parameters.
Definition: vtkblackoilmodule.hpp:53
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkblackoilmodule.hpp:141
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkblackoilmodule.hpp:88
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkblackoilmodule.hpp:246
VtkBlackOilModule(const Simulator &simulator)
Definition: vtkblackoilmodule.hpp:78
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkblackoilmodule.hpp:97
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:66
Declare the properties used by the infrastructure code of the finite volume discretizations.
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:235
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Struct holding the parameters for VtkBlackoilOutputModule.
Definition: vtkblackoilparams.hpp:53
bool gasSaturationPressureOutput_
Definition: vtkblackoilparams.hpp:66
bool oilFormationVolumeFactorOutput_
Definition: vtkblackoilparams.hpp:62
static void registerParameters()
Registers the parameters in parameter system.
bool gasFormationVolumeFactorOutput_
Definition: vtkblackoilparams.hpp:63
bool oilVaporizationFactorOutput_
Definition: vtkblackoilparams.hpp:61
bool saturatedGasOilVaporizationFactorOutput_
Definition: vtkblackoilparams.hpp:68
bool primaryVarsMeaningOutput_
Definition: vtkblackoilparams.hpp:70
bool gasDissolutionFactorOutput_
Definition: vtkblackoilparams.hpp:60
bool waterFormationVolumeFactorOutput_
Definition: vtkblackoilparams.hpp:64
bool saturatedOilGasDissolutionFactorOutput_
Definition: vtkblackoilparams.hpp:67
bool saturationRatiosOutput_
Definition: vtkblackoilparams.hpp:69
bool oilSaturationPressureOutput_
Definition: vtkblackoilparams.hpp:65
void read()
Reads the parameter values from the parameter system.