My Project
Loading...
Searching...
No Matches
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
28#include <opm/grid/utility/platform_dependent/disable_warnings.h>
29
30#include <opm/grid/CpGrid.hpp>
31#include <dune/common/iteratorfacades.hh>
32
33#include <opm/grid/utility/platform_dependent/reenable_warnings.h>
34
35namespace Dune
36{
37namespace cpgrid
38{
41{
42public:
46 FaceCellsProxy(const Dune::CpGrid* grid, int cell_index)
47 : grid_(grid), cell_index_(cell_index)
48 {}
50 int operator[](int local_index)
51 {
52 return grid_->faceCell(cell_index_, local_index);
53 }
54private:
55 const Dune::CpGrid* grid_;
56 int cell_index_;
57};
58
62{
63public:
65
69 : grid_(grid)
70 {}
73 FaceCellsProxy operator[](int cell_index) const
74 {
75 return FaceCellsProxy(grid_, cell_index);
76 }
82 int operator()(int cell_index, int local_index) const
83 {
84 return grid_->faceCell(cell_index, local_index);
85 }
86private:
87 const Dune::CpGrid* grid_;
88};
89
90
92 {
93 public:
94 explicit IndexIterator(int index)
95 : index_(index)
96 {}
97
98 void increment()
99 {
100 ++index_;
101 }
102 void decrement()
103 {
104 --index_;
105 }
106 void advance(int n)
107 {
108 index_+=n;
109 }
110 int distanceTo(const IndexIterator& o)const
111 {
112 return o.index_-index_;
113 }
114 bool equals(const IndexIterator& o) const
115 {
116 return index_==o.index_;
117 }
118 protected:
119 int index_;
120 };
121
122
128template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
129 int (Dune::CpGrid::*SizeMethod)(int)const>
131{
132public:
134 : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
135 public IndexIterator
136 {
137 public:
138 iterator(const Dune::CpGrid* grid, int outer_index, int inner_index)
139 : IndexIterator(inner_index), grid_(grid), outer_index_(outer_index)
140 {}
141 int dereference() const
142 {
143 return std::mem_fn(AccessMethod)(*grid_, outer_index_, this->index_);
144 }
145 int elementAt(int n) const
146 {
147 return std::mem_fn(AccessMethod)(*grid_, outer_index_, n);
148 }
149 private:
150 const Dune::CpGrid* grid_;
151 int outer_index_;
152 };
153
154 typedef iterator const_iterator;
155
159 LocalIndexProxy(const Dune::CpGrid* grid, int cell_index)
160 : grid_(grid), cell_index_(cell_index)
161 {}
163 int operator[](int local_index)
164 {
165 return std::mem_fn(AccessMethod)(*grid_, cell_index_, local_index);
166 }
167 const_iterator begin()
168 {
169 return const_iterator(grid_, cell_index_, 0);
170 }
171 const_iterator end()
172 {
173 return const_iterator(grid_, cell_index_,
174 std::mem_fn(SizeMethod)(*grid_, cell_index_));
175 }
176private:
177 const Dune::CpGrid* grid_;
178 int cell_index_;
179};
180
186template<int (Dune::CpGrid::*AccessMethod)(int,int)const,
187 int (Dune::CpGrid::*SizeMethod)(int)const>
189{
190public:
195 : grid_(grid)
196 {}
199 row_type operator[](int cell_index) const
200 {
201 return row_type(grid_, cell_index);
202 }
208 int operator()(int cell_index, int local_index) const
209 {
210 return std::mem_fn(AccessMethod)(*grid_, cell_index, local_index);
211 }
212private:
213 const Dune::CpGrid* grid_;
214};
215
228
230{
231public:
233 : public Dune::RandomAccessIteratorFacade<iterator,int, int, int>,
234 public IndexIterator
235 {
236 public:
238 int index, int cell_index)
239 : IndexIterator(index), row_(row), cell_index_(cell_index)
240 {}
241 int dereference() const
242 {
243 return row_[this->index_].index();
244 }
245 int elementAt(int n) const
246 {
247 return row_[n].index();
248 }
249 int getCellIndex()const
250 {
251 return cell_index_;
252 }
253
254 private:
255 // Note that row_ is a Opm::iterator_range returned by Opm::SparseTable.
256 // and stored in CellFacesRow. It is a temporary object that only lives
257 // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
258 // Therefore it needs to be copied. As it is rather light weight this
259 // should be fast anyway. A const reference would mean that the iterator
260 // is not assignable.
262 int cell_index_;
263 };
264
265 typedef iterator const_iterator;
266
268 const int cell_index)
269 : row_(row), cell_index_(cell_index)
270 {}
271
272 const_iterator begin() const
273 {
274 return const_iterator(row_, 0, cell_index_);
275 }
276
277 const_iterator end() const
278 {
279 return const_iterator(row_, row_.size(), cell_index_);
280 }
281
282private:
283 // Note that row_ is a Opm::iterator_range returned by Opm::SparseTable.
284 // and stored in CellFacesRow. It is a temporary object that only lives
285 // as long as there is a reference to it e.g. by storing the Cell2FacesRow.
286 // Therefore it needs to be copied. As it is rather light weight this
287 // should be fast anyway. A const reference would mean that the row
288 // is not assignable.
290 const int cell_index_;
291};
292
293
295{
296public:
297 typedef Cell2FacesRow row_type;
298
299 explicit Cell2FacesContainer(const Dune::CpGrid* grid)
300 : grid_(grid)
301 {};
302
303 Cell2FacesRow operator[](int cell_index) const
304 {
305 auto& row=grid_->cellFaceRow(cell_index);
306 return Cell2FacesRow(row, cell_index);
307 }
308
310 std::size_t noEntries() const
311 {
312 return grid_->numCellFaces();
313 }
314private:
315 const Dune::CpGrid* grid_;
316};
317
318} // end namespace cpgrid
319} // end namespace Dune
320
321namespace Opm
322{
323namespace UgGridHelpers
324{
325template<>
326struct Cell2FacesTraits<Dune::CpGrid>
327{
329};
331template<const Dune::FieldVector<double, 3>& (Dune::CpGrid::*Method)(int)const>
332class CpGridCentroidIterator
333 : public Dune::RandomAccessIteratorFacade<CpGridCentroidIterator<Method>, Dune::FieldVector<double, 3>,
334 const Dune::FieldVector<double, 3>&, int>
335{
336public:
340 CpGridCentroidIterator(const Dune::CpGrid& grid, int cell_index)
341 : grid_(&grid), cell_index_(cell_index)
342 {}
343
344 const Dune::FieldVector<double, 3>& dereference() const
345 {
346 return std::mem_fn(Method)(*grid_, cell_index_);
347 }
348 void increment()
349 {
350 ++cell_index_;
351 }
352 const Dune::FieldVector<double, 3>& elementAt(int n) const
353 {
354 return std::mem_fn(Method)(*grid_, n);
355 }
356 void advance(int n)
357 {
358 cell_index_+=n;
359 }
360 void decrement()
361 {
362 --cell_index_;
363 }
364 int distanceTo(const CpGridCentroidIterator& o) const
365 {
366 return o.cell_index_-cell_index_;
367 }
368 bool equals(const CpGridCentroidIterator& o) const
369 {
370 return o.grid_==grid_ && o.cell_index_==cell_index_;
371 }
372
373private:
374 const Dune::CpGrid* grid_;
375 int cell_index_;
376};
377
378template<>
379struct CellCentroidTraits<Dune::CpGrid>
380{
381 typedef CpGridCentroidIterator<&Dune::CpGrid::cellCentroid> IteratorType;
382 typedef const double* ValueType;
383};
384
385typedef Dune::FieldVector<double, 3> Vector;
386
388int numCells(const Dune::CpGrid& grid);
389
391int numFaces(const Dune::CpGrid& grid);
392
394int dimensions(const Dune::CpGrid& grid);
395
397int numCellFaces(const Dune::CpGrid& grid);
398
400const int* cartDims(const Dune::CpGrid& grid);
401
406const int* globalCell(const Dune::CpGrid&);
407
408#if HAVE_ECL_INPUT
413std::vector<int> createACTNUM(const Dune::CpGrid& grid);
414
417EclipseGrid createEclipseGrid(const Dune::CpGrid& grid, const EclipseGrid& inputGrid);
418#endif
419
420CellCentroidTraits<Dune::CpGrid>::IteratorType
421beginCellCentroids(const Dune::CpGrid& grid);
422
427double cellCentroidCoordinate(const Dune::CpGrid& grid, int cell_index,
428 int coordinate);
429
433const double* cellCentroid(const Dune::CpGrid& grid, int cell_index);
434
438double cellCenterDepth(const Dune::CpGrid& grid, int cell_index);
439
440
446Vector faceCenterEcl(const Dune::CpGrid& grid, int cell_index, int face_tag);
447
454Vector faceAreaNormalEcl(const Dune::CpGrid& grid, int face_index);
455
459double cellVolume(const Dune::CpGrid& grid, int cell_index);
460
462class CellVolumeIterator
463 : public Dune::RandomAccessIteratorFacade<CellVolumeIterator, double, double, int>
464{
465public:
469 CellVolumeIterator(const Dune::CpGrid& grid, int cell_index)
470 : grid_(&grid), cell_index_(cell_index)
471 {}
472
473 double dereference() const
474 {
475 return grid_->cellVolume(cell_index_);
476 }
477 void increment()
478 {
479 ++cell_index_;
480 }
481 double elementAt(int n) const
482 {
483 return grid_->cellVolume(n);
484 }
485 void advance(int n)
486 {
487 cell_index_+=n;
488 }
489 void decrement()
490 {
491 --cell_index_;
492 }
493 int distanceTo(const CellVolumeIterator& o) const
494 {
495 return o.cell_index_-cell_index_;
496 }
497 bool equals(const CellVolumeIterator& o) const
498 {
499 return o.grid_==grid_ && o.cell_index_==cell_index_;
500 }
501
502private:
503 const Dune::CpGrid* grid_;
504 int cell_index_;
505};
506
507template<>
508struct CellVolumeIteratorTraits<Dune::CpGrid>
509{
510 typedef CellVolumeIterator IteratorType;
511};
512
514CellVolumeIterator beginCellVolumes(const Dune::CpGrid& grid);
515
517CellVolumeIterator endCellVolumes(const Dune::CpGrid& grid);
518
519template<>
520struct FaceCentroidTraits<Dune::CpGrid>
521{
522 typedef CpGridCentroidIterator<&Dune::CpGrid::faceCentroid> IteratorType;
523 typedef const Dune::CpGrid::Vector ValueType;
524};
525
527FaceCentroidTraits<Dune::CpGrid>::IteratorType
528beginFaceCentroids(const Dune::CpGrid& grid);
529
534const FaceCentroidTraits<Dune::CpGrid>::ValueType&
535faceCentroid(const Dune::CpGrid& grid, int face_index);
536
537template<>
538struct FaceCellTraits<Dune::CpGrid>
539{
541};
543Dune::cpgrid::Cell2FacesContainer cell2Faces(const Dune::CpGrid& grid);
544
546FaceCellTraits<Dune::CpGrid>::Type
547faceCells(const Dune::CpGrid& grid);
548
549template<>
550struct Face2VerticesTraits<Dune::CpGrid>
551{
553};
554
556Face2VerticesTraits<Dune::CpGrid>::Type
557face2Vertices(const Dune::CpGrid& grid);
558
562const double* vertexCoordinates(const Dune::CpGrid& grid, int index);
563
564const double* faceNormal(const Dune::CpGrid& grid, int face_index);
565
566double faceArea(const Dune::CpGrid& grid, int face_index);
567
572int faceTag(const Dune::CpGrid& grid,
574
575
576} // end namespace UgGridHelpers
577
578} // end namespace Opm
579
580#endif
[ provides Dune::Grid ]
Definition CpGrid.hpp:201
int faceCell(int face, int local_index) const
Get the index identifying a cell attached to a face.
Definition CpGrid.cpp:1113
Definition GridHelpers.hpp:295
std::size_t noEntries() const
Get the number of non-zero entries.
Definition GridHelpers.hpp:310
Definition GridHelpers.hpp:235
Definition GridHelpers.hpp:230
A class representing the face to cells mapping similar to the way done in UnstructuredGrid.
Definition GridHelpers.hpp:62
FaceCellsContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition GridHelpers.hpp:68
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition GridHelpers.hpp:82
FaceCellsProxy operator[](int cell_index) const
Get the mapping for a cell.
Definition GridHelpers.hpp:73
A proxy class representing a row of FaceCellsContainer.
Definition GridHelpers.hpp:41
FaceCellsProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition GridHelpers.hpp:46
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition GridHelpers.hpp:50
A class representing the face to vertices mapping similar to the way done in UnstructuredGrid.
Definition GridHelpers.hpp:220
FaceVerticesContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition GridHelpers.hpp:224
Definition GridHelpers.hpp:92
A class representing the sparse mapping of entity relations (e.g.
Definition GridHelpers.hpp:189
LocalIndexContainerProxy(const Dune::CpGrid *grid)
Constructor.
Definition GridHelpers.hpp:194
row_type operator[](int cell_index) const
Get the mapping for a cell.
Definition GridHelpers.hpp:199
int operator()(int cell_index, int local_index) const
Get a face associated with a cell.
Definition GridHelpers.hpp:208
Definition GridHelpers.hpp:136
A proxy class representing a row of LocalIndexContainerProxy.
Definition GridHelpers.hpp:131
LocalIndexProxy(const Dune::CpGrid *grid, int cell_index)
Constructor.
Definition GridHelpers.hpp:159
int operator[](int local_index)
Get the index of the cell associated with a local_index.
Definition GridHelpers.hpp:163
A class used as a row type for OrientedEntityTable.
Definition OrientedEntityTable.hpp:55
Copyright 2019 Equinor AS.
Definition CartesianIndexMapper.hpp:10
Holds the implementation of the CpGrid as a pimple.
Definition CellQuadrature.cpp:68
face_tag
Connection taxonomy.
Definition preprocess.h:66