27#ifndef EWOMS_VTK_MULTI_PHASE_MODULE_HH
28#define EWOMS_VTK_MULTI_PHASE_MODULE_HH
30#include <dune/common/fvector.hh>
32#include <opm/material/common/MathToolbox.hpp>
33#include <opm/material/common/Valgrind.hpp>
83template<
class TypeTag>
96 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
99 enum { dimWorld = GridView::dimensionworld };
100 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
107 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
121 Parameters::Register<Parameters::VtkWriteExtrusionFactor>
122 (
"Include the extrusion factor of the degrees of freedom into the VTK output files");
123 Parameters::Register<Parameters::VtkWritePressures>
124 (
"Include the phase pressures in the VTK output files");
125 Parameters::Register<Parameters::VtkWriteDensities>
126 (
"Include the phase densities in the VTK output files");
127 Parameters::Register<Parameters::VtkWriteSaturations>
128 (
"Include the phase saturations in the VTK output files");
129 Parameters::Register<Parameters::VtkWriteMobilities>
130 (
"Include the phase mobilities in the VTK output files");
131 Parameters::Register<Parameters::VtkWriteRelativePermeabilities>
132 (
"Include the phase relative permeabilities in the VTK output files");
133 Parameters::Register<Parameters::VtkWriteViscosities>
134 (
"Include component phase viscosities in the VTK output files");
135 Parameters::Register<Parameters::VtkWriteAverageMolarMasses>
136 (
"Include the average phase mass in the VTK output files");
137 Parameters::Register<Parameters::VtkWritePorosity>
138 (
"Include the porosity in the VTK output files");
139 Parameters::Register<Parameters::VtkWriteIntrinsicPermeabilities>
140 (
"Include the intrinsic permeability in the VTK output files");
141 Parameters::Register<Parameters::VtkWriteFilterVelocities>
142 (
"Include in the filter velocities of the phases the VTK output files");
143 Parameters::Register<Parameters::VtkWritePotentialGradients>
144 (
"Include the phase pressure potential gradients in the VTK output files");
165 if (velocityOutput_()) {
166 size_t nDof = this->
simulator_.model().numGridDof();
167 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
168 velocity_[phaseIdx].resize(nDof);
169 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
170 velocity_[phaseIdx][dofIdx].resize(dimWorld);
171 velocity_[phaseIdx][dofIdx] = 0.0;
177 if (potentialGradientOutput_()) {
178 size_t nDof = this->
simulator_.model().numGridDof();
179 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
180 potentialGradient_[phaseIdx].resize(nDof);
181 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
182 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
183 potentialGradient_[phaseIdx][dofIdx] = 0.0;
197 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
201 const auto& problem = elemCtx.problem();
202 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
203 unsigned I = elemCtx.globalSpaceIndex(i, 0);
204 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
205 const auto& fs = intQuants.fluidState();
207 if (extrusionFactorOutput_()) extrusionFactor_[I] = intQuants.extrusionFactor();
208 if (porosityOutput_()) porosity_[I] = getValue(intQuants.porosity());
210 if (intrinsicPermeabilityOutput_()) {
211 const auto& K = problem.intrinsicPermeability(elemCtx, i, 0);
212 for (
unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx)
213 for (
unsigned colIdx = 0; colIdx < K.cols; ++colIdx)
214 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
217 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
218 if (!FluidSystem::phaseIsActive(phaseIdx)) {
221 if (pressureOutput_())
222 pressure_[phaseIdx][I] = getValue(fs.pressure(phaseIdx));
223 if (densityOutput_())
224 density_[phaseIdx][I] = getValue(fs.density(phaseIdx));
225 if (saturationOutput_())
226 saturation_[phaseIdx][I] = getValue(fs.saturation(phaseIdx));
227 if (mobilityOutput_())
228 mobility_[phaseIdx][I] = getValue(intQuants.mobility(phaseIdx));
229 if (relativePermeabilityOutput_())
230 relativePermeability_[phaseIdx][I] = getValue(intQuants.relativePermeability(phaseIdx));
231 if (viscosityOutput_())
232 viscosity_[phaseIdx][I] = getValue(fs.viscosity(phaseIdx));
233 if (averageMolarMassOutput_())
234 averageMolarMass_[phaseIdx][I] = getValue(fs.averageMolarMass(phaseIdx));
238 if (potentialGradientOutput_()) {
240 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
241 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
243 unsigned i = extQuants.interiorIndex();
244 unsigned I = elemCtx.globalSpaceIndex(i, 0);
246 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
247 Scalar weight = extQuants.extrusionFactor();
249 potentialWeight_[phaseIdx][I] += weight;
251 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
253 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx)
254 pGrad[dimIdx] = getValue(inputPGrad[dimIdx])*weight;
255 potentialGradient_[phaseIdx][I] += pGrad;
260 if (velocityOutput_()) {
262 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
263 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
265 unsigned i = extQuants.interiorIndex();
266 unsigned I = elemCtx.globalSpaceIndex(i, 0);
268 unsigned j = extQuants.exteriorIndex();
269 unsigned J = elemCtx.globalSpaceIndex(j, 0);
271 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
272 Scalar weight = std::max<Scalar>(1e-16,
273 std::abs(getValue(extQuants.volumeFlux(phaseIdx))));
274 Valgrind::CheckDefined(extQuants.extrusionFactor());
275 assert(extQuants.extrusionFactor() > 0);
276 weight *= extQuants.extrusionFactor();
278 const auto& inputV = extQuants.filterVelocity(phaseIdx);
280 for (
unsigned k = 0; k < dimWorld; ++k)
281 v[k] = getValue(inputV[k]);
282 if (v.two_norm() > 1e-20)
283 weight /= v.two_norm();
286 velocity_[phaseIdx][I] += v;
287 velocity_[phaseIdx][J] += v;
289 velocityWeight_[phaseIdx][I] += weight;
290 velocityWeight_[phaseIdx][J] += weight;
305 if (extrusionFactorOutput_())
307 if (pressureOutput_())
309 if (densityOutput_())
311 if (saturationOutput_())
313 if (mobilityOutput_())
315 if (relativePermeabilityOutput_())
317 if (viscosityOutput_())
319 if (averageMolarMassOutput_())
322 if (porosityOutput_())
324 if (intrinsicPermeabilityOutput_())
327 if (velocityOutput_()) {
328 size_t numDof = this->
simulator_.model().numGridDof();
330 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
333 for (
unsigned i = 0; i < numDof; ++i)
334 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
337 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
339 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
343 if (potentialGradientOutput_()) {
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 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
353 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx).data());
355 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
356 potentialGradient_[phaseIdx],
372 return velocityOutput_() || potentialGradientOutput_();
376 static bool extrusionFactorOutput_()
378 static bool val = Parameters::Get<Parameters::VtkWriteExtrusionFactor>();
382 static bool pressureOutput_()
384 static bool val = Parameters::Get<Parameters::VtkWritePressures>();
388 static bool densityOutput_()
390 static bool val = Parameters::Get<Parameters::VtkWriteDensities>();
394 static bool saturationOutput_()
396 static bool val = Parameters::Get<Parameters::VtkWriteSaturations>();
400 static bool mobilityOutput_()
402 static bool val = Parameters::Get<Parameters::VtkWriteMobilities>();
406 static bool relativePermeabilityOutput_()
408 static bool val = Parameters::Get<Parameters::VtkWriteRelativePermeabilities>();
412 static bool viscosityOutput_()
414 static bool val = Parameters::Get<Parameters::VtkWriteViscosities>();
418 static bool averageMolarMassOutput_()
420 static bool val = Parameters::Get<Parameters::VtkWriteAverageMolarMasses>();
424 static bool porosityOutput_()
426 static bool val = Parameters::Get<Parameters::VtkWritePorosity>();
430 static bool intrinsicPermeabilityOutput_()
432 static bool val = Parameters::Get<Parameters::VtkWriteIntrinsicPermeabilities>();
436 static bool velocityOutput_()
438 static bool val = Parameters::Get<Parameters::VtkWriteFilterVelocities>();
442 static bool potentialGradientOutput_()
444 static bool val = Parameters::Get<Parameters::VtkWritePotentialGradients>();
448 ScalarBuffer extrusionFactor_;
449 PhaseBuffer pressure_;
450 PhaseBuffer density_;
451 PhaseBuffer saturation_;
452 PhaseBuffer mobility_;
453 PhaseBuffer relativePermeability_;
454 PhaseBuffer viscosity_;
455 PhaseBuffer averageMolarMass_;
457 ScalarBuffer porosity_;
458 TensorBuffer intrinsicPermeability_;
460 PhaseVectorBuffer velocity_;
461 PhaseBuffer velocityWeight_;
463 PhaseVectorBuffer potentialGradient_;
464 PhaseBuffer potentialWeight_;
The base class for writer modules.
Definition: baseoutputmodule.hh:67
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:85
BaseOutputWriter::TensorBuffer TensorBuffer
Definition: baseoutputmodule.hh:87
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:345
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:201
const Simulator & simulator_
Definition: baseoutputmodule.hh:434
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:94
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:288
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:90
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: baseoutputmodule.hh:86
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition: baseoutputmodule.hh:166
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 quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hh:85
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hh:195
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hh:119
VtkMultiPhaseModule(const Simulator &simulator)
Definition: vtkmultiphasemodule.hh:112
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hh:151
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:370
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hh:299
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:66
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilnewtonmethodparameters.hh:31
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.
Definition: vtkmultiphasemodule.hh:55
static constexpr bool value
Definition: vtkmultiphasemodule.hh:55
Definition: vtkmultiphasemodule.hh:50
static constexpr bool value
Definition: vtkmultiphasemodule.hh:50
Definition: vtkmultiphasemodule.hh:48
static constexpr bool value
Definition: vtkmultiphasemodule.hh:48
Definition: vtkmultiphasemodule.hh:59
static constexpr bool value
Definition: vtkmultiphasemodule.hh:59
Definition: vtkmultiphasemodule.hh:57
static constexpr bool value
Definition: vtkmultiphasemodule.hh:57
Definition: vtkmultiphasemodule.hh:52
static constexpr bool value
Definition: vtkmultiphasemodule.hh:52
Definition: vtkmultiphasemodule.hh:56
static constexpr bool value
Definition: vtkmultiphasemodule.hh:56
Definition: vtkmultiphasemodule.hh:58
static constexpr bool value
Definition: vtkmultiphasemodule.hh:58
Definition: vtkmultiphasemodule.hh:49
static constexpr bool value
Definition: vtkmultiphasemodule.hh:49
Definition: vtkmultiphasemodule.hh:53
static constexpr bool value
Definition: vtkmultiphasemodule.hh:53
Definition: vtkmultiphasemodule.hh:51
static constexpr bool value
Definition: vtkmultiphasemodule.hh:51
Definition: vtkmultiphasemodule.hh:54
static constexpr bool value
Definition: vtkmultiphasemodule.hh:54