My Project
Loading...
Searching...
No Matches
baseoutputmodule.hh
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 EWOMS_BASE_OUTPUT_MODULE_HH
28#define EWOMS_BASE_OUTPUT_MODULE_HH
29
30#include <dune/common/fvector.hh>
31
32#include <dune/istl/bvector.hh>
33
35
37
39
43
44#include <array>
45#include <cstdio>
46#include <sstream>
47#include <string>
48
49namespace Opm::Properties {
50
51template <class TypeTag, class MyTypeTag>
52struct FluidSystem;
53
54} // namespace Opm::Properties
55
56namespace Opm {
57
65template<class TypeTag>
67{
75
77 enum { numComponents = getPropValue<TypeTag, Properties::NumComponents>() };
79 enum { dim = GridView::dimension };
80 enum { dimWorld = GridView::dimensionworld };
81 using Vector = BaseOutputWriter::Vector;
82 using Tensor = BaseOutputWriter::Tensor;
83
84public:
85 using ScalarBuffer = BaseOutputWriter::ScalarBuffer;
86 using VectorBuffer = BaseOutputWriter::VectorBuffer;
87 using TensorBuffer = BaseOutputWriter::TensorBuffer;
88
89 using EqBuffer = std::array<ScalarBuffer, numEq>;
90 using PhaseBuffer = std::array<ScalarBuffer, numPhases>;
91 using ComponentBuffer = std::array<ScalarBuffer, numComponents>;
92 using PhaseComponentBuffer = std::array<std::array<ScalarBuffer, numComponents>, numPhases>;
93
94 using PhaseVectorBuffer = std::array<VectorBuffer, numPhases>;
95
96 BaseOutputModule(const Simulator& simulator)
97 : simulator_(simulator)
98 {}
99
100 virtual ~BaseOutputModule()
101 {}
102
111 virtual void allocBuffers() = 0;
112
121 virtual void processElement(const ElementContext& elemCtx) = 0;
122
127
138 virtual bool needExtensiveQuantities() const
139 { return false; }
140
141protected:
152
156 void resizeScalarBuffer_(ScalarBuffer& buffer,
158 {
159 buffer.resize(this->getBufferSize(bufferType));
160 std::fill(buffer.begin(), buffer.end(), 0.0);
161 }
162
166 void resizeTensorBuffer_(TensorBuffer& buffer,
168 {
169 buffer.resize(this->getBufferSize(bufferType));
170 Tensor nullMatrix(dimWorld, dimWorld, 0.0);
171 std::fill(buffer.begin(), buffer.end(), nullMatrix);
172 }
173
174 void resizeVectorBuffer_(VectorBuffer& buffer,
176 {
177 buffer.resize(this->getBufferSize(bufferType));
178 Vector zerovector(dimWorld,0.0);
179 zerovector = 0.0;
180 std::fill(buffer.begin(), buffer.end(), zerovector);
181 }
182
187 void resizeEqBuffer_(EqBuffer& buffer,
189 {
190 const std::size_t n = this->getBufferSize(bufferType);
191 for (unsigned i = 0; i < numEq; ++i) {
192 buffer[i].resize(n);
193 std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
194 }
195 }
196
201 void resizePhaseBuffer_(PhaseBuffer& buffer,
203 {
204 const std::size_t n = this->getBufferSize(bufferType);
205 for (unsigned i = 0; i < numPhases; ++i) {
206 buffer[i].resize(n);
207 std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
208 }
209 }
210
215 void resizeComponentBuffer_(ComponentBuffer& buffer,
217 {
218 const std::size_t n = this->getBufferSize(bufferType);
219 for (unsigned i = 0; i < numComponents; ++i) {
220 buffer[i].resize(n);
221 std::fill(buffer[i].begin(), buffer[i].end(), 0.0);
222 }
223 }
224
229 void resizePhaseComponentBuffer_(PhaseComponentBuffer& buffer,
231 {
232 const std::size_t n = this->getBufferSize(bufferType);
233 for (unsigned i = 0; i < numPhases; ++i) {
234 for (unsigned j = 0; j < numComponents; ++j) {
235 buffer[i][j].resize(n);
236 std::fill(buffer[i][j].begin(), buffer[i][j].end(), 0.0);
237 }
238 }
239 }
240
245 const char *name,
246 ScalarBuffer& buffer,
248 {
249 switch (bufferType) {
250 case DofBuffer:
251 DiscBaseOutputModule::attachScalarDofData_(baseWriter, buffer, name);
252 break;
253 case VertexBuffer:
254 attachScalarVertexData_(baseWriter, buffer, name);
255 break;
256 case ElementBuffer:
257 attachScalarElementData_(baseWriter, buffer, name);
258 break;
259 default: throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
260 }
261 }
262
267 const char *name,
268 VectorBuffer& buffer,
270 {
271 switch (bufferType) {
272 case DofBuffer:
273 DiscBaseOutputModule::attachVectorDofData_(baseWriter, buffer, name);
274 break;
275 case VertexBuffer:
276 attachVectorVertexData_(baseWriter, buffer, name);
277 break;
278 case ElementBuffer:
279 attachVectorElementData_(baseWriter, buffer, name);
280 break;
281 default: throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
282 }
283 }
284
289 const char *name,
290 TensorBuffer& buffer,
292 {
293 switch (bufferType) {
294 case DofBuffer:
295 DiscBaseOutputModule::attachTensorDofData_(baseWriter, buffer, name);
296 break;
297 case VertexBuffer:
298 attachTensorVertexData_(baseWriter, buffer, name);
299 break;
300 case ElementBuffer:
301 attachTensorElementData_(baseWriter, buffer, name);
302 break;
303 default: throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
304 }
305 }
306
311 const char *pattern,
312 EqBuffer& buffer,
314 {
315 char name[512];
316 for (unsigned i = 0; i < numEq; ++i) {
317 std::string eqName = simulator_.model().primaryVarName(i);
318 snprintf(name, 512, pattern, eqName.c_str());
319
320 this->commitScalarBuffer_(baseWriter, name, buffer[i], bufferType);
321 }
322 }
323
328 const char *pattern,
329 EqBuffer& buffer,
331 {
332 char name[512];
333 for (unsigned i = 0; i < numEq; ++i) {
334 std::ostringstream oss;
335 oss << i;
336 snprintf(name, 512, pattern, oss.str().c_str());
337
338 this->commitScalarBuffer_(baseWriter, name, buffer[i], bufferType);
339 }
340 }
341
346 const char *pattern,
347 PhaseBuffer& buffer,
349 {
350 char name[512];
351 for (unsigned i = 0; i < numPhases; ++i) {
352 snprintf(name, 512, pattern, FluidSystem::phaseName(i).data());
353
354 this->commitScalarBuffer_(baseWriter, name, buffer[i], bufferType);
355 }
356 }
357
362 const char *pattern,
363 ComponentBuffer& buffer,
365 {
366 char name[512];
367 for (unsigned i = 0; i < numComponents; ++i) {
368 snprintf(name, 512, pattern, FluidSystem::componentName(i).data());
369
370 this->commitScalarBuffer_(baseWriter, name, buffer[i], bufferType);
371 }
372 }
373
378 const char *pattern,
379 PhaseComponentBuffer& buffer,
381 {
382 char name[512];
383 for (unsigned i= 0; i < numPhases; ++i) {
384 for (unsigned j = 0; j < numComponents; ++j) {
385 snprintf(name, 512, pattern,
386 FluidSystem::phaseName(i).data(),
387 FluidSystem::componentName(j).data());
388
389 this->commitScalarBuffer_(baseWriter, name, buffer[i][j], bufferType);
390 }
391 }
392 }
393
394 void attachScalarElementData_(BaseOutputWriter& baseWriter,
395 ScalarBuffer& buffer,
396 const char *name)
397 { baseWriter.attachScalarElementData(buffer, name); }
398
399 void attachScalarVertexData_(BaseOutputWriter& baseWriter,
400 ScalarBuffer& buffer,
401 const char *name)
402 { baseWriter.attachScalarVertexData(buffer, name); }
403
404 void attachVectorElementData_(BaseOutputWriter& baseWriter,
405 VectorBuffer& buffer,
406 const char *name)
407 { baseWriter.attachVectorElementData(buffer, name); }
408
409 void attachVectorVertexData_(BaseOutputWriter& baseWriter,
410 VectorBuffer& buffer,
411 const char *name)
412 { baseWriter.attachVectorVertexData(buffer, name); }
413
414 void attachTensorElementData_(BaseOutputWriter& baseWriter,
415 TensorBuffer& buffer,
416 const char *name)
417 { baseWriter.attachTensorElementData(buffer, name); }
418
419 void attachTensorVertexData_(BaseOutputWriter& baseWriter,
420 TensorBuffer& buffer,
421 const char *name)
422 { baseWriter.attachTensorVertexData(buffer, name); }
423
424 std::size_t getBufferSize(BufferType bufferType) const
425 {
426 switch (bufferType) {
427 case VertexBuffer: return simulator_.gridView().size(dim);
428 case ElementBuffer: return simulator_.gridView().size(0);
429 case DofBuffer: return simulator_.model().numGridDof();
430 default: throw std::logic_error("bufferType must be one of Dof, Vertex or Element");
431 }
432 }
433
434 const Simulator& simulator_;
435};
436
437} // namespace Opm
438
439#endif
The base class for all output writers.
Defines a type tags and some fundamental properties all models.
The base class for writer modules.
Definition baseoutputmodule.hh:67
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
void commitComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a component-specific buffer to the result file.
Definition baseoutputmodule.hh:361
virtual bool needExtensiveQuantities() const
Returns true iff the module needs to access the extensive quantities of a context to do its job.
Definition baseoutputmodule.hh:138
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
void resizeComponentBuffer_(ComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a component specific quantity.
Definition baseoutputmodule.hh:215
void commitEqBuffer_(BaseOutputWriter &baseWriter, const char *pattern, EqBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer with as many variables as PDEs to the result file.
Definition baseoutputmodule.hh:327
void resizePhaseComponentBuffer_(PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a phase and component specific buffer.
Definition baseoutputmodule.hh:229
BufferType
Definition baseoutputmodule.hh:142
@ DofBuffer
Buffer contains data associated with the degrees of freedom.
Definition baseoutputmodule.hh:144
@ VertexBuffer
Buffer contains data associated with the grid's vertices.
Definition baseoutputmodule.hh:147
@ ElementBuffer
Buffer contains data associated with the grid's elements.
Definition baseoutputmodule.hh:150
void commitPriVarsBuffer_(BaseOutputWriter &baseWriter, const char *pattern, EqBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer with as many variables as PDEs to the result file.
Definition baseoutputmodule.hh:310
void resizeEqBuffer_(EqBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a equation specific quantity.
Definition baseoutputmodule.hh:187
void commitPhaseComponentBuffer_(BaseOutputWriter &baseWriter, const char *pattern, PhaseComponentBuffer &buffer, BufferType bufferType=DofBuffer)
Add a phase and component specific quantities to the output.
Definition baseoutputmodule.hh:377
virtual void processElement(const ElementContext &elemCtx)=0
Modify the internal buffers according to the intensive quanties relevant for an element.
void resizeTensorBuffer_(TensorBuffer &buffer, BufferType bufferType=DofBuffer)
Allocate the space for a buffer storing a tensorial quantity.
Definition baseoutputmodule.hh:166
virtual void commitBuffers(BaseOutputWriter &writer)=0
Add all buffers to the VTK output writer.
void commitVectorBuffer_(BaseOutputWriter &baseWriter, const char *name, VectorBuffer &buffer, BufferType bufferType=DofBuffer)
Add a buffer containing vectorial quantities to the result file.
Definition baseoutputmodule.hh:266
virtual void allocBuffers()=0
Allocate memory for the scalar fields we would like to write to disk.
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
Declare the properties used by the infrastructure code of the finite volume discretizations.
Defines the common properties required by the porous medium multi-phase models.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
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.