25 #ifndef EWOMS_VTK_MULTI_PHASE_MODULE_HH
26 #define EWOMS_VTK_MULTI_PHASE_MODULE_HH
34 #include <opm/material/common/MathToolbox.hpp>
36 #include <dune/common/fvector.hh>
41 namespace Properties {
65 SET_BOOL_PROP(VtkMultiPhase, VtkWriteRelativePermeabilities,
true);
67 SET_BOOL_PROP(VtkMultiPhase, VtkWriteAverageMolarMasses,
false);
69 SET_BOOL_PROP(VtkMultiPhase, VtkWriteIntrinsicPermeabilities,
false);
70 SET_BOOL_PROP(VtkMultiPhase, VtkWritePotentialGradients,
false);
71 SET_BOOL_PROP(VtkMultiPhase, VtkWriteFilterVelocities,
false);
92 template<
class TypeTag>
98 typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
99 typedef typename GET_PROP_TYPE(TypeTag, Evaluation) Evaluation;
100 typedef typename GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
102 typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
103 typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
104 typedef typename GET_PROP_TYPE(TypeTag, DiscBaseOutputModule) DiscBaseOutputModule;
106 static const int vtkFormat =
GET_PROP_VALUE(TypeTag, VtkOutputFormat);
109 typedef Opm::MathToolbox<Evaluation> Toolbox;
111 enum { dimWorld = GridView::dimensionworld };
119 typedef Dune::FieldVector<Scalar, dimWorld> DimVector;
125 : ParentType(simulator)
134 "Include the phase pressures in the VTK output files");
136 "Include the phase densities in the VTK output files");
138 "Include the phase saturations in the VTK output files");
140 "Include the phase mobilities in the VTK output files");
142 "Include the phase relative permeabilities in the VTK output files");
144 "Include component phase viscosities in the VTK output files");
146 "Include the average phase mass in the VTK output files");
148 "Include the porosity in the VTK output files");
150 "Include the intrinsic permeability in the VTK output files");
152 "Include in the filter velocities of the phases the VTK output files");
154 "Include the phase pressure potential gradients in the VTK output files");
174 if (velocityOutput_()) {
176 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
177 velocity_[phaseIdx].resize(nDof);
178 for (
int dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
179 velocity_[phaseIdx][dofIdx].resize(dimWorld);
180 velocity_[phaseIdx][dofIdx] = 0.0;
186 if (potentialGradientOutput_()) {
188 for (
int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
189 potentialGradient_[phaseIdx].resize(nDof);
190 for (
int dofIdx = 0; dofIdx < nDof; ++ dofIdx) {
191 potentialGradient_[phaseIdx][dofIdx].resize(dimWorld);
192 potentialGradient_[phaseIdx][dofIdx] = 0.0;
206 typedef Opm::MathToolbox<Evaluation> Toolbox;
208 for (
int i = 0; i < elemCtx.numPrimaryDof(0); ++i) {
209 int I = elemCtx.globalSpaceIndex(i, 0);
210 const auto &intQuants = elemCtx.intensiveQuantities(i, 0);
211 const auto &fs = intQuants.fluidState();
213 if (porosityOutput_()) porosity_[I] = Toolbox::value(intQuants.porosity());
214 if (intrinsicPermeabilityOutput_()) {
215 const auto& K = intQuants.intrinsicPermeability();
216 intrinsicPermeability_[I].resize(K.rows, K.cols);
217 for (
int rowIdx = 0; rowIdx < K.rows; ++rowIdx)
218 for (
int colIdx = 0; colIdx < K.cols; ++colIdx)
219 intrinsicPermeability_[I][rowIdx][colIdx] = K[rowIdx][colIdx];
222 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
223 if (pressureOutput_())
224 pressure_[phaseIdx][I] = Toolbox::value(fs.pressure(phaseIdx));
225 if (densityOutput_())
226 density_[phaseIdx][I] = Toolbox::value(fs.density(phaseIdx));
227 if (saturationOutput_())
228 saturation_[phaseIdx][I] = Toolbox::value(fs.saturation(phaseIdx));
229 if (mobilityOutput_())
230 mobility_[phaseIdx][I] = Toolbox::value(intQuants.mobility(phaseIdx));
231 if (relativePermeabilityOutput_())
232 relativePermeability_[phaseIdx][I] = Toolbox::value(intQuants.relativePermeability(phaseIdx));
233 if (viscosityOutput_())
234 viscosity_[phaseIdx][I] = Toolbox::value(fs.viscosity(phaseIdx));
235 if (averageMolarMassOutput_())
236 averageMolarMass_[phaseIdx][I] = Toolbox::value(fs.averageMolarMass(phaseIdx));
240 if (potentialGradientOutput_()) {
242 for (
int faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
243 const auto &extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
245 int i = extQuants.interiorIndex();
246 int I = elemCtx.globalSpaceIndex(i, 0);
248 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
249 Scalar weight = extQuants.extrusionFactor();
251 potentialWeight_[phaseIdx][I] += weight;
253 const auto &inputPGrad = extQuants.potentialGrad(phaseIdx);
255 for (
int j = 0; j < numPhases; ++j)
256 pGrad[j] = Toolbox::value(inputPGrad[j])*weight;
257 potentialGradient_[phaseIdx][I] += pGrad;
262 if (velocityOutput_()) {
264 for (
int faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(0); ++ faceIdx) {
265 const auto &extQuants = elemCtx.extensiveQuantities(faceIdx, 0);
267 int i = extQuants.interiorIndex();
268 int I = elemCtx.globalSpaceIndex(i, 0);
270 int j = extQuants.exteriorIndex();
271 int J = elemCtx.globalSpaceIndex(j, 0);
273 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
274 Scalar weight = std::max<Scalar>(1e-16,
275 std::abs(Toolbox::value(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 (
int k = 0; k < dimWorld; ++k)
283 v[k] = Toolbox::value(inputV[k]);
284 if (v.two_norm() > 1e-20)
285 weight /= v.two_norm();
288 velocity_[phaseIdx][I] += v;
289 velocity_[phaseIdx][J] += v;
291 velocityWeight_[phaseIdx][I] += weight;
292 velocityWeight_[phaseIdx][J] += weight;
303 VtkMultiWriter *vtkWriter =
dynamic_cast<VtkMultiWriter*
>(&baseWriter);
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_()) {
330 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
333 for (
int i = 0; i < nDof; ++i)
334 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
337 snprintf(name, 512,
"filterVelocity_%s", FluidSystem::phaseName(phaseIdx));
339 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
343 if (potentialGradientOutput_()) {
346 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
349 for (
int i = 0; i < nDof; ++i)
350 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
353 snprintf(name, 512,
"gradP_%s", FluidSystem::phaseName(phaseIdx));
355 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
356 potentialGradient_[phaseIdx],
363 static bool pressureOutput_()
366 static bool densityOutput_()
369 static bool saturationOutput_()
372 static bool mobilityOutput_()
375 static bool relativePermeabilityOutput_()
376 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteRelativePermeabilities); }
378 static bool viscosityOutput_()
381 static bool averageMolarMassOutput_()
384 static bool porosityOutput_()
387 static bool intrinsicPermeabilityOutput_()
388 {
return EWOMS_GET_PARAM(TypeTag,
bool, VtkWriteIntrinsicPermeabilities); }
390 static bool velocityOutput_()
393 static bool potentialGradientOutput_()
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_;
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:319
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkmultiphasemodule.hh:131
void processElement(const ElementContext &elemCtx)
Modify the internal buffers according to the intensive quantities seen on an element.
Definition: vtkmultiphasemodule.hh:204
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition: baseoutputmodule.hh:168
void resizeScalarBuffer_(ScalarBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a scalar quantity.
Definition: baseoutputmodule.hh:148
void resizePhaseBuffer_(PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase-specific quantity.
Definition: baseoutputmodule.hh:213
VtkMultiPhaseModule(const Simulator &simulator)
Definition: vtkmultiphasemodule.hh:124
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:468
BaseOutputWriter::TensorBuffer TensorBuffer
Definition: baseoutputmodule.hh:92
void commitBuffers(BaseOutputWriter &baseWriter)
Add all buffers to the VTK output writer.
Definition: vtkmultiphasemodule.hh:301
std::array< VectorBuffer, numPhases > PhaseVectorBuffer
Definition: baseoutputmodule.hh:99
VTK output module for quantities which make sense for all models which deal with multiple fluid phase...
Definition: vtkmultiphasemodule.hh:93
NEW_PROP_TAG(Grid)
The type of the DUNE grid.
This file provides the infrastructure to retrieve run-time parameters.
The base class for all output writers.
Definition: baseoutputwriter.hh:41
Simplifies writing multi-file VTK datasets.
BaseOutputWriter::VectorBuffer VectorBuffer
Definition: baseoutputmodule.hh:91
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:73
Definition: baseauxiliarymodule.hh:35
#define EWOMS_REGISTER_PARAM(TypeTag, ParamType, ParamName, Description)
Register a run-time parameter.
Definition: parametersystem.hh:64
Model & model()
Return the physical model used in the simulation.
Definition: simulator.hh:176
std::array< ScalarBuffer, numPhases > PhaseBuffer
Definition: baseoutputmodule.hh:95
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:283
Provides the magic behind the eWoms property system.
The base class for writer modules.
void commitPhaseBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase-specific buffer to the result file.
Definition: baseoutputmodule.hh:386
SET_BOOL_PROP(FvBaseDiscretization, EnableVtkOutput, true)
Enable the VTK output by default.
The base class for writer modules.
Definition: baseoutputmodule.hh:71
const Simulator & simulator_
Definition: baseoutputmodule.hh:487
Simplifies writing multi-file VTK datasets.
Definition: vtkmultiwriter.hh:61
#define EWOMS_GET_PARAM(TypeTag, ParamType, ParamName)
Retrieve a runtime parameter.
Definition: parametersystem.hh:95
void allocBuffers()
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkmultiphasemodule.hh:161
BaseOutputWriter::ScalarBuffer ScalarBuffer
Definition: baseoutputmodule.hh:90