3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
griddata.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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program 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 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program 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 this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_IO_GRID_DATA_HH
25#define DUMUX_IO_GRID_DATA_HH
26
27#include <vector>
28#include <memory>
29#include <type_traits>
30
31#include <dune/common/exceptions.hh>
32#include <dune/common/parallel/communication.hh>
33#include <dune/common/parallel/mpihelper.hh>
34#include <dune/grid/common/gridfactory.hh>
35#include <dune/grid/io/file/dgfparser/parser.hh>
36#include <dune/grid/io/file/dgfparser/gridptr.hh>
38
39// UGGrid specific includes
40#if HAVE_UG
41#include <dune/grid/uggrid.hh>
42#endif
43
44#include "gmshgriddatahandle.hh"
45
46namespace Dumux {
47
48namespace Detail {
49
50template<class Grid>
51struct isUG : public std::false_type {};
52
53#if HAVE_UG
54template<int dim>
55struct isUG<Dune::UGGrid<dim>> : public std::true_type {};
56#endif
57
58} // end namespace Details
59
64template <class Grid>
66{
67 using Intersection = typename Grid::LeafIntersection;
68 using Element = typename Grid::template Codim<0>::Entity;
69 using Vertex = typename Grid::template Codim<Grid::dimension>::Entity;
71
72 enum DataSourceType { dgf, gmsh, vtk };
73
74public:
76 GridData(std::shared_ptr<Grid> grid, std::shared_ptr<Dune::GridFactory<Grid>> factory,
77 std::vector<int>&& elementMarkers, std::vector<int>&& boundaryMarkers, std::vector<int>&& faceMarkers = std::vector<int>{})
78 : gridPtr_(grid)
79 , gridFactory_(factory)
80 , elementMarkers_(elementMarkers)
81 , boundaryMarkers_(boundaryMarkers)
82 , faceMarkers_(faceMarkers)
83 , dataSourceType_(DataSourceType::gmsh)
84 {}
85
87 GridData(Dune::GridPtr<Grid> grid)
88 : dgfGrid_(grid)
89 , dataSourceType_(DataSourceType::dgf)
90 {}
91
93 GridData(std::shared_ptr<Grid> grid, std::shared_ptr<Dune::GridFactory<Grid>> factory,
94 VTKReader::Data&& cellData, VTKReader::Data&& pointData)
95 : gridPtr_(grid)
96 , gridFactory_(factory)
97 , cellData_(cellData)
98 , pointData_(pointData)
99 , dataSourceType_(DataSourceType::vtk)
100 {}
101
102
106 // \{
107
112 const std::vector<double>& parameters(const Vertex& vertex) const
113 {
114 if (dataSourceType_ == DataSourceType::dgf)
115 {
116 if (vertex.level() != 0)
117 DUNE_THROW(Dune::IOError, "You can only obtain parameters for level 0 vertices!");
118
119 return dgfGrid_.parameters(vertex);
120 }
121 else
122 DUNE_THROW(Dune::InvalidStateException, "The parameters method is only available if the grid was constructed with a DGF file.");
123 }
124
128 const std::vector<double>& parameters(const Element& element) const
129 {
130 if (dataSourceType_ == DataSourceType::dgf)
131 {
132 if (element.hasFather())
133 {
134 auto level0Element = element;
135 while(level0Element.hasFather())
136 level0Element = level0Element.father();
137
138 return dgfGrid_.parameters(level0Element);
139 }
140 else
141 {
142 return dgfGrid_.parameters(element);
143 }
144 }
145 else
146 DUNE_THROW(Dune::InvalidStateException, "The parameters method is only available if the grid was constructed with a DGF file.");
147 }
148
152 template <class GridImp, class IntersectionImp>
153 const Dune::DGFBoundaryParameter::type& parameters(const Dune::Intersection<GridImp, IntersectionImp>& intersection) const
154 {
155 if (dataSourceType_ == DataSourceType::dgf)
156 return dgfGrid_.parameters(intersection);
157 else
158 DUNE_THROW(Dune::InvalidStateException, "The parameters method is only available if the grid was constructed with a DGF file.");
159 }
160
161 // \}
162
166 // \{
167
173 int getBoundaryDomainMarker(int boundarySegmentIndex) const
174 {
175 if (dataSourceType_ != DataSourceType::gmsh)
176 DUNE_THROW(Dune::InvalidStateException, "Domain markers are only available for gmsh grids.");
177 if (boundarySegmentIndex >= boundaryMarkers_.size())
178 DUNE_THROW(Dune::RangeError, "Boundary segment index "<< boundarySegmentIndex << " bigger than number of boundary segments in grid.\n"
179 "Make sure to call this function only for boundaries that were defined as physical entities in gmsh.");
180 return boundaryMarkers_[boundarySegmentIndex];
181 }
182
188 int getBoundaryDomainMarker(const Intersection& intersection) const
189 { return getBoundaryDomainMarker(intersection.boundarySegmentIndex()); }
190
194 bool wasInserted(const Intersection& intersection) const
195 { return gridFactory_->wasInserted(intersection); }
196
202 int getElementDomainMarker(const Element& element) const
203 {
204 if (dataSourceType_ != DataSourceType::gmsh)
205 DUNE_THROW(Dune::InvalidStateException, "Domain markers are only available for gmsh grids.");
206
207 // parameters are only given for level 0 elements
208 auto level0element = element;
209 while (level0element.hasFather())
210 level0element = level0element.father();
211
212 // in the parallel case the data is load balanced and then accessed with indices of the index set
213 // for UGGrid element data is read on all processes since UGGrid can't communicate element data (yet)
214 if (gridPtr_->comm().size() > 1 && !Detail::isUG<Grid>::value)
215 return elementMarkers_[gridPtr_->levelGridView(0).indexSet().index(level0element)];
216 else
217 return elementMarkers_[gridFactory_->insertionIndex(level0element)];
218 }
219
224 template<bool ug = Detail::isUG<Grid>::value, typename std::enable_if_t<!ug, int> = 0>
226 {
227 return DataHandle(*gridPtr_, *gridFactory_, elementMarkers_, boundaryMarkers_, faceMarkers_);
228 }
229
234 template<bool ug = Detail::isUG<Grid>::value, typename std::enable_if_t<ug, int> = 0>
236 {
237 return DataHandle(*gridPtr_, *gridFactory_, elementMarkers_, boundaryMarkers_);
238 }
239
240 // \}
241
245 // \{
246
252 double getParameter(const Element& element, const std::string& fieldName) const
253 {
254 if (dataSourceType_ != DataSourceType::vtk)
255 DUNE_THROW(Dune::InvalidStateException, "This access function is only available for data from VTK files.");
256
257 if (cellData_.count(fieldName) == 0)
258 DUNE_THROW(Dune::IOError, "No field with the name " << fieldName << " found in cell data");
259
260 // parameters are only given for level 0 elements
261 auto level0element = element;
262 while (level0element.hasFather())
263 level0element = level0element.father();
264
265 return cellData_.at(fieldName)[gridFactory_->insertionIndex(level0element)];
266 }
267
274 double getParameter(const Vertex& vertex, const std::string& fieldName) const
275 {
276 if (dataSourceType_ != DataSourceType::vtk)
277 DUNE_THROW(Dune::InvalidStateException, "This access function is only available for data from VTK files.");
278
279 if (vertex.level() != 0)
280 DUNE_THROW(Dune::IOError, "You can only obtain parameters for level 0 vertices!");
281
282 if (pointData_.count(fieldName) == 0)
283 DUNE_THROW(Dune::IOError, "No field with the name " << fieldName << " found in point data");
284
285 return pointData_.at(fieldName)[gridFactory_->insertionIndex(vertex)];
286 }
287
288 // \}
289
290private:
291 // grid and grid factor for gmsh grid data / vtk grid data
292 std::shared_ptr<Grid> gridPtr_;
293 std::shared_ptr<Dune::GridFactory<Grid>> gridFactory_;
294
299 std::vector<int> elementMarkers_;
300 std::vector<int> boundaryMarkers_;
301 std::vector<int> faceMarkers_;
302
306 VTKReader::Data cellData_, pointData_;
307
308 // dgf grid data
309 Dune::GridPtr<Grid> dgfGrid_;
310
311 // specify which type of data we have
312 // TODO unfortunately all grid readers provide different data types, should be streamlined (changes in Dune)
313 DataSourceType dataSourceType_;
314};
315
316} // namespace Dumux
317
318#endif
A data handle for commucating grid data for gmsh grids.
A vtk file reader using tinyxml2 as xml backend.
Definition: adapt.hh:29
Definition: common/pdesolver.hh:35
A data handle for commucating grid data for gmsh grids.
Definition: gmshgriddatahandle.hh:49
Definition: griddata.hh:51
Class for grid data attached to dgf or gmsh grid files.
Definition: griddata.hh:66
GridData(std::shared_ptr< Grid > grid, std::shared_ptr< Dune::GridFactory< Grid > > factory, std::vector< int > &&elementMarkers, std::vector< int > &&boundaryMarkers, std::vector< int > &&faceMarkers=std::vector< int >{})
constructor for gmsh grid data
Definition: griddata.hh:76
double getParameter(const Vertex &vertex, const std::string &fieldName) const
Call the parameters function of the DGF grid pointer if available for vertex data.
Definition: griddata.hh:274
int getBoundaryDomainMarker(const Intersection &intersection) const
Return the boundary domain marker (Gmsh physical entity number) of an intersection Only available whe...
Definition: griddata.hh:188
double getParameter(const Element &element, const std::string &fieldName) const
Get a element parameter.
Definition: griddata.hh:252
GridData(std::shared_ptr< Grid > grid, std::shared_ptr< Dune::GridFactory< Grid > > factory, VTKReader::Data &&cellData, VTKReader::Data &&pointData)
constructor for gmsh grid data
Definition: griddata.hh:93
const Dune::DGFBoundaryParameter::type & parameters(const Dune::Intersection< GridImp, IntersectionImp > &intersection) const
Call the parameters function of the DGF grid pointer if available.
Definition: griddata.hh:153
GridData(Dune::GridPtr< Grid > grid)
constructor for dgf grid data
Definition: griddata.hh:87
DataHandle createGmshDataHandle()
Create a data handle for communication of the data in parallel simulations.
Definition: griddata.hh:225
const std::vector< double > & parameters(const Vertex &vertex) const
Call the parameters function of the DGF grid pointer if available for vertex data.
Definition: griddata.hh:112
int getBoundaryDomainMarker(int boundarySegmentIndex) const
Return the boundary domain marker (Gmsh physical entity number) of an intersection Only available whe...
Definition: griddata.hh:173
int getElementDomainMarker(const Element &element) const
Return the element domain marker (Gmsh physical entity number) of an element. Only available when usi...
Definition: griddata.hh:202
const std::vector< double > & parameters(const Element &element) const
Call the parameters function of the DGF grid pointer if available for element data.
Definition: griddata.hh:128
bool wasInserted(const Intersection &intersection) const
Returns true if an intersection was inserted during grid creation.
Definition: griddata.hh:194
std::unordered_map< std::string, std::vector< double > > Data
the cell / point data type for point data read from a grid file
Definition: vtkreader.hh:59