27 #ifndef EWOMS_VTK_MULTI_PHASE_MODULE_HH 28 #define EWOMS_VTK_MULTI_PHASE_MODULE_HH 36 #include <opm/material/common/MathToolbox.hpp> 37 #include <opm/common/Valgrind.hpp> 39 #include <dune/common/fvector.hh> 44 namespace Properties {
71 SET_BOOL_PROP(VtkMultiPhase, VtkWriteRelativePermeabilities,
true);
73 SET_BOOL_PROP(VtkMultiPhase, VtkWriteAverageMolarMasses,
false);
75 SET_BOOL_PROP(VtkMultiPhase, VtkWriteIntrinsicPermeabilities,
false);
76 SET_BOOL_PROP(VtkMultiPhase, VtkWritePotentialGradients,
false);
77 SET_BOOL_PROP(VtkMultiPhase, VtkWriteFilterVelocities,
false);
98 template<
class TypeTag>
104 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
105 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
106 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
108 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
109 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
110 typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
112 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
115 typedef Opm::MathToolbox<Evaluation> Toolbox;
117 enum { dimWorld = GridView::dimensionworld };
120 typedef typename ParentType::ScalarBuffer ScalarBuffer;
121 typedef typename ParentType::VectorBuffer VectorBuffer;
122 typedef typename ParentType::TensorBuffer TensorBuffer;
123 typedef typename ParentType::PhaseBuffer PhaseBuffer;
125 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
127 typedef std::array<VectorBuffer, numPhases> PhaseVectorBuffer;
140 "Include the extrusion factor of the degrees of freedom into the VTK output files");
142 "Include the phase pressures in the VTK output files");
144 "Include the phase densities in the VTK output files");
146 "Include the phase saturations in the VTK output files");
148 "Include the phase mobilities in the VTK output files");
150 "Include the phase relative permeabilities in the VTK output files");
152 "Include component phase viscosities in the VTK output files");
154 "Include the average phase mass in the VTK output files");
156 "Include the porosity in the VTK output files");
158 "Include the intrinsic permeability in the VTK output files");
160 "Include in the filter velocities of the phases the VTK output files");
162 "Include the phase pressure potential gradients in the VTK output files");
183 if (velocityOutput_()) {
184 size_t nDof = this->simulator_.
model().numGridDof();
185 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
186 velocity_[phaseIdx].resize(nDof);
187 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
188 velocity_[phaseIdx][dofIdx].resize(dimWorld);
189 velocity_[phaseIdx][dofIdx] = 0.0;
195 if (potentialGradientOutput_()) {
196 size_t nDof = this->simulator_.
model().numGridDof();
197 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
198 potentialGradient_[phaseIdx].resize(nDof);
199 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
200 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
201 potentialGradient_[phaseIdx][dofIdx] = 0.0;
215 typedef Opm::MathToolbox<Evaluation> Toolbox;
220 const auto& problem = elemCtx.problem();
221 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
222 unsigned I = elemCtx.globalSpaceIndex(i, 0);
223 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
224 const auto& fs = intQuants.fluidState();
226 if (extrusionFactorOutput_()) extrusionFactor_[I] = intQuants.extrusionFactor();
227 if (porosityOutput_()) porosity_[I] = Toolbox::value(intQuants.porosity());
228 if (intrinsicPermeabilityOutput_()) {
229 const auto& K = problem.intrinsicPermeability(elemCtx, i, 0);
230 intrinsicPermeability_[I].resize(K.rows, K.cols);
231 for (
unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx)
232 for (
unsigned colIdx = 0; colIdx < K.cols; ++colIdx)
233 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
236 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
237 if (pressureOutput_())
238 pressure_[phaseIdx][I] = Toolbox::value(fs.pressure(phaseIdx));
239 if (densityOutput_())
240 density_[phaseIdx][I] = Toolbox::value(fs.density(phaseIdx));
241 if (saturationOutput_())
242 saturation_[phaseIdx][I] = Toolbox::value(fs.saturation(phaseIdx));
243 if (mobilityOutput_())
244 mobility_[phaseIdx][I] = Toolbox::value(intQuants.mobility(phaseIdx));
245 if (relativePermeabilityOutput_())
246 relativePermeability_[phaseIdx][I] = Toolbox::value(intQuants.relativePermeability(phaseIdx));
247 if (viscosityOutput_())
248 viscosity_[phaseIdx][I] = Toolbox::value(fs.viscosity(phaseIdx));
249 if (averageMolarMassOutput_())
250 averageMolarMass_[phaseIdx][I] = Toolbox::value(fs.averageMolarMass(phaseIdx));
254 if (potentialGradientOutput_()) {
256 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
257 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
259 unsigned i = extQuants.interiorIndex();
260 unsigned I = elemCtx.globalSpaceIndex(i, 0);
262 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
263 Scalar weight = extQuants.extrusionFactor();
265 potentialWeight_[phaseIdx][I] += weight;
267 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
269 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
270 pGrad[dimIdx] = Toolbox::value(inputPGrad[dimIdx])*weight;
271 potentialGradient_[phaseIdx][I] += pGrad;
276 if (velocityOutput_()) {
278 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
279 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
281 unsigned i = extQuants.interiorIndex();
282 unsigned I = elemCtx.globalSpaceIndex(i, 0);
284 unsigned j = extQuants.exteriorIndex();
285 unsigned J = elemCtx.globalSpaceIndex(j, 0);
287 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
288 Scalar weight = std::max<Scalar>(1e-16,
289 std::abs(Toolbox::value(extQuants.volumeFlux(phaseIdx))));
290 Opm::Valgrind::CheckDefined(extQuants.extrusionFactor());
291 assert(extQuants.extrusionFactor() > 0);
292 weight *= extQuants.extrusionFactor();
294 const auto& inputV = extQuants.filterVelocity(phaseIdx);
296 for (
unsigned k = 0; k < dimWorld; ++k)
297 v[k] = Toolbox::value(inputV[k]);
298 if (v.two_norm() > 1e-20)
299 weight /= v.two_norm();
302 velocity_[phaseIdx][I] += v;
303 velocity_[phaseIdx][J] += v;
305 velocityWeight_[phaseIdx][I] += weight;
306 velocityWeight_[phaseIdx][J] += weight;
321 if (extrusionFactorOutput_())
323 if (pressureOutput_())
325 if (densityOutput_())
327 if (saturationOutput_())
329 if (mobilityOutput_())
331 if (relativePermeabilityOutput_())
333 if (viscosityOutput_())
335 if (averageMolarMassOutput_())
338 if (porosityOutput_())
340 if (intrinsicPermeabilityOutput_())
343 if (velocityOutput_()) {
344 size_t numDof = this->simulator_.
model().numGridDof();
346 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
349 for (
unsigned i = 0; i < numDof; ++i)
350 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
353 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx));
355 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
359 if (potentialGradientOutput_()) {
360 size_t numDof = this->simulator_.
model().numGridDof();
362 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
365 for (
unsigned i = 0; i < numDof; ++i)
366 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
369 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx));
371 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
372 potentialGradient_[phaseIdx],
388 return velocityOutput_() || potentialGradientOutput_();
392 static bool extrusionFactorOutput_()
395 static bool pressureOutput_()
398 static bool densityOutput_()
401 static bool saturationOutput_()
404 static bool mobilityOutput_()
407 static bool relativePermeabilityOutput_()
408 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteRelativePermeabilities); }
410 static bool viscosityOutput_()
413 static bool averageMolarMassOutput_()
416 static bool porosityOutput_()
419 static bool intrinsicPermeabilityOutput_()
420 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteIntrinsicPermeabilities); }
422 static bool velocityOutput_()
425 static bool potentialGradientOutput_()
428 ScalarBuffer extrusionFactor_;
429 PhaseBuffer pressure_;
430 PhaseBuffer density_;
431 PhaseBuffer saturation_;
432 PhaseBuffer mobility_;
433 PhaseBuffer relativePermeability_;
434 PhaseBuffer viscosity_;
435 PhaseBuffer averageMolarMass_;
437 ScalarBuffer porosity_;
438 TensorBuffer intrinsicPermeability_;
440 PhaseVectorBuffer velocity_;
441 PhaseBuffer velocityWeight_;
443 PhaseVectorBuffer potentialGradient_;
444 PhaseBuffer potentialWeight_;
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hh:213
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
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:189
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:408
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hh:99
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:63
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:305
#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 commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hh:315
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition: baseoutputmodule.hh:190
void commitTensorBuffer_(BaseOutputWriter &baseWriter, const char *name, TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer containing tensorial quantities to the result file.
Definition: baseoutputmodule.hh:341
#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.
The base class for writer modules.
Definition: baseoutputmodule.hh:80
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hh:169
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:99
virtual bool needExtensiveQuantities() const final
Returns true iff the module needs to access the extensive quantities of a context to do its job...
Definition: vtkmultiphasemodule.hh:386
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:235
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hh:137
Provides the magic behind the eWoms property system.
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:170
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