dune-grid  2.7.0
utility/gridinfo.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 
4 #ifndef DUNE_GRID_UTILITY_GRIDINFO_HH
5 #define DUNE_GRID_UTILITY_GRIDINFO_HH
6 
7 #include <algorithm>
8 #include <cstddef>
9 #include <functional>
10 #include <limits>
11 #include <map>
12 #include <ostream>
13 #include <string>
14 #include <vector>
15 
16 #include <dune/common/classname.hh>
17 #include <dune/common/exceptions.hh>
18 #include <dune/common/fvector.hh>
19 #include <dune/common/hybridutilities.hh>
20 #include <dune/common/std/utility.hh>
21 
22 #include <dune/geometry/multilineargeometry.hh>
23 #include <dune/geometry/referenceelements.hh>
24 #include <dune/geometry/type.hh>
25 
27 
28 namespace Dune {
29 
31  template<class ctype>
32  struct EntityInfo {
34  std::size_t count;
36 
41  ctype volumeMin;
43 
48  ctype volumeMax;
49 
51 
55  ctype volumeSum;
56 
58 
64  count(0), volumeMin(std::numeric_limits<ctype>::infinity()),
65  volumeMax(-std::numeric_limits<ctype>::infinity()), volumeSum(0)
66  { }
67  };
68 
70 
77  public std::binary_function<GeometryType, GeometryType, bool>
78  {
80  inline bool operator()(const GeometryType &a, const GeometryType &b) const
81  {
82  return a.dim() < b.dim() ||
83  (a.dim() == b.dim() && (a.isNone() < b.isNone() ||
84  (a.isNone() == b.isNone() && (a.id() >> 1) < (b.id() >> 1))));
85  // topologyId is set to 0 for None, so no harm im comparing them even if
86  // isNone()==true
87  }
88  };
89 
91 
96  template<class ctype>
97  struct GridViewInfo :
98  public std::map<GeometryType, EntityInfo<ctype>, GridViewInfoGTCompare>
99  {
101  std::string gridName;
103  std::string gridViewName;
105 
109  std::string partitionName;
110 
112 
125  void print(std::ostream &stream, std::string prefix) const {
126  if(!gridName.empty()) {
127  stream << prefix << gridName << ":\n";
128  prefix += " ";
129  }
130  if(!gridViewName.empty()) {
131  stream << prefix << gridViewName << ":\n";
132  prefix += " ";
133  }
134  if(!partitionName.empty()) {
135  stream << prefix << partitionName << ":\n";
136  prefix += " ";
137  }
138 
139  typedef typename GridViewInfo::const_iterator Iterator;
140  std::size_t dim = ~0;
141  const Iterator &end = this->end();
142  for(Iterator it = this->begin(); it != end; ++it) {
143  if(it->first.dim() != dim) {
144  dim = it->first.dim();
145  stream << prefix << "Dim = " << dim << ":\n";
146  }
147  stream << prefix << " " << it->first << ": Count = "
148  << it->second.count << ", Volume range = "
149  << "(" << it->second.volumeMin << ".."
150  << it->second.volumeMax << "), Total volume = "
151  << it->second.volumeSum << "\n";
152  }
153  }
154  };
155 
157 
162  template<class ctype>
163  std::ostream &operator<<(std::ostream &stream,
164  const GridViewInfo<ctype> &info)
165  {
166  info.print(stream, "");
167  return stream;
168  }
169 
170 #ifndef DOXYGEN
171  template<int codim>
173  struct FillGridInfoOperation {
174  template<class Entity, class Mapper, class Visited, class RefElem>
175  static void apply(const Entity &e, const Mapper &mapper, Visited &visited,
176  const typename Entity::Geometry &geo,
177  RefElem refelem,
179  {
180  typedef typename Entity::Geometry::ctype ctype;
181  static const std::size_t dimw = Entity::Geometry::coorddimension;
182  static const std::size_t dim = Entity::dimension;
183  std::vector<FieldVector<ctype, dimw> > coords;
184  for(int i = 0; i < refelem.size(codim); ++i) {
185  int index = mapper.map(e, i, codim);
186  if(visited[index])
187  continue;
188  visited[index] = true;
189 
190  GeometryType gt = refelem.type(i, codim);
191  coords.clear();
192  coords.resize( refelem.size(i, codim, dim) );
193  for(std::size_t corner = 0; corner < coords.size(); ++corner)
194  coords[ corner ] = geo.corner( refelem.subEntity( i, codim, corner, dim ) );
195  MultiLinearGeometry<ctype, dim-codim, dimw> mygeo(gt, coords);
196 
197  ctype volume = mygeo.volume();
198  EntityInfo<ctype> &ei = gridViewInfo[mygeo.type()];
199  ei.volumeMin = std::min(ei.volumeMin, volume);
200  ei.volumeMax = std::max(ei.volumeMax, volume);
201  ei.volumeSum += volume;
202  }
203  }
204  };
205 
206  template<int dimgrid>
207  struct MCMGNonElementLayout {
208  bool contains(GeometryType gt) const { return gt.dim() < dimgrid; }
209  };
210 #endif // !DOXYGEN
211 
213 
217  template<class GV>
218  void fillGridViewInfoSerial(const GV &gv,
219  GridViewInfo<typename GV::ctype> &gridViewInfo)
220  {
221  typedef typename GV::ctype ctype;
222  static const std::size_t dim = GV::dimension;
223  typedef typename GV::template Codim<0>::Iterator EIterator;
224  typedef typename GV::template Codim<0>::Geometry EGeometry;
225  typedef typename GV::IndexSet IndexSet;
226 
227  typedef typename GridViewInfo<ctype>::iterator InfoIterator;
228 
229  typedef ReferenceElements<ctype, dim> RefElems;
230 
232  std::vector<bool> visited(mapper.size(), false);
233 
234  gridViewInfo.gridName = className<typename GV::Grid>();
235  gridViewInfo.gridViewName = className<GV>();
236  gridViewInfo.partitionName = "";
237  gridViewInfo.clear();
238 
239  const EIterator &eend = gv.template end<0>();
240  for(EIterator eit = gv.template begin<0>(); eit != eend; ++eit) {
241  ctype volume = eit->geometry().volume();
242  EntityInfo<ctype> &ei = gridViewInfo[eit->type()];
243  ei.volumeMin = std::min(ei.volumeMin, volume);
244  ei.volumeMax = std::max(ei.volumeMax, volume);
245  ei.volumeSum += volume;
246 
247  if(!eit->type().isNone()) {
248  const EGeometry &geo = eit->geometry();
249  Hybrid::forEach(Std::make_index_sequence< dim >{},
250  [ & ](auto i){ FillGridInfoOperation< i+1 >::apply(*eit, mapper, visited, geo, RefElems::general(eit->type()), gridViewInfo); } );
251  }
252  }
253 
254  GeometryType gt;
255  gt.makeNone(dim);
256  if(gridViewInfo.count(gt) > 0) {
257  for(std::size_t codim = 0; codim < dim; ++codim) {
258  gt.makeNone(dim-codim);
259  EntityInfo<ctype> & ei = gridViewInfo[gt];
260  ei.volumeMin = ei.volumeMax = ei.volumeSum =
261  std::numeric_limits<ctype>::quiet_NaN();
262  }
263  gt.makeNone(0);
264  EntityInfo<ctype> & ei = gridViewInfo[gt];
265  ei.volumeMin = ei.volumeMax = ei.volumeSum = 0;
266  }
267 
268  const InfoIterator &end = gridViewInfo.end();
269  const IndexSet &is = gv.indexSet();
270  for(InfoIterator it = gridViewInfo.begin(); it != end; ++it) {
271  it->second.count = is.size(it->first);
272  if(it->second.count == 0)
273  DUNE_THROW(Exception, "Found Entities of geomentry type " <<
274  it->first << " while iterating through the grid, but "
275  "indexSet.size() == 0 for that geometry type");
276  }
277 
278  }
279 
280 } // namespace Dune
281 
282 
283 #endif // DUNE_GRID_UTILITY_GRIDINFO_HH
Dune::GridViewInfoGTCompare::operator()
bool operator()(const GeometryType &a, const GeometryType &b) const
compare two GeometryTypes
Definition: utility/gridinfo.hh:80
Dune::EntityInfo::EntityInfo
EntityInfo()
initialize the structure
Definition: utility/gridinfo.hh:63
Dune::GridViewInfo::gridName
std::string gridName
name of the grid class this information was extracted from
Definition: utility/gridinfo.hh:101
Dune::GridViewInfo::partitionName
std::string partitionName
name of the partition this information was extracted from
Definition: utility/gridinfo.hh:109
Dune::GridViewInfoGTCompare
Comparison object to sort GeometryType by majorly dimension.
Definition: utility/gridinfo.hh:76
Dune::GridViewInfo::print
void print(std::ostream &stream, std::string prefix) const
print the information contained in this object
Definition: utility/gridinfo.hh:125
Dune::Mapper
Mapper interface.
Definition: mapper.hh:107
Dune::GridViewInfo::gridViewName
std::string gridViewName
name of the class of the GridView this information was extracted from
Definition: utility/gridinfo.hh:103
Dune::Entity::Geometry
GridImp::template Codim< cd >::Geometry Geometry
The corresponding geometry type.
Definition: common/entity.hh:98
Dune::fillGridViewInfoSerial
void fillGridViewInfoSerial(const GV &gv, GridViewInfo< typename GV::ctype > &gridViewInfo)
fill a GridViewInfo structure from a serial grid
Definition: utility/gridinfo.hh:218
Dune::MultipleCodimMultipleGeomTypeMapper::size
size_type size() const
Return total number of entities in the entity set managed by the mapper.
Definition: mcmgmapper.hh:287
Dune::Entity::dimension
@ dimension
Know the grid dimension.
Definition: common/entity.hh:109
Dune::GridViewInfo::operator<<
std::ostream & operator<<(std::ostream &stream, const GridViewInfo< ctype > &info)
write a GridViewInfo object
Definition: utility/gridinfo.hh:163
Dune::EntityInfo::volumeMin
ctype volumeMin
minimum volume of all entities in the set.
Definition: utility/gridinfo.hh:41
Dune::EntityInfo::count
std::size_t count
number of entities in the set
Definition: utility/gridinfo.hh:34
Dune::GridViewInfo
structure to hold information about a certain GridView.
Definition: utility/gridinfo.hh:97
Dune::EntityInfo
Structure to hold statistical information about one type of entity.
Definition: utility/gridinfo.hh:32
mcmgmapper.hh
Mapper for multiple codim and multiple geometry types.
Dune::EntityInfo::volumeMax
ctype volumeMax
maximum volume of all entities in the set.
Definition: utility/gridinfo.hh:48
Dune::IndexSet::size
IndexType size(GeometryType type) const
Return total number of entities of given geometry type in entity set .
Definition: indexidset.hh:220
Dune::IndexSet
Index Set Interface base class.
Definition: common/grid.hh:346
Dune::MultipleCodimMultipleGeomTypeMapper
Implementation class for a multiple codim and multiple geometry type mapper.
Definition: mcmgmapper.hh:198
Dune
Include standard header files.
Definition: agrid.hh:58
Dune::VTK::GeometryType
GeometryType
Type representing VTK's entity geometry types.
Definition: common.hh:180
Dune::Entity
Wrapper class for entities.
Definition: common/entity.hh:63
Dune::EntityInfo::volumeSum
ctype volumeSum
sum of volumes of all entities in the set.
Definition: utility/gridinfo.hh:55