vtkdiscretefracturemodule.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_DISCRETE_FRACTURE_MODULE_HPP
28#define OPM_VTK_DISCRETE_FRACTURE_MODULE_HPP
29
30#include <dune/common/fvector.hh>
31
32#include <opm/material/common/Valgrind.hpp>
33
35
39
42
43#include <algorithm>
44#include <cassert>
45#include <cmath>
46#include <cstdio>
47#include <cstddef>
48
49namespace Opm {
50
64template <class TypeTag>
66{
68
72
76
78
79 static constexpr auto vtkFormat = getPropValue<TypeTag, Properties::VtkOutputFormat>();
81
82 enum { dim = GridView::dimension };
83 enum { dimWorld = GridView::dimensionworld };
84 enum { numPhases = getPropValue<TypeTag, Properties::NumPhases>() };
85
86 using BufferType = typename ParentType::BufferType;
87 using PhaseBuffer = typename ParentType::PhaseBuffer;
90
91public:
92 explicit VtkDiscreteFractureModule(const Simulator& simulator)
93 : ParentType(simulator)
94 {
95 params_.read();
96 }
97
101 static void registerParameters()
102 {
104 }
105
110 void allocBuffers() override
111 {
112 if (params_.saturationOutput_) {
113 this->resizePhaseBuffer_(fractureSaturation_, BufferType::Dof);
114 }
115 if (params_.mobilityOutput_) {
116 this->resizePhaseBuffer_(fractureMobility_, BufferType::Dof);
117 }
118 if (params_.relativePermeabilityOutput_) {
119 this->resizePhaseBuffer_(fractureRelativePermeability_, BufferType::Dof);
120 }
121
122 if (params_.porosityOutput_) {
123 this->resizeScalarBuffer_(fracturePorosity_, BufferType::Dof);
124 }
125 if (params_.intrinsicPermeabilityOutput_) {
126 this->resizeScalarBuffer_(fractureIntrinsicPermeability_, BufferType::Dof);
127 }
128 if (params_.volumeFractionOutput_) {
129 this->resizeScalarBuffer_(fractureVolumeFraction_, BufferType::Dof);
130 }
131
132 if (params_.velocityOutput_) {
133 const std::size_t nDof = this->simulator_.model().numGridDof();
134 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
135 fractureVelocity_[phaseIdx].resize(nDof);
136 for (unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
137 fractureVelocity_[phaseIdx][dofIdx].resize(dimWorld);
138 fractureVelocity_[phaseIdx][dofIdx] = 0.0;
139 }
140 }
141 this->resizePhaseBuffer_(fractureVelocityWeight_, BufferType::Dof);
142 }
143 }
144
149 void processElement(const ElementContext& elemCtx) override
150 {
151 if (!Parameters::Get<Parameters::EnableVtkOutput>()) {
152 return;
153 }
154
155 const auto& fractureMapper = elemCtx.simulator().vanguard().fractureMapper();
156
157 for (unsigned i = 0; i < elemCtx.numPrimaryDof(/*timeIdx=*/0); ++i) {
158 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
159 if (!fractureMapper.isFractureVertex(I)) {
160 continue;
161 }
162
163 const auto& intQuants = elemCtx.intensiveQuantities(i, /*timeIdx=*/0);
164 const auto& fs = intQuants.fractureFluidState();
165
166 if (params_.porosityOutput_) {
167 Valgrind::CheckDefined(intQuants.fracturePorosity());
168 fracturePorosity_[I] = intQuants.fracturePorosity();
169 }
170 if (params_.intrinsicPermeabilityOutput_) {
171 const auto& K = intQuants.fractureIntrinsicPermeability();
172 fractureIntrinsicPermeability_[I] = K[0][0];
173 }
174
175 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
176 if (params_.saturationOutput_) {
177 Valgrind::CheckDefined(fs.saturation(phaseIdx));
178 fractureSaturation_[phaseIdx][I] = fs.saturation(phaseIdx);
179 }
180 if (params_.mobilityOutput_) {
181 Valgrind::CheckDefined(intQuants.fractureMobility(phaseIdx));
182 fractureMobility_[phaseIdx][I] = intQuants.fractureMobility(phaseIdx);
183 }
184 if (params_.relativePermeabilityOutput_) {
185 Valgrind::CheckDefined(intQuants.fractureRelativePermeability(phaseIdx));
186 fractureRelativePermeability_[phaseIdx][I] =
187 intQuants.fractureRelativePermeability(phaseIdx);
188 }
189 if (params_.volumeFractionOutput_) {
190 Valgrind::CheckDefined(intQuants.fractureVolume());
191 fractureVolumeFraction_[I] += intQuants.fractureVolume();
192 }
193 }
194 }
195
196 if (params_.velocityOutput_) {
197 // calculate velocities if requested by the simulator
198 for (unsigned scvfIdx = 0; scvfIdx < elemCtx.numInteriorFaces(/*timeIdx=*/0); ++scvfIdx) {
199 const auto& extQuants = elemCtx.extensiveQuantities(scvfIdx, /*timeIdx=*/0);
200
201 const unsigned i = extQuants.interiorIndex();
202 const unsigned I = elemCtx.globalSpaceIndex(i, /*timeIdx=*/0);
203
204 const unsigned j = extQuants.exteriorIndex();
205 const unsigned J = elemCtx.globalSpaceIndex(j, /*timeIdx=*/0);
206
207 if (!fractureMapper.isFractureEdge(I, J)) {
208 continue;
209 }
210
211 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
212 Scalar weight = std::max(Scalar{1e-16},
213 std::abs(extQuants.fractureVolumeFlux(phaseIdx)));
214 Valgrind::CheckDefined(extQuants.extrusionFactor());
215 assert(extQuants.extrusionFactor() > 0);
216 weight *= extQuants.extrusionFactor();
217
218 Dune::FieldVector<Scalar, dim> v(extQuants.fractureFilterVelocity(phaseIdx));
219 v *= weight;
220
221 for (unsigned dimIdx = 0; dimIdx < dimWorld; ++dimIdx) {
222 fractureVelocity_[phaseIdx][I][dimIdx] += v[dimIdx];
223 fractureVelocity_[phaseIdx][J][dimIdx] += v[dimIdx];
224 }
225
226 fractureVelocityWeight_[phaseIdx][I] += weight;
227 fractureVelocityWeight_[phaseIdx][J] += weight;
228 }
229 }
230 }
231 }
232
236 void commitBuffers(BaseOutputWriter& baseWriter) override
237 {
238 if (!dynamic_cast<VtkMultiWriter*>(&baseWriter)) {
239 return;
240 }
241
242 if (params_.saturationOutput_) {
243 this->commitPhaseBuffer_(baseWriter, "fractureSaturation_%s",
244 fractureSaturation_, BufferType::Dof);
245 }
246 if (params_.mobilityOutput_) {
247 this->commitPhaseBuffer_(baseWriter, "fractureMobility_%s",
248 fractureMobility_, BufferType::Dof);
249 }
250 if (params_.relativePermeabilityOutput_) {
251 this->commitPhaseBuffer_(baseWriter, "fractureRelativePerm_%s",
252 fractureRelativePermeability_, BufferType::Dof);
253 }
254
255 if (params_.porosityOutput_) {
256 this->commitScalarBuffer_(baseWriter, "fracturePorosity",
257 fracturePorosity_, BufferType::Dof);
258 }
259 if (params_.intrinsicPermeabilityOutput_) {
260 this->commitScalarBuffer_(baseWriter, "fractureIntrinsicPerm",
261 fractureIntrinsicPermeability_, BufferType::Dof);
262 }
263 if (params_.volumeFractionOutput_) {
264 // divide the fracture volume by the total volume of the finite volumes
265 for (unsigned I = 0; I < fractureVolumeFraction_.size(); ++I) {
266 fractureVolumeFraction_[I] /= this->simulator_.model().dofTotalVolume(I);
267 }
268 this->commitScalarBuffer_(baseWriter, "fractureVolumeFraction",
269 fractureVolumeFraction_, BufferType::Dof);
270 }
271
272 if (params_.velocityOutput_) {
273 const std::size_t nDof = this->simulator_.model().numGridDof();
274 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
275 // first, divide the velocity field by the
276 // respective finite volume's surface area
277 for (unsigned dofIdx = 0; dofIdx < nDof; ++dofIdx) {
278 fractureVelocity_[phaseIdx][dofIdx] /=
279 std::max<Scalar>(1e-20, fractureVelocityWeight_[phaseIdx][dofIdx]);
280 }
281 // commit the phase velocity
282 char name[512];
283 snprintf(name, 512, "fractureFilterVelocity_%s", FluidSystem::phaseName(phaseIdx).data());
284
285 DiscBaseOutputModule::attachVectorDofData_(baseWriter, fractureVelocity_[phaseIdx], name);
286 }
287 }
288 }
289
290private:
292 PhaseBuffer fractureSaturation_{};
293 PhaseBuffer fractureMobility_{};
294 PhaseBuffer fractureRelativePermeability_{};
295
296 ScalarBuffer fracturePorosity_{};
297 ScalarBuffer fractureVolumeFraction_{};
298 ScalarBuffer fractureIntrinsicPermeability_{};
299
300 PhaseVectorBuffer fractureVelocity_{};
301 PhaseBuffer fractureVelocityWeight_{};
302
303 PhaseVectorBuffer potentialGradient_{};
304 PhaseBuffer potentialWeight_{};
305};
306
307} // namespace Opm
308
309#endif // OPM_VTK_DISCRETE_FRACTURE_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
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 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 discrete fractures i...
Definition: vtkdiscretefracturemodule.hpp:66
void processElement(const ElementContext &elemCtx) override
Modify the internal buffers according to the intensive quantities relevant for an element.
Definition: vtkdiscretefracturemodule.hpp:149
static void registerParameters()
Register all run-time parameters for the multi-phase VTK output module.
Definition: vtkdiscretefracturemodule.hpp:101
void commitBuffers(BaseOutputWriter &baseWriter) override
Add all buffers to the VTK output writer.
Definition: vtkdiscretefracturemodule.hpp:236
VtkDiscreteFractureModule(const Simulator &simulator)
Definition: vtkdiscretefracturemodule.hpp:92
void allocBuffers() override
Allocate memory for the scalar fields we would like to write to the VTK file.
Definition: vtkdiscretefracturemodule.hpp:110
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 VtkDiscreteFractureModule.
Definition: vtkdiscretefractureparams.hpp:49
void read()
Reads the parameter values from the parameter system.
bool velocityOutput_
Definition: vtkdiscretefractureparams.hpp:62
bool relativePermeabilityOutput_
Definition: vtkdiscretefractureparams.hpp:58
bool porosityOutput_
Definition: vtkdiscretefractureparams.hpp:59
static void registerParameters()
Registers the parameters in parameter system.
bool intrinsicPermeabilityOutput_
Definition: vtkdiscretefractureparams.hpp:60
bool volumeFractionOutput_
Definition: vtkdiscretefractureparams.hpp:61
bool saturationOutput_
Definition: vtkdiscretefractureparams.hpp:56
bool mobilityOutput_
Definition: vtkdiscretefractureparams.hpp:57