All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Groups Pages
GridHelpers.hpp
1 /*
2  Copyright 2014, 2015 Dr. Markus Blatt - HPC-Simulation-Software & Services.
3  Copyright 2014 Statoil AS
4  Copyright 2015 NTNU
5 
6  This file is part of the Open Porous Media project (OPM).
7 
8  OPM is free software: you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation, either version 3 of the License, or
11  (at your option) any later version.
12 
13  OPM is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with OPM. If not, see <http://www.gnu.org/licenses/>.
20 */
21 #ifndef DUNE_CORNERPOINT_GRIDHELPERS_HEADER_INCLUDED
22 #define DUNE_CORNERPOINT_GRIDHELPERS_HEADER_INCLUDED
23 
24 
25 #include <opm/core/grid/GridHelpers.hpp>
26 #include <opm/grid/utility/OpmParserIncludes.hpp>
27 
28 
29 #include <opm/grid/utility/platform_dependent/disable_warnings.h>
30 
31 #include <dune/grid/CpGrid.hpp>
32 #include <dune/common/iteratorfacades.hh>
33 
34 #include <opm/grid/utility/platform_dependent/reenable_warnings.h>
35 
36 namespace Dune
37 {
38 namespace cpgrid
39 {
42 {
43 public:
47  FaceCellsProxy(const Dune::CpGrid* grid, int cell_index)
48  : grid_(grid), cell_index_(cell_index)
49  {}
51  int operator[](int local_index)
52  {
53  return grid_->faceCell(cell_index_, local_index);
54  }
55 private:
56  const Dune::CpGrid* grid_;
57  int cell_index_;
58 };
59 
63 {
64 public:
65  typedef FaceCellsProxy row_type;
66 
69  explicit FaceCellsContainerProxy(const Dune::CpGrid* grid)
70  : grid_(grid)
71  {}
74  FaceCellsProxy operator[](int cell_index) const
75  {
76  return FaceCellsProxy(grid_, cell_index);
77  }
83  int operator()(int cell_index, int local_index) const
84  {
85  return grid_->faceCell(cell_index, local_index);
86  }
87 private:
88  const Dune::CpGrid* grid_;
89 };
90 
91 
93  {
94  public:
95  explicit IndexIterator(int index)
96  : index_(index)
97  {}
98 
99  void increment()
100  {
101  ++index_;
102  }
103  void decrement()
104  {
105  --index_;
106  }
107  void advance(int n)
108  {
109  index_+=n;
110  }
111  int distanceTo(const IndexIterator& o)const
112  {
113  return o.index_-index_;
114  }
115  bool equals(const IndexIterator& o) const
116  {
117  return index_==o.index_;
118  }
119  protected:
120  int index_;
121  };
122 
123 
129 template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
130  int (Dune::CpGrid::*SizeMethod)(int)const>
132 {
133 public:
134  class iterator
135  : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
136  public IndexIterator
137  {
138  public:
139  iterator(const Dune::CpGrid* grid, int outer_index, int inner_index)
140  : IndexIterator(inner_index), grid_(grid), outer_index_(outer_index)
141  {}
142  int dereference() const
143  {
144  return std::mem_fn(AccessMethod)(*grid_, outer_index_, this->index_);
145  }
146  int elementAt(int n) const
147  {
148  return std::mem_fn(AccessMethod)(*grid_, outer_index_, n);
149  }
150  private:
151  const Dune::CpGrid* grid_;
152  int outer_index_;
153  };
154 
155  typedef iterator const_iterator;
156 
160  LocalIndexProxy(const Dune::CpGrid* grid, int cell_index)
161  : grid_(grid), cell_index_(cell_index)
162  {}
164  int operator[](int local_index)
165  {
166  return std::mem_fn(AccessMethod)(*grid_, cell_index_, local_index);
167  }
168  const_iterator begin()
169  {
170  return const_iterator(grid_, cell_index_, 0);
171  }
172  const_iterator end()
173  {
174  return const_iterator(grid_, cell_index_,
175  std::mem_fn(SizeMethod)(*grid_, cell_index_));
176  }
177 private:
178  const Dune::CpGrid* grid_;
179  int cell_index_;
180 };
181 
187 template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
188  int (Dune::CpGrid::*SizeMethod)(int)const>
190 {
191 public:
196  : grid_(grid)
197  {}
200  row_type operator[](int cell_index) const
201  {
202  return row_type(grid_, cell_index);
203  }
209  int operator()(int cell_index, int local_index) const
210  {
211  return std::mem_fn(AccessMethod)(*grid_, cell_index, local_index);
212  }
213 private:
214  const Dune::CpGrid* grid_;
215 };
216 
221 {
222 public:
226  : LocalIndexContainerProxy<&Dune::CpGrid::faceVertex, &Dune::CpGrid::numFaceVertices>(grid)
227  {}
228 };
229 
231 {
232 public:
233  class iterator
234  : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
235  public IndexIterator
236  {
237  public:
239  int index, int cell_index)
240  : IndexIterator(index), row_(row), cell_index_(cell_index)
241  {}
242  int dereference() const
243  {
244  return row_[this->index_].index();
245  }
246  int elementAt(int n) const
247  {
248  return row_[n].index();
249  }
250  int getCellIndex()const
251  {
252  return cell_index_;
253  }
254 
255  private:
256  // Note that row_ is a boost::iterator_range returned by Opm::SparseTable.
257  // and stored in CellFacesRow. It is a temporary object that only lives
258  // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
259  // Therefore it needs to be copied. As it is rather light weight this
260  // should be fast anyway. A const reference would mean that the iterator
261  // is not assignable.
263  int cell_index_;
264  };
265 
266  typedef iterator const_iterator;
267 
269  const int cell_index)
270  : row_(row), cell_index_(cell_index)
271  {}
272 
273  const_iterator begin() const
274  {
275  return const_iterator(row_, 0, cell_index_);
276  }
277 
278  const_iterator end() const
279  {
280  return const_iterator(row_, row_.size(), cell_index_);
281  }
282 
283 private:
284  // Note that row_ is a boost::iterator_range returned by Opm::SparseTable.
285  // and stored in CellFacesRow. It is a temporary object that only lives
286  // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
287  // Therefore it needs to be copied. As it is rather light weight this
288  // should be fast anyway. A const reference would mean that the row
289  // is not assignable.
291  const int cell_index_;
292 };
293 
294 
296 {
297 public:
298  typedef Cell2FacesRow row_type;
299 
300  explicit Cell2FacesContainer(const Dune::CpGrid* grid)
301  : grid_(grid)
302  {};
303 
304  Cell2FacesRow operator[](int cell_index) const
305  {
306  auto& row=grid_->cellFaceRow(cell_index);
307  return Cell2FacesRow(row, cell_index);
308  }
309 
311  std::size_t noEntries() const
312  {
313  return grid_->numCellFaces();
314  }
315 private:
316  const Dune::CpGrid* grid_;
317 };
318 
319 } // end namespace cpgrid
320 } // end namespace Dune
321 
322 namespace Opm
323 {
324 namespace UgGridHelpers
325 {
326 template<>
328 {
330 };
332 template<const Dune::FieldVector<double, 3>& (Dune::CpGrid::*Method)(int)const>
334  : public Dune::RandomAccessIteratorFacade<CpGridCentroidIterator<Method>, Dune::FieldVector<double, 3>,
335  const Dune::FieldVector<double, 3>&, int>
336 {
337 public:
341  CpGridCentroidIterator(const Dune::CpGrid& grid, int cell_index)
342  : grid_(&grid), cell_index_(cell_index)
343  {}
344 
345  const Dune::FieldVector<double, 3>& dereference() const
346  {
347  return std::mem_fn(Method)(*grid_, cell_index_);
348  }
349  void increment()
350  {
351  ++cell_index_;
352  }
353  const Dune::FieldVector<double, 3>& elementAt(int n) const
354  {
355  return std::mem_fn(Method)(*grid_, n);
356  }
357  void advance(int n)
358  {
359  cell_index_+=n;
360  }
361  void decrement()
362  {
363  --cell_index_;
364  }
365  int distanceTo(const CpGridCentroidIterator& o) const
366  {
367  return o.cell_index_-cell_index_;
368  }
369  bool equals(const CpGridCentroidIterator& o) const
370  {
371  return o.grid_==grid_ && o.cell_index_==cell_index_;
372  }
373 
374 private:
375  const Dune::CpGrid* grid_;
376  int cell_index_;
377 };
378 
379 template<>
380 struct CellCentroidTraits<Dune::CpGrid>
381 {
383  typedef const double* ValueType;
384 };
385 
386 typedef Dune::FieldVector<double, 3> Vector;
387 
389 int numCells(const Dune::CpGrid& grid);
390 
392 int numFaces(const Dune::CpGrid& grid);
393 
395 int dimensions(const Dune::CpGrid& grid);
396 
398 int numCellFaces(const Dune::CpGrid& grid);
399 
401 const int* cartDims(const Dune::CpGrid& grid);
402 
407 const int* globalCell(const Dune::CpGrid&);
408 
409 #if HAVE_OPM_PARSER
410 EclipseGrid createEclipseGrid(const Dune::CpGrid& grid, const EclipseGrid& inputGrid);
413 #endif
414 
416 beginCellCentroids(const Dune::CpGrid& grid);
417 
422 double cellCentroidCoordinate(const Dune::CpGrid& grid, int cell_index,
423  int coordinate);
424 
428 const double* cellCentroid(const Dune::CpGrid& grid, int cell_index);
429 
433 double cellCenterDepth(const Dune::CpGrid& grid, int cell_index);
434 
435 
441 Vector faceCenterEcl(const Dune::CpGrid& grid, int cell_index, int face_tag);
442 
449 Vector faceAreaNormalEcl(const Dune::CpGrid& grid, int face_index);
450 
454 double cellVolume(const Dune::CpGrid& grid, int cell_index);
455 
458  : public Dune::RandomAccessIteratorFacade<CellVolumeIterator, double, double, int>
459 {
460 public:
464  CellVolumeIterator(const Dune::CpGrid& grid, int cell_index)
465  : grid_(&grid), cell_index_(cell_index)
466  {}
467 
468  double dereference() const
469  {
470  return grid_->cellVolume(cell_index_);
471  }
472  void increment()
473  {
474  ++cell_index_;
475  }
476  double elementAt(int n) const
477  {
478  return grid_->cellVolume(n);
479  }
480  void advance(int n)
481  {
482  cell_index_+=n;
483  }
484  void decrement()
485  {
486  --cell_index_;
487  }
488  int distanceTo(const CellVolumeIterator& o) const
489  {
490  return o.cell_index_-cell_index_;
491  }
492  bool equals(const CellVolumeIterator& o) const
493  {
494  return o.grid_==grid_ && o.cell_index_==cell_index_;
495  }
496 
497 private:
498  const Dune::CpGrid* grid_;
499  int cell_index_;
500 };
501 
502 template<>
503 struct CellVolumeIteratorTraits<Dune::CpGrid>
504 {
506 };
507 
509 CellVolumeIterator beginCellVolumes(const Dune::CpGrid& grid);
510 
512 CellVolumeIterator endCellVolumes(const Dune::CpGrid& grid);
513 
514 template<>
515 struct FaceCentroidTraits<Dune::CpGrid>
516 {
518  typedef const Dune::CpGrid::Vector ValueType;
519 };
520 
523 beginFaceCentroids(const Dune::CpGrid& grid);
524 
530 faceCentroid(const Dune::CpGrid& grid, int face_index);
531 
532 template<>
533 struct FaceCellTraits<Dune::CpGrid>
534 {
536 };
538 Dune::cpgrid::Cell2FacesContainer cell2Faces(const Dune::CpGrid& grid);
539 
542 faceCells(const Dune::CpGrid& grid);
543 
544 template<>
545 struct Face2VerticesTraits<Dune::CpGrid>
546 {
548 };
549 
552 face2Vertices(const Dune::CpGrid& grid);
553 
557 const double* vertexCoordinates(const Dune::CpGrid& grid, int index);
558 
559 const double* faceNormal(const Dune::CpGrid& grid, int face_index);
560 
561 double faceArea(const Dune::CpGrid& grid, int face_index);
562 
567 int faceTag(const Dune::CpGrid& grid,
568  const Dune::cpgrid::Cell2FacesRow::iterator& cell_face);
569 
570 
571 } // end namespace UgGridHelpers
572 
573 } // end namespace Opm
574 
575 #endif
FaceCellsProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition: GridHelpers.hpp:47
CpGridCentroidIterator(const Dune::CpGrid &grid, int cell_index)
Creates an iterator.
Definition: GridHelpers.hpp:341
row_type operator[](int cell_index) const
Get the mapping for a cell.
Definition: GridHelpers.hpp:200
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.hpp:769
FaceVerticesContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:225
A proxy class representing a row of FaceCellsContainer.
Definition: GridHelpers.hpp:41
FaceCellsContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:69
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition: GridHelpers.hpp:209
A proxy class representing a row of LocalIndexContainerProxy.
Definition: GridHelpers.hpp:131
Definition: GridHelpers.hpp:233
[ provides Dune::Grid ]
Definition: CpGrid.hpp:213
CellVolumeIterator(const Dune::CpGrid &grid, int cell_index)
Creates an iterator.
Definition: GridHelpers.hpp:464
A class representing the sparse mapping of entity relations (e.g.
Definition: GridHelpers.hpp:189
An iterator over the cell volumes.
Definition: GridHelpers.hpp:457
Definition: GridHelpers.hpp:295
FaceCellsProxy operator[](int cell_index) const
Get the mapping for a cell.
Definition: GridHelpers.hpp:74
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition: GridHelpers.hpp:83
Definition: GridHelpers.hpp:230
Traits of the cell centroids of a grid.
Definition: GridHelpers.hpp:131
A class representing the face to cells mapping similar to the way done in UnstructuredGrid.
Definition: GridHelpers.hpp:62
LocalIndexProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition: GridHelpers.hpp:160
Traits of the face to attached cell mappping of a grid.
Definition: GridHelpers.hpp:345
std::size_t noEntries() const
Get the number of non-zero entries.
Definition: GridHelpers.hpp:311
A class representing the face to vertices mapping similar to the way done in UnstructuredGrid.
Definition: GridHelpers.hpp:219
An iterator over the cell volumes.
Definition: GridHelpers.hpp:333
Maps the grid type to the associated type of the face to vertices mapping.
Definition: GridHelpers.hpp:303
face_tag
Connection taxonomy.
Definition: preprocess.h:67
Definition: GridHelpers.hpp:92
Traits of the face centroids of a grid.
Definition: GridHelpers.hpp:246
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition: GridHelpers.hpp:51
LocalIndexContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:195
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
Definition: GridHelpers.hpp:134
The mapping of the grid type to type of the iterator over the cell volumes.
Definition: GridHelpers.hpp:196
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition: GridHelpers.hpp:164