My Project
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#include <functional>
25
26#include <opm/grid/GridHelpers.hpp>
27#include <opm/grid/utility/OpmParserIncludes.hpp>
28
29
30#include <opm/grid/utility/platform_dependent/disable_warnings.h>
31
32#include <opm/grid/CpGrid.hpp>
33#include <dune/common/iteratorfacades.hh>
34
35#include <opm/grid/utility/platform_dependent/reenable_warnings.h>
36
37namespace Dune
38{
39namespace cpgrid
40{
43{
44public:
48 FaceCellsProxy(const Dune::CpGrid* grid, int cell_index)
49 : grid_(grid), cell_index_(cell_index)
50 {}
52 int operator[](int local_index)
53 {
54 return grid_->faceCell(cell_index_, local_index);
55 }
56private:
57 const Dune::CpGrid* grid_;
58 int cell_index_;
59};
60
64{
65public:
67
71 : grid_(grid)
72 {}
75 FaceCellsProxy operator[](int cell_index) const
76 {
77 return FaceCellsProxy(grid_, cell_index);
78 }
84 int operator()(int cell_index, int local_index) const
85 {
86 return grid_->faceCell(cell_index, local_index);
87 }
88private:
89 const Dune::CpGrid* grid_;
90};
91
92
94 {
95 public:
96 explicit IndexIterator(int index)
97 : index_(index)
98 {}
99
100 void increment()
101 {
102 ++index_;
103 }
104 void decrement()
105 {
106 --index_;
107 }
108 void advance(int n)
109 {
110 index_+=n;
111 }
112 int distanceTo(const IndexIterator& o)const
113 {
114 return o.index_-index_;
115 }
116 bool equals(const IndexIterator& o) const
117 {
118 return index_==o.index_;
119 }
120 protected:
121 int index_;
122 };
123
124
130template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
131 int (Dune::CpGrid::*SizeMethod)(int)const>
133{
134public:
136 : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
137 public IndexIterator
138 {
139 public:
140 iterator(const Dune::CpGrid* grid, int outer_index, int inner_index)
141 : IndexIterator(inner_index), grid_(grid), outer_index_(outer_index)
142 {}
143 int dereference() const
144 {
145 return std::mem_fn(AccessMethod)(*grid_, outer_index_, this->index_);
146 }
147 int elementAt(int n) const
148 {
149 return std::mem_fn(AccessMethod)(*grid_, outer_index_, n);
150 }
151 private:
152 const Dune::CpGrid* grid_;
153 int outer_index_;
154 };
155
156 typedef iterator const_iterator;
157
161 LocalIndexProxy(const Dune::CpGrid* grid, int cell_index)
162 : grid_(grid), cell_index_(cell_index)
163 {}
165 int operator[](int local_index)
166 {
167 return std::mem_fn(AccessMethod)(*grid_, cell_index_, local_index);
168 }
169 const_iterator begin()
170 {
171 return const_iterator(grid_, cell_index_, 0);
172 }
173 const_iterator end()
174 {
175 return const_iterator(grid_, cell_index_,
176 std::mem_fn(SizeMethod)(*grid_, cell_index_));
177 }
178private:
179 const Dune::CpGrid* grid_;
180 int cell_index_;
181};
182
188template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
189 int (Dune::CpGrid::*SizeMethod)(int)const>
191{
192public:
197 : grid_(grid)
198 {}
201 row_type operator[](int cell_index) const
202 {
203 return row_type(grid_, cell_index);
204 }
210 int operator()(int cell_index, int local_index) const
211 {
212 return std::mem_fn(AccessMethod)(*grid_, cell_index, local_index);
213 }
214private:
215 const Dune::CpGrid* grid_;
216};
217
222{
223public:
227 : LocalIndexContainerProxy<&Dune::CpGrid::faceVertex, &Dune::CpGrid::numFaceVertices>(grid)
228 {}
229};
230
232{
233public:
235 : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
236 public IndexIterator
237 {
238 public:
240 int index, int cell_index)
241 : IndexIterator(index), row_(row), cell_index_(cell_index)
242 {}
243 int dereference() const
244 {
245 return row_[this->index_].index();
246 }
247 int elementAt(int n) const
248 {
249 return row_[n].index();
250 }
251 int getCellIndex()const
252 {
253 return cell_index_;
254 }
255
256 private:
257 // Note that row_ is a Opm::iterator_range returned by Opm::SparseTable.
258 // and stored in CellFacesRow. It is a temporary object that only lives
259 // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
260 // Therefore it needs to be copied. As it is rather light weight this
261 // should be fast anyway. A const reference would mean that the iterator
262 // is not assignable.
264 int cell_index_;
265 };
266
267 typedef iterator const_iterator;
268
270 const int cell_index)
271 : row_(row), cell_index_(cell_index)
272 {}
273
274 const_iterator begin() const
275 {
276 return const_iterator(row_, 0, cell_index_);
277 }
278
279 const_iterator end() const
280 {
281 return const_iterator(row_, row_.size(), cell_index_);
282 }
283
284private:
285 // Note that row_ is a Opm::iterator_range returned by Opm::SparseTable.
286 // and stored in CellFacesRow. It is a temporary object that only lives
287 // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
288 // Therefore it needs to be copied. As it is rather light weight this
289 // should be fast anyway. A const reference would mean that the row
290 // is not assignable.
292 const int cell_index_;
293};
294
295
297{
298public:
299 typedef Cell2FacesRow row_type;
300
301 explicit Cell2FacesContainer(const Dune::CpGrid* grid)
302 : grid_(grid)
303 {};
304
305 Cell2FacesRow operator[](int cell_index) const
306 {
307 auto& row=grid_->cellFaceRow(cell_index);
308 return Cell2FacesRow(row, cell_index);
309 }
310
312 std::size_t noEntries() const
313 {
314 return grid_->numCellFaces();
315 }
316private:
317 const Dune::CpGrid* grid_;
318};
319
320} // end namespace cpgrid
321} // end namespace Dune
322
323namespace Opm
324{
325namespace UgGridHelpers
326{
327template<>
329{
331};
333template<const Dune::FieldVector<double, 3>& (Dune::CpGrid::*Method)(int)const>
335 : public Dune::RandomAccessIteratorFacade<CpGridCentroidIterator<Method>, Dune::FieldVector<double, 3>,
336 const Dune::FieldVector<double, 3>&, int>
337{
338public:
342 CpGridCentroidIterator(const Dune::CpGrid& grid, int cell_index)
343 : grid_(&grid), cell_index_(cell_index)
344 {}
345
346 const Dune::FieldVector<double, 3>& dereference() const
347 {
348 return std::mem_fn(Method)(*grid_, cell_index_);
349 }
350 void increment()
351 {
352 ++cell_index_;
353 }
354 const Dune::FieldVector<double, 3>& elementAt(int n) const
355 {
356 return std::mem_fn(Method)(*grid_, n);
357 }
358 void advance(int n)
359 {
360 cell_index_+=n;
361 }
362 void decrement()
363 {
364 --cell_index_;
365 }
366 int distanceTo(const CpGridCentroidIterator& o) const
367 {
368 return o.cell_index_-cell_index_;
369 }
370 bool equals(const CpGridCentroidIterator& o) const
371 {
372 return o.grid_==grid_ && o.cell_index_==cell_index_;
373 }
374
375private:
376 const Dune::CpGrid* grid_;
377 int cell_index_;
378};
379
380template<>
382{
384 typedef const double* ValueType;
385};
386
387typedef Dune::FieldVector<double, 3> Vector;
388
390int numCells(const Dune::CpGrid& grid);
391
393int numFaces(const Dune::CpGrid& grid);
394
396int dimensions(const Dune::CpGrid& grid);
397
399int numCellFaces(const Dune::CpGrid& grid);
400
402const int* cartDims(const Dune::CpGrid& grid);
403
408const int* globalCell(const Dune::CpGrid&);
409
414std::vector<int> createACTNUM(const Dune::CpGrid& grid);
415
416#if HAVE_ECL_INPUT
419EclipseGrid createEclipseGrid(const Dune::CpGrid& grid, const EclipseGrid& inputGrid);
420#endif
421
423beginCellCentroids(const Dune::CpGrid& grid);
424
429double cellCentroidCoordinate(const Dune::CpGrid& grid, int cell_index,
430 int coordinate);
431
435const double* cellCentroid(const Dune::CpGrid& grid, int cell_index);
436
440double cellCenterDepth(const Dune::CpGrid& grid, int cell_index);
441
442
448Vector faceCenterEcl(const Dune::CpGrid& grid, int cell_index, int face_tag);
449
456Vector faceAreaNormalEcl(const Dune::CpGrid& grid, int face_index);
457
461double cellVolume(const Dune::CpGrid& grid, int cell_index);
462
465 : public Dune::RandomAccessIteratorFacade<CellVolumeIterator, double, double, int>
466{
467public:
471 CellVolumeIterator(const Dune::CpGrid& grid, int cell_index)
472 : grid_(&grid), cell_index_(cell_index)
473 {}
474
475 double dereference() const
476 {
477 return grid_->cellVolume(cell_index_);
478 }
479 void increment()
480 {
481 ++cell_index_;
482 }
483 double elementAt(int n) const
484 {
485 return grid_->cellVolume(n);
486 }
487 void advance(int n)
488 {
489 cell_index_+=n;
490 }
491 void decrement()
492 {
493 --cell_index_;
494 }
495 int distanceTo(const CellVolumeIterator& o) const
496 {
497 return o.cell_index_-cell_index_;
498 }
499 bool equals(const CellVolumeIterator& o) const
500 {
501 return o.grid_==grid_ && o.cell_index_==cell_index_;
502 }
503
504private:
505 const Dune::CpGrid* grid_;
506 int cell_index_;
507};
508
509template<>
511{
513};
514
516CellVolumeIterator beginCellVolumes(const Dune::CpGrid& grid);
517
519CellVolumeIterator endCellVolumes(const Dune::CpGrid& grid);
520
521template<>
523{
525 typedef const Dune::CpGrid::Vector ValueType;
526};
527
530beginFaceCentroids(const Dune::CpGrid& grid);
531
537faceCentroid(const Dune::CpGrid& grid, int face_index);
538
539template<>
540struct FaceCellTraits<Dune::CpGrid>
541{
543};
545Dune::cpgrid::Cell2FacesContainer cell2Faces(const Dune::CpGrid& grid);
546
549faceCells(const Dune::CpGrid& grid);
550
551template<>
553{
555};
556
559face2Vertices(const Dune::CpGrid& grid);
560
564const double* vertexCoordinates(const Dune::CpGrid& grid, int index);
565
566const double* faceNormal(const Dune::CpGrid& grid, int face_index);
567
568double faceArea(const Dune::CpGrid& grid, int face_index);
569
574int faceTag(const Dune::CpGrid& grid,
576
577
578} // end namespace UgGridHelpers
579
580} // end namespace Opm
581
582#endif
[ provides Dune::Grid ]
Definition: CpGrid.hpp:210
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition: CpGrid.hpp:1011
Definition: GridHelpers.hpp:297
std::size_t noEntries() const
Get the number of non-zero entries.
Definition: GridHelpers.hpp:312
Definition: GridHelpers.hpp:237
Definition: GridHelpers.hpp:232
A class representing the face to cells mapping similar to the way done in UnstructuredGrid.
Definition: GridHelpers.hpp:64
FaceCellsContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:70
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition: GridHelpers.hpp:84
FaceCellsProxy operator[](int cell_index) const
Get the mapping for a cell.
Definition: GridHelpers.hpp:75
A proxy class representing a row of FaceCellsContainer.
Definition: GridHelpers.hpp:43
FaceCellsProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition: GridHelpers.hpp:48
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition: GridHelpers.hpp:52
A class representing the face to vertices mapping similar to the way done in UnstructuredGrid.
Definition: GridHelpers.hpp:222
FaceVerticesContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:226
Definition: GridHelpers.hpp:94
A class representing the sparse mapping of entity relations (e.g.
Definition: GridHelpers.hpp:191
LocalIndexContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition: GridHelpers.hpp:196
row_type operator[](int cell_index) const
Get the mapping for a cell.
Definition: GridHelpers.hpp:201
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition: GridHelpers.hpp:210
Definition: GridHelpers.hpp:138
A proxy class representing a row of LocalIndexContainerProxy.
Definition: GridHelpers.hpp:133
LocalIndexProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition: GridHelpers.hpp:161
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition: GridHelpers.hpp:165
A class used as a row type for OrientedEntityTable.
Definition: OrientedEntityTable.hpp:55
An iterator over the cell volumes.
Definition: GridHelpers.hpp:466
CellVolumeIterator(const Dune::CpGrid &grid, int cell_index)
Creates an iterator.
Definition: GridHelpers.hpp:471
An iterator over the cell volumes.
Definition: GridHelpers.hpp:337
CpGridCentroidIterator(const Dune::CpGrid &grid, int cell_index)
Creates an iterator.
Definition: GridHelpers.hpp:342
Copyright 2019 Equinor AS.
Definition: CartesianIndexMapper.hpp:10
Holds the implementation of the CpGrid as a pimple.
Definition: CellQuadrature.hpp:29
face_tag
Connection taxonomy.
Definition: preprocess.h:66
Traits of the cell centroids of a grid.
Definition: GridHelpers.hpp:126
The mapping of the grid type to type of the iterator over the cell volumes.
Definition: GridHelpers.hpp:191
Maps the grid type to the associated type of the face to vertices mapping.
Definition: GridHelpers.hpp:292
Traits of the face to attached cell mappping of a grid.
Definition: GridHelpers.hpp:334
Traits of the face centroids of a grid.
Definition: GridHelpers.hpp:241