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/material/common/Valgrind.hpp>
39#include <dune/common/fvector.hh>
54template<
class TypeTag,
class MyTypeTag>
56template<
class TypeTag,
class MyTypeTag>
58template<
class TypeTag,
class MyTypeTag>
60template<
class TypeTag,
class MyTypeTag>
62template<
class TypeTag,
class MyTypeTag>
64template<
class TypeTag,
class MyTypeTag>
66template<
class TypeTag,
class MyTypeTag>
68template<
class TypeTag,
class MyTypeTag>
70template<
class TypeTag,
class MyTypeTag>
72template<
class TypeTag,
class MyTypeTag>
74template<
class TypeTag,
class MyTypeTag>
76template<
class TypeTag,
class MyTypeTag>
80template<
class TypeTag>
82template<
class TypeTag>
83struct VtkWritePressures<TypeTag, TTag::VtkMultiPhase> {
static constexpr bool value =
true; };
84template<
class TypeTag>
85struct VtkWriteDensities<TypeTag, TTag::VtkMultiPhase> {
static constexpr bool value =
true; };
86template<
class TypeTag>
88template<
class TypeTag>
90template<
class TypeTag>
92template<
class TypeTag>
94template<
class TypeTag>
96template<
class TypeTag>
97struct VtkWritePorosity<TypeTag, TTag::VtkMultiPhase> {
static constexpr bool value =
true; };
98template<
class TypeTag>
100template<
class TypeTag>
102template<
class TypeTag>
127template<
class TypeTag>
140 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
143 enum { dimWorld = GridView::dimensionworld };
144 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
151 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
165 Parameters::registerParam<TypeTag, Properties::VtkWriteExtrusionFactor>
166 (
"Include the extrusion factor of the degrees of freedom into the VTK output files");
167 Parameters::registerParam<TypeTag, Properties::VtkWritePressures>
168 (
"Include the phase pressures in the VTK output files");
169 Parameters::registerParam<TypeTag, Properties::VtkWriteDensities>
170 (
"Include the phase densities in the VTK output files");
171 Parameters::registerParam<TypeTag, Properties::VtkWriteSaturations>
172 (
"Include the phase saturations in the VTK output files");
173 Parameters::registerParam<TypeTag, Properties::VtkWriteMobilities>
174 (
"Include the phase mobilities in the VTK output files");
175 Parameters::registerParam<TypeTag, Properties::VtkWriteRelativePermeabilities>
176 (
"Include the phase relative permeabilities in the VTK output files");
177 Parameters::registerParam<TypeTag, Properties::VtkWriteViscosities>
178 (
"Include component phase viscosities in the VTK output files");
179 Parameters::registerParam<TypeTag, Properties::VtkWriteAverageMolarMasses>
180 (
"Include the average phase mass in the VTK output files");
181 Parameters::registerParam<TypeTag, Properties::VtkWritePorosity>
182 (
"Include the porosity in the VTK output files");
183 Parameters::registerParam<TypeTag, Properties::VtkWriteIntrinsicPermeabilities>
184 (
"Include the intrinsic permeability in the VTK output files");
185 Parameters::registerParam<TypeTag, Properties::VtkWriteFilterVelocities>
186 (
"Include in the filter velocities of the phases the VTK output files");
187 Parameters::registerParam<TypeTag, Properties::VtkWritePotentialGradients>
188 (
"Include the phase pressure potential gradients in the VTK output files");
209 if (velocityOutput_()) {
210 size_t nDof = this->
simulator_.model().numGridDof();
211 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
212 velocity_[phaseIdx].resize(nDof);
213 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
214 velocity_[phaseIdx][dofIdx].resize(dimWorld);
215 velocity_[phaseIdx][dofIdx] = 0.0;
221 if (potentialGradientOutput_()) {
222 size_t nDof = this->
simulator_.model().numGridDof();
223 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
224 potentialGradient_[phaseIdx].resize(nDof);
225 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
226 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
227 potentialGradient_[phaseIdx][dofIdx] = 0.0;
241 if (!Parameters::get<TypeTag, Properties::EnableVtkOutput>())
244 const auto& problem = elemCtx.problem();
245 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
246 unsigned I = elemCtx.globalSpaceIndex(i, 0);
247 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
248 const auto& fs = intQuants.fluidState();
250 if (extrusionFactorOutput_()) extrusionFactor_[I] = intQuants.extrusionFactor();
251 if (porosityOutput_()) porosity_[I] = getValue(intQuants.porosity());
253 if (intrinsicPermeabilityOutput_()) {
254 const auto& K = problem.intrinsicPermeability(elemCtx, i, 0);
255 for (
unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx)
256 for (
unsigned colIdx = 0; colIdx < K.cols; ++colIdx)
257 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
260 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
261 if (!FluidSystem::phaseIsActive(phaseIdx)) {
264 if (pressureOutput_())
265 pressure_[phaseIdx][I] = getValue(fs.pressure(phaseIdx));
266 if (densityOutput_())
267 density_[phaseIdx][I] = getValue(fs.density(phaseIdx));
268 if (saturationOutput_())
269 saturation_[phaseIdx][I] = getValue(fs.saturation(phaseIdx));
270 if (mobilityOutput_())
271 mobility_[phaseIdx][I] = getValue(intQuants.mobility(phaseIdx));
272 if (relativePermeabilityOutput_())
273 relativePermeability_[phaseIdx][I] = getValue(intQuants.relativePermeability(phaseIdx));
274 if (viscosityOutput_())
275 viscosity_[phaseIdx][I] = getValue(fs.viscosity(phaseIdx));
276 if (averageMolarMassOutput_())
277 averageMolarMass_[phaseIdx][I] = getValue(fs.averageMolarMass(phaseIdx));
281 if (potentialGradientOutput_()) {
283 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
284 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
286 unsigned i = extQuants.interiorIndex();
287 unsigned I = elemCtx.globalSpaceIndex(i, 0);
289 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
290 Scalar weight = extQuants.extrusionFactor();
292 potentialWeight_[phaseIdx][I] += weight;
294 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
296 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
297 pGrad[dimIdx] = getValue(inputPGrad[dimIdx])*weight;
298 potentialGradient_[phaseIdx][I] += pGrad;
303 if (velocityOutput_()) {
305 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
306 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
308 unsigned i = extQuants.interiorIndex();
309 unsigned I = elemCtx.globalSpaceIndex(i, 0);
311 unsigned j = extQuants.exteriorIndex();
312 unsigned J = elemCtx.globalSpaceIndex(j, 0);
314 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
315 Scalar weight = std::max<Scalar>(1e-16,
316 std::abs(getValue(extQuants.volumeFlux(phaseIdx))));
317 Valgrind::CheckDefined(extQuants.extrusionFactor());
318 assert(extQuants.extrusionFactor() > 0);
319 weight *= extQuants.extrusionFactor();
321 const auto& inputV = extQuants.filterVelocity(phaseIdx);
323 for (
unsigned k = 0; k < dimWorld; ++k)
324 v[k] = getValue(inputV[k]);
325 if (v.two_norm() > 1e-20)
326 weight /= v.two_norm();
329 velocity_[phaseIdx][I] += v;
330 velocity_[phaseIdx][J] += v;
332 velocityWeight_[phaseIdx][I] += weight;
333 velocityWeight_[phaseIdx][J] += weight;
348 if (extrusionFactorOutput_())
350 if (pressureOutput_())
352 if (densityOutput_())
354 if (saturationOutput_())
356 if (mobilityOutput_())
358 if (relativePermeabilityOutput_())
360 if (viscosityOutput_())
362 if (averageMolarMassOutput_())
365 if (porosityOutput_())
367 if (intrinsicPermeabilityOutput_())
370 if (velocityOutput_()) {
371 size_t numDof = this->
simulator_.model().numGridDof();
373 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
376 for (
unsigned i = 0; i < numDof; ++i)
377 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
380 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
382 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
386 if (potentialGradientOutput_()) {
387 size_t numDof = this->
simulator_.model().numGridDof();
389 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
392 for (
unsigned i = 0; i < numDof; ++i)
393 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
396 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx).data());
398 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
399 potentialGradient_[phaseIdx],
415 return velocityOutput_() || potentialGradientOutput_();
419 static bool extrusionFactorOutput_()
421 static bool val = Parameters::get<TypeTag, Properties::VtkWriteExtrusionFactor>();
425 static bool pressureOutput_()
427 static bool val = Parameters::get<TypeTag, Properties::VtkWritePressures>();
431 static bool densityOutput_()
433 static bool val = Parameters::get<TypeTag, Properties::VtkWriteDensities>();
437 static bool saturationOutput_()
439 static bool val = Parameters::get<TypeTag, Properties::VtkWriteSaturations>();
443 static bool mobilityOutput_()
445 static bool val = Parameters::get<TypeTag, Properties::VtkWriteMobilities>();
449 static bool relativePermeabilityOutput_()
451 static bool val = Parameters::get<TypeTag, Properties::VtkWriteRelativePermeabilities>();
455 static bool viscosityOutput_()
457 static bool val = Parameters::get<TypeTag, Properties::VtkWriteViscosities>();
461 static bool averageMolarMassOutput_()
463 static bool val = Parameters::get<TypeTag, Properties::VtkWriteAverageMolarMasses>();
467 static bool porosityOutput_()
469 static bool val = Parameters::get<TypeTag, Properties::VtkWritePorosity>();
473 static bool intrinsicPermeabilityOutput_()
475 static bool val = Parameters::get<TypeTag, Properties::VtkWriteIntrinsicPermeabilities>();
479 static bool velocityOutput_()
481 static bool val = Parameters::get<TypeTag, Properties::VtkWriteFilterVelocities>();
485 static bool potentialGradientOutput_()
487 static bool val = Parameters::get<TypeTag, Properties::VtkWritePotentialGradients>();
491 ScalarBuffer extrusionFactor_;
492 PhaseBuffer pressure_;
493 PhaseBuffer density_;
494 PhaseBuffer saturation_;
495 PhaseBuffer mobility_;
496 PhaseBuffer relativePermeability_;
497 PhaseBuffer viscosity_;
498 PhaseBuffer averageMolarMass_;
500 ScalarBuffer porosity_;
501 TensorBuffer intrinsicPermeability_;
503 PhaseVectorBuffer velocity_;
504 PhaseBuffer velocityWeight_;
506 PhaseVectorBuffer potentialGradient_;
507 PhaseBuffer potentialWeight_;
The base class for writer modules.
Definition: baseoutputmodule.hh:73
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:91
BaseOutputWriter::TensorBuffer TensorBuffer
Definition: baseoutputmodule.hh:93
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:419
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:246
const Simulator & simulator_
Definition: baseoutputmodule.hh:519
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:100
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:352
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:96
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: baseoutputmodule.hh:92
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition: baseoutputmodule.hh:182
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:316
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:162
The base class for all output writers.
Definition: baseoutputwriter.hh:44
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hh:129
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hh:239
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hh:163
VtkMultiPhaseModule(const Simulator &simulator)
Definition: vtkmultiphasemodule.hh:156
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hh:195
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:413
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hh:342
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: vtkmultiphasemodule.hh:49
a tag to mark properties as undefined
Definition: propertysystem.hh:40
Definition: vtkmultiphasemodule.hh:69
Definition: vtkmultiphasemodule.hh:59
Definition: vtkmultiphasemodule.hh:55
Definition: vtkmultiphasemodule.hh:77
Definition: vtkmultiphasemodule.hh:73
Definition: vtkmultiphasemodule.hh:63
Definition: vtkmultiphasemodule.hh:71
Definition: vtkmultiphasemodule.hh:75
Definition: vtkmultiphasemodule.hh:57
Definition: vtkmultiphasemodule.hh:65
Definition: vtkmultiphasemodule.hh:61
Definition: vtkmultiphasemodule.hh:67