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>
66template<
class TypeTag>
79 static const int vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
82 enum { dimWorld = GridView::dimensionworld };
83 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
90 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
148 size_t nDof = this->
simulator_.model().numGridDof();
149 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
150 velocity_[phaseIdx].resize(nDof);
151 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
152 velocity_[phaseIdx][dofIdx].resize(dimWorld);
153 velocity_[phaseIdx][dofIdx] = 0.0;
160 size_t nDof = this->
simulator_.model().numGridDof();
161 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
162 potentialGradient_[phaseIdx].resize(nDof);
163 for (
unsigned dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
164 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
165 potentialGradient_[phaseIdx][dofIdx] = 0.0;
179 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
183 const auto& problem = elemCtx.problem();
184 for (
unsigned i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
185 unsigned I = elemCtx.globalSpaceIndex(i, 0);
186 const auto& intQuants = elemCtx.intensiveQuantities(i, 0);
187 const auto& fs = intQuants.fluidState();
190 extrusionFactor_[I] = intQuants.extrusionFactor();
193 porosity_[I] = getValue(intQuants.porosity());
197 const auto& K = problem.intrinsicPermeability(elemCtx, i, 0);
198 for (
unsigned rowIdx = 0; rowIdx < K.rows; ++rowIdx) {
199 for (
unsigned colIdx = 0; colIdx < K.cols; ++colIdx) {
200 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
205 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
206 if (!FluidSystem::phaseIsActive(phaseIdx)) {
210 pressure_[phaseIdx][I] = getValue(fs.pressure(phaseIdx));
213 density_[phaseIdx][I] = getValue(fs.density(phaseIdx));
216 saturation_[phaseIdx][I] = getValue(fs.saturation(phaseIdx));
219 mobility_[phaseIdx][I] = getValue(intQuants.mobility(phaseIdx));
222 relativePermeability_[phaseIdx][I] = getValue(intQuants.relativePermeability(phaseIdx));
225 viscosity_[phaseIdx][I] = getValue(fs.viscosity(phaseIdx));
228 averageMolarMass_[phaseIdx][I] = getValue(fs.averageMolarMass(phaseIdx));
235 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
236 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
238 unsigned i = extQuants.interiorIndex();
239 unsigned I = elemCtx.globalSpaceIndex(i, 0);
241 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
242 Scalar weight = extQuants.extrusionFactor();
244 potentialWeight_[phaseIdx][I] += weight;
246 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
248 for (
unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
249 pGrad[dimIdx] = getValue(inputPGrad[dimIdx])*weight;
251 potentialGradient_[phaseIdx][I] += pGrad;
258 for (
unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
259 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
261 unsigned i = extQuants.interiorIndex();
262 unsigned I = elemCtx.globalSpaceIndex(i, 0);
264 unsigned j = extQuants.exteriorIndex();
265 unsigned J = elemCtx.globalSpaceIndex(j, 0);
267 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
268 Scalar weight = std::max<Scalar>(1e-16,
269 std::abs(getValue(extQuants.volumeFlux(phaseIdx))));
270 Valgrind::CheckDefined(extQuants.extrusionFactor());
271 assert(extQuants.extrusionFactor() > 0);
272 weight *= extQuants.extrusionFactor();
274 const auto& inputV = extQuants.filterVelocity(phaseIdx);
276 for (
unsigned k = 0; k < dimWorld; ++k) {
277 v[k] = getValue(inputV[k]);
279 if (v.two_norm() > 1e-20) {
280 weight /= v.two_norm();
284 velocity_[phaseIdx][I] += v;
285 velocity_[phaseIdx][J] += v;
287 velocityWeight_[phaseIdx][I] += weight;
288 velocityWeight_[phaseIdx][J] += weight;
337 size_t numDof = this->
simulator_.model().numGridDof();
339 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
342 for (
unsigned i = 0; i < numDof; ++i) {
343 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
347 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
349 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
354 size_t numDof = this->
simulator_.model().numGridDof();
356 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
359 for (
unsigned i = 0; i < numDof; ++i) {
360 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
364 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx).data());
366 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
367 potentialGradient_[phaseIdx],
388 ScalarBuffer extrusionFactor_{};
389 PhaseBuffer pressure_{};
390 PhaseBuffer density_{};
391 PhaseBuffer saturation_{};
392 PhaseBuffer mobility_{};
393 PhaseBuffer relativePermeability_{};
394 PhaseBuffer viscosity_{};
395 PhaseBuffer averageMolarMass_{};
397 ScalarBuffer porosity_{};
398 TensorBuffer intrinsicPermeability_{};
400 PhaseVectorBuffer velocity_{};
401 PhaseBuffer velocityWeight_{};
403 PhaseVectorBuffer potentialGradient_{};
404 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.hpp:68
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hpp:177
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hpp:104
VtkMultiPhaseModule(const Simulator &simulator)
Definition: vtkmultiphasemodule.hpp:95
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hpp:113
bool needExtensiveQuantities() const final
Returns true iff the module needs to access the extensive quantities of a context to do its job.
Definition: vtkmultiphasemodule.hpp:381
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hpp:297
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 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