27#ifndef OPM_VTK_MULTI_PHASE_MODULE_HPP
28#define OPM_VTK_MULTI_PHASE_MODULE_HPP
30#include <dune/common/fvector.hh>
32#include <opm/material/common/MathToolbox.hpp>
33#include <opm/material/common/Valgrind.hpp>
71template<
class TypeTag>
84 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
87 enum { dimWorld = GridView::dimensionworld };
88 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
96 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
154 const std::size_t nDof = this->
simulator_.model().numGridDof();
155 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
156 velocity_[phaseIdx].resize(nDof);
157 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
158 velocity_[phaseIdx][dofIdx].resize(dimWorld);
159 velocity_[phaseIdx][dofIdx] = 0.0;
166 const std::size_t nDof = this->
simulator_.model().numGridDof();
167 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
168 potentialGradient_[phaseIdx].resize(nDof);
169 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
170 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
171 potentialGradient_[phaseIdx][dofIdx] = 0.0;
185 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
189 const auto& problem = elemCtx.problem();
190 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
191 const unsigned I = elemCtx.globalSpaceIndex(i, 0);
192 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
193 const auto& fs = intQuants.fluidState();
196 extrusionFactor_[I] = intQuants.extrusionFactor();
199 porosity_[I] = getValue(intQuants.porosity());
203 const auto& K = problem.intrinsicPermeability(elemCtx, i, 0);
204 for (
unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx) {
205 for (
unsigned colIdx = 0; colIdx < K.cols; ++colIdx) {
206 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
211 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
212 if (!FluidSystem::phaseIsActive(phaseIdx)) {
216 pressure_[phaseIdx][I] = getValue(fs.pressure(phaseIdx));
219 density_[phaseIdx][I] = getValue(fs.density(phaseIdx));
222 saturation_[phaseIdx][I] = getValue(fs.saturation(phaseIdx));
225 mobility_[phaseIdx][I] = getValue(intQuants.mobility(phaseIdx));
228 relativePermeability_[phaseIdx][I] = getValue(intQuants.relativePermeability(phaseIdx));
231 viscosity_[phaseIdx][I] = getValue(fs.viscosity(phaseIdx));
234 averageMolarMass_[phaseIdx][I] = getValue(fs.averageMolarMass(phaseIdx));
241 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++faceIdx) {
242 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
244 const unsigned i = extQuants.interiorIndex();
245 const unsigned I = elemCtx.globalSpaceIndex(i, 0);
247 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
248 const Scalar weight = extQuants.extrusionFactor();
250 potentialWeight_[phaseIdx][I] += weight;
252 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
254 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
255 pGrad[dimIdx] = getValue(inputPGrad[dimIdx]) * weight;
257 potentialGradient_[phaseIdx][I] += pGrad;
264 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++faceIdx) {
265 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
267 const unsigned i = extQuants.interiorIndex();
268 const unsigned I = elemCtx.globalSpaceIndex(i, 0);
270 const unsigned j = extQuants.exteriorIndex();
271 const unsigned J = elemCtx.globalSpaceIndex(j, 0);
273 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
274 Scalar weight = std::max(Scalar{1e-16},
275 std::abs(getValue(extQuants.volumeFlux(phaseIdx))));
276 Valgrind::CheckDefined(extQuants.extrusionFactor());
277 assert(extQuants.extrusionFactor() > 0);
278 weight *= extQuants.extrusionFactor();
280 const auto& inputV = extQuants.filterVelocity(phaseIdx);
282 for (
unsigned k = 0; k < dimWorld; ++k) {
283 v[k] = getValue(inputV[k]);
285 if (v.two_norm() > 1e-20) {
286 weight /= v.two_norm();
290 velocity_[phaseIdx][I] += v;
291 velocity_[phaseIdx][J] += v;
293 velocityWeight_[phaseIdx][I] += weight;
294 velocityWeight_[phaseIdx][J] += weight;
311 extrusionFactor_, BufferType::Dof);
327 relativePermeability_, BufferType::Dof);
334 averageMolarMass_, BufferType::Dof);
342 intrinsicPermeability_, BufferType::Dof);
346 const std::size_t numDof = this->
simulator_.model().numGridDof();
347 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
350 for (
unsigned i = 0; i < numDof; ++i) {
351 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
355 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
357 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
362 const std::size_t numDof = this->
simulator_.model().numGridDof();
363 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
366 for (
unsigned i = 0; i < numDof; ++i) {
367 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
371 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx).data());
373 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
374 potentialGradient_[phaseIdx],
395 ScalarBuffer extrusionFactor_{};
396 PhaseBuffer pressure_{};
397 PhaseBuffer density_{};
398 PhaseBuffer saturation_{};
399 PhaseBuffer mobility_{};
400 PhaseBuffer relativePermeability_{};
401 PhaseBuffer viscosity_{};
402 PhaseBuffer averageMolarMass_{};
404 ScalarBuffer porosity_{};
405 TensorBuffer intrinsicPermeability_{};
407 PhaseVectorBuffer velocity_{};
408 PhaseBuffer velocityWeight_{};
410 PhaseVectorBuffer potentialGradient_{};
411 PhaseBuffer potentialWeight_{};
The base class for writer modules.
Definition: baseoutputmodule.hh:68
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:86
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:157
BaseOutputWriter::TensorBuffer TensorBuffer
Definition: baseoutputmodule.hh:88
const Simulator & simulator_
Definition: baseoutputmodule.hh:426
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:95
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:91
BufferType
Definition: baseoutputmodule.hh:143
void commitTensorBuffer_(BaseOutputWriter &baseWriter, const char *name, TensorBuffer &buffer, BufferType bufferType)
Add a buffer containing tensorial quantities to the result file.
Definition: baseoutputmodule.hh:282
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: baseoutputmodule.hh:87
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a tensorial quantity.
Definition: baseoutputmodule.hh:166
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:337
void commitScalarBuffer_(BaseOutputWriter &baseWriter, const char *name, ScalarBuffer &buffer, BufferType bufferType)
Add a buffer containing scalar quantities to the result file.
Definition: baseoutputmodule.hh:238
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:198
The base class for all output writers.
Definition: baseoutputwriter.hh:46
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hpp:73
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hpp:110
VtkMultiPhaseModule(const Simulator &simulator)
Definition: vtkmultiphasemodule.hpp:101
bool needExtensiveQuantities() const override
Returns true iff the module needs to access the extensive quantities of a context to do its job.
Definition: vtkmultiphasemodule.hpp:388
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hpp:303
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hpp:119
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hpp:183
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:65
Declare the properties used by the infrastructure code of the finite volume discretizations.
Definition: blackoilboundaryratevector.hh:39
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:233
This file provides the infrastructure to retrieve run-time parameters.
The Opm property system, traits with inheritance.
Struct holding the parameters for VtkMultiPhaseModule.
Definition: vtkmultiphaseparams.hpp:54
bool averageMolarMassOutput_
Definition: vtkmultiphaseparams.hpp:68
bool pressureOutput_
Definition: vtkmultiphaseparams.hpp:62
bool extrusionFactorOutput_
Definition: vtkmultiphaseparams.hpp:61
bool porosityOutput_
Definition: vtkmultiphaseparams.hpp:69
bool velocityOutput_
Definition: vtkmultiphaseparams.hpp:71
bool densityOutput_
Definition: vtkmultiphaseparams.hpp:63
static void registerParameters()
Registers the parameters in parameter system.
bool saturationOutput_
Definition: vtkmultiphaseparams.hpp:64
bool potentialGradientOutput_
Definition: vtkmultiphaseparams.hpp:72
bool intrinsicPermeabilityOutput_
Definition: vtkmultiphaseparams.hpp:70
void read()
Reads the parameter values from the parameter system.
bool viscosityOutput_
Definition: vtkmultiphaseparams.hpp:67
bool relativePermeabilityOutput_
Definition: vtkmultiphaseparams.hpp:66
bool mobilityOutput_
Definition: vtkmultiphaseparams.hpp:65