vtkmultiphasemodule.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 2 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18
19 Consult the COPYING file in the top-level source directory of this
20 module for the precise wording of the license and the list of
21 copyright holders.
22*/
27#ifndef OPM_VTK_MULTI_PHASE_MODULE_HPP
28#define OPM_VTK_MULTI_PHASE_MODULE_HPP
29
30#include <dune/common/fvector.hh>
31
32#include <opm/material/common/MathToolbox.hpp>
33#include <opm/material/common/Valgrind.hpp>
34
36
40
43
44#include <algorithm>
45#include <array>
46#include <cassert>
47#include <cmath>
48#include <cstddef>
49#include <cstdio>
50
51namespace Opm {
52
71template<class TypeTag>
73{
75
79
83
84 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
86
87 enum { dimWorld = GridView::dimensionworld };
88 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
89
90 using BufferType = typename ParentType::BufferType;
92 using VectorBuffer = typename ParentType::VectorBuffer;
94 using PhaseBuffer = typename ParentType::PhaseBuffer;
95
96 using DimVector = Dune::FieldVector<Scalar, dimWorld>;
97
98 using PhaseVectorBuffer = std::array<VectorBuffer, numPhases>;
99
100public:
101 explicit VtkMultiPhaseModule(const Simulator& simulator)
102 : ParentType(simulator)
103 {
104 params_.read();
105 }
106
110 static void registerParameters()
111 {
113 }
114
119 void allocBuffers() override
120 {
121 if (params_.extrusionFactorOutput_) {
122 this->resizeScalarBuffer_(extrusionFactor_, BufferType::Dof);
123 }
124 if (params_.pressureOutput_) {
125 this->resizePhaseBuffer_(pressure_, BufferType::Dof);
126 }
127 if (params_.densityOutput_) {
128 this->resizePhaseBuffer_(density_, BufferType::Dof);
129 }
130 if (params_.saturationOutput_) {
131 this->resizePhaseBuffer_(saturation_, BufferType::Dof);
132 }
133 if (params_.mobilityOutput_) {
134 this->resizePhaseBuffer_(mobility_, BufferType::Dof);
135 }
136 if (params_.relativePermeabilityOutput_) {
137 this->resizePhaseBuffer_(relativePermeability_, BufferType::Dof);
138 }
139 if (params_.viscosityOutput_) {
140 this->resizePhaseBuffer_(viscosity_, BufferType::Dof);
141 }
142 if (params_.averageMolarMassOutput_) {
143 this->resizePhaseBuffer_(averageMolarMass_, BufferType::Dof);
144 }
145
146 if (params_.porosityOutput_) {
147 this->resizeScalarBuffer_(porosity_, BufferType::Dof);
148 }
149 if (params_.intrinsicPermeabilityOutput_) {
150 this->resizeTensorBuffer_(intrinsicPermeability_, BufferType::Dof);
151 }
152
153 if (params_.velocityOutput_) {
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;
160 }
161 }
162 this->resizePhaseBuffer_(velocityWeight_, BufferType::Dof);
163 }
164
165 if (params_.potentialGradientOutput_) {
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;
172 }
173 }
174
175 this->resizePhaseBuffer_(potentialWeight_, BufferType::Dof);
176 }
177 }
178
183 void processElement(const ElementContext& elemCtx) override
184 {
185 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
186 return;
187 }
188
189 const auto& problem = elemCtx.problem();
190 for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
191 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
192 const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
193 const auto& fs = intQuants.fluidState();
194
195 if (params_.extrusionFactorOutput_) {
196 extrusionFactor_[I] = intQuants.extrusionFactor();
197 }
198 if (params_.porosityOutput_) {
199 porosity_[I] = getValue(intQuants.porosity());
200 }
201
202 if (params_.intrinsicPermeabilityOutput_) {
203 const auto& K = problem.intrinsicPermeability(elemCtx, i, /*timeIdx=*/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];
207 }
208 }
209 }
210
211 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
212 if (!FluidSystem::phaseIsActive(phaseIdx)) {
213 continue;
214 }
215 if (params_.pressureOutput_) {
216 pressure_[phaseIdx][I] = getValue(fs.pressure(phaseIdx));
217 }
218 if (params_.densityOutput_) {
219 density_[phaseIdx][I] = getValue(fs.density(phaseIdx));
220 }
221 if (params_.saturationOutput_) {
222 saturation_[phaseIdx][I] = getValue(fs.saturation(phaseIdx));
223 }
224 if (params_.mobilityOutput_) {
225 mobility_[phaseIdx][I] = getValue(intQuants.mobility(phaseIdx));
226 }
227 if (params_.relativePermeabilityOutput_) {
228 relativePermeability_[phaseIdx][I] = getValue(intQuants.relativePermeability(phaseIdx));
229 }
230 if (params_.viscosityOutput_) {
231 viscosity_[phaseIdx][I] = getValue(fs.viscosity(phaseIdx));
232 }
233 if (params_.averageMolarMassOutput_) {
234 averageMolarMass_[phaseIdx][I] = getValue(fs.averageMolarMass(phaseIdx));
235 }
236 }
237 }
238
239 if (params_.potentialGradientOutput_) {
240 // calculate velocities if requested
241 for (unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++faceIdx) {
242 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, /*timeIdx=*/0);
243
244 const unsigned i = extQuants.interiorIndex();
245 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
246
247 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
248 const Scalar weight = extQuants.extrusionFactor();
249
250 potentialWeight_[phaseIdx][I] += weight;
251
252 const auto& inputPGrad = extQuants.potentialGrad(phaseIdx);
253 DimVector pGrad;
254 for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
255 pGrad[dimIdx] = getValue(inputPGrad[dimIdx]) * weight;
256 }
257 potentialGradient_[phaseIdx][I] += pGrad;
258 } // end for all phases
259 } // end for all faces
260 }
261
262 if (params_.velocityOutput_) {
263 // calculate velocities if requested
264 for (unsigned faceIdx = 0; faceIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++faceIdx) {
265 const auto& extQuants = elemCtx.extensiveQuantities(faceIdx, /*timeIdx=*/0);
266
267 const unsigned i = extQuants.interiorIndex();
268 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
269
270 const unsigned j = extQuants.exteriorIndex();
271 const unsigned J = elemCtx.globalSpaceIndex(j, /*timeIdx=*/0);
272
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();
279
280 const auto& inputV = extQuants.filterVelocity(phaseIdx);
281 DimVector v;
282 for (unsigned k = 0; k < dimWorld; ++k) {
283 v[k] = getValue(inputV[k]);
284 }
285 if (v.two_norm() > 1e-20) {
286 weight /= v.two_norm();
287 }
288 v *= weight;
289
290 velocity_[phaseIdx][I] += v;
291 velocity_[phaseIdx][J] += v;
292
293 velocityWeight_[phaseIdx][I] += weight;
294 velocityWeight_[phaseIdx][J] += weight;
295 } // end for all phases
296 } // end for all faces
297 }
298 }
299
303 void commitBuffers(BaseOutputWriter& baseWriter) override
304 {
305 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
306 return;
307 }
308
309 if (params_.extrusionFactorOutput_) {
310 this->commitScalarBuffer_(baseWriter, "extrusionFactor",
311 extrusionFactor_, BufferType::Dof);
312 }
313 if (params_.pressureOutput_) {
314 this->commitPhaseBuffer_(baseWriter, "pressure_%s", pressure_, BufferType::Dof);
315 }
316 if (params_.densityOutput_) {
317 this->commitPhaseBuffer_(baseWriter, "density_%s", density_, BufferType::Dof);
318 }
319 if (params_.saturationOutput_) {
320 this->commitPhaseBuffer_(baseWriter, "saturation_%s", saturation_, BufferType::Dof);
321 }
322 if (params_.mobilityOutput_) {
323 this->commitPhaseBuffer_(baseWriter, "mobility_%s", mobility_, BufferType::Dof);
324 }
325 if (params_.relativePermeabilityOutput_) {
326 this->commitPhaseBuffer_(baseWriter, "relativePerm_%s",
327 relativePermeability_, BufferType::Dof);
328 }
329 if (params_.viscosityOutput_) {
330 this->commitPhaseBuffer_(baseWriter, "viscosity_%s", viscosity_, BufferType::Dof);
331 }
332 if (params_.averageMolarMassOutput_) {
333 this->commitPhaseBuffer_(baseWriter, "averageMolarMass_%s",
334 averageMolarMass_, BufferType::Dof);
335 }
336
337 if (params_.porosityOutput_) {
338 this->commitScalarBuffer_(baseWriter, "porosity", porosity_, BufferType::Dof);
339 }
340 if (params_.intrinsicPermeabilityOutput_) {
341 this->commitTensorBuffer_(baseWriter, "intrinsicPerm",
342 intrinsicPermeability_, BufferType::Dof);
343 }
344
345 if (params_.velocityOutput_) {
346 const std::size_t numDof = this->simulator_.model().numGridDof();
347 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
348 // first, divide the velocity field by the
349 // respective finite volume's surface area
350 for (unsigned i = 0; i < numDof; ++i) {
351 velocity_[phaseIdx][i] /= velocityWeight_[phaseIdx][i];
352 }
353 // commit the phase velocity
354 char name[512];
355 snprintf(name, 512, "filterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
356
357 DiscBaseOutputModule::attachVectorDofData_(baseWriter, velocity_[phaseIdx], name);
358 }
359 }
360
361 if (params_.potentialGradientOutput_) {
362 const std::size_t numDof = this->simulator_.model().numGridDof();
363 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
364 // first, divide the velocity field by the
365 // respective finite volume's surface area
366 for (unsigned i = 0; i < numDof; ++i) {
367 potentialGradient_[phaseIdx][i] /= potentialWeight_[phaseIdx][i];
368 }
369 // commit the phase velocity
370 char name[512];
371 snprintf(name, 512, "gradP_%s", FluidSystem::phaseName(phaseIdx).data());
372
373 DiscBaseOutputModule::attachVectorDofData_(baseWriter,
374 potentialGradient_[phaseIdx],
375 name);
376 }
377 }
378 }
379
388 bool needExtensiveQuantities() const override
389 {
390 return params_.velocityOutput_ || params_.potentialGradientOutput_;
391 }
392
393private:
394 VtkMultiPhaseParams params_{};
395 ScalarBuffer extrusionFactor_{};
396 PhaseBuffer pressure_{};
397 PhaseBuffer density_{};
398 PhaseBuffer saturation_{};
399 PhaseBuffer mobility_{};
400 PhaseBuffer relativePermeability_{};
401 PhaseBuffer viscosity_{};
402 PhaseBuffer averageMolarMass_{};
403
404 ScalarBuffer porosity_{};
405 TensorBuffer intrinsicPermeability_{};
406
407 PhaseVectorBuffer velocity_{};
408 PhaseBuffer velocityWeight_{};
409
410 PhaseVectorBuffer potentialGradient_{};
411 PhaseBuffer potentialWeight_{};
412};
413
414} // namespace Opm
415
416#endif // OPM_VTK_MULTI_PHASE_MODULE_HPP
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