My Project
GridHelpers.hpp
1 /*
2  Copyright 2014, 2015 Dr. Markus Blatt - HPC-Simulation-Software & Services
3  Copyright 2014 Statoil AS
4  Copyright 2015
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 OPM_CORE_GRIDHELPERS_HEADER_INCLUDED
22 #define OPM_CORE_GRIDHELPERS_HEADER_INCLUDED
23 
24 #include <dune/common/fvector.hh>
25 
27 
28 #include <opm/grid/utility/IteratorRange.hpp>
29 #include <opm/grid/utility/OpmParserIncludes.hpp>
30 
31 
32 namespace Opm
33 {
34 namespace UgGridHelpers
35 {
36 
43 {
44 public:
47 
53  SparseTableView(int* data, int *offset, std::size_t size_arg)
54  : data_(data), offset_(offset), size_(size_arg)
55  {}
56 
60  row_type operator[](std::size_t row) const
61  {
62  assert(row<=size());
63  return row_type{data_ + offset_[row],
64  data_ + offset_[row+1]};
65  }
66 
69  std::size_t size() const
70  {
71  return size_;
72  }
73 
75  std::size_t noEntries() const
76  {
77  return offset_[size_];
78  }
79 
80 private:
82  const int* data_;
86  const int* offset_;
88  std::size_t size_;
89 };
90 
92 int numCells(const UnstructuredGrid& grid);
93 
95 int numFaces(const UnstructuredGrid& grid);
96 
98 int dimensions(const UnstructuredGrid& grid);
99 
101 int numCellFaces(const UnstructuredGrid& grid);
102 
104 const int* cartDims(const UnstructuredGrid& grid);
105 
110 const int* globalCell(const UnstructuredGrid& grid);
111 
116 std::vector<int> createACTNUM(const UnstructuredGrid& grid);
117 
118 
124 template<class G>
126 {
127 };
128 
129 template<>
131 {
132  typedef const double* IteratorType;
133  typedef const double* ValueType;
134 };
135 
141 beginCellCentroids(const UnstructuredGrid& grid);
142 
143 
147 double cellCenterDepth(const UnstructuredGrid& grid, int cell_index);
148 
154 Dune::FieldVector<double,3> faceCenterEcl(const UnstructuredGrid& grid, int cell_index, int face_tag);
155 
162 Dune::FieldVector<double,3> faceAreaNormalEcl(const UnstructuredGrid& grid, int face_index);
163 
164 
169 double cellCentroidCoordinate(const UnstructuredGrid& grid, int cell_index,
170  int coordinate);
171 
172 
176 const double* cellCentroid(const UnstructuredGrid& grid, int cell_index);
177 
178 
182 double cellVolume(const UnstructuredGrid& grid, int cell_index);
183 
189 template<class T>
191 {
192 };
193 
194 template<>
196 {
197  typedef const double* IteratorType;
198 };
199 
220 #if HAVE_ECL_INPUT
223 Opm::EclipseGrid createEclipseGrid(const UnstructuredGrid& grid, const Opm::EclipseGrid& inputGrid );
224 #endif
225 
227 const double* beginCellVolumes(const UnstructuredGrid& grid);
228 
230 const double* endCellVolumes(const UnstructuredGrid& grid);
231 
232 
233 
239 template<class G>
241 {
242 };
243 
244 template<>
246 {
247  typedef const double* IteratorType;
248  typedef const double* ValueType;
249 };
250 
253 beginFaceCentroids(const UnstructuredGrid& grid);
254 
259 faceCentroid(const UnstructuredGrid& grid, int face_index);
260 
264 const double* faceNormal(const UnstructuredGrid& grid, int face_index);
265 
269 double faceArea(const UnstructuredGrid& grid, int face_index);
270 
275 template<class T>
277 {
278 };
279 
280 template<>
282 {
283  typedef SparseTableView Type;
284 };
285 
290 template<class T>
292 {
293 };
294 
295 template<>
297 {
298  typedef SparseTableView Type;
299 };
300 
303 cell2Faces(const UnstructuredGrid& grid);
304 
307 face2Vertices(const UnstructuredGrid& grid);
308 
312 const double* vertexCoordinates(const UnstructuredGrid& grid, int index);
313 
315 {
316 public:
317  FaceCellsProxy(const UnstructuredGrid& grid)
318  : face_cells_(grid.face_cells)
319  {}
320  int operator()(int face_index, int local_index) const
321  {
322  return face_cells_[2*face_index+local_index];
323  }
324 private:
325  const int* face_cells_;
326 };
327 
332 template<class T>
334 {};
335 
336 template<>
338 {
339  typedef FaceCellsProxy Type;
340 };
341 
344 
350 template<class T>
351 T* increment(T* cc, int i, int dim)
352 {
353  return cc+(i*dim);
354 }
359 template<class T>
360 T increment(const T& t, int i, int)
361 {
362  return t+i;
363 }
364 
369 template<class T>
370 double getCoordinate(T* cc, int i)
371 {
372  return cc[i];
373 }
374 
380 template<class T>
381 double getCoordinate(T t, int i)
382 {
383  return (*t)[i];
384 }
385 
386 
387 
388 
389 
390 } // end namespace UGGridHelpers
391 } // end namespace OPM
392 #endif
Main OPM-Core grid data structure along with helper functions for construction, destruction and readi...
Definition: GridHelpers.hpp:315
Allows viewing a sparse table consisting out of C-array.
Definition: GridHelpers.hpp:43
std::size_t noEntries() const
Get the number of non-zero entries.
Definition: GridHelpers.hpp:75
SparseTableView(int *data, int *offset, std::size_t size_arg)
Creates a sparse table view.
Definition: GridHelpers.hpp:53
std::size_t size() const
Get the size of the table.
Definition: GridHelpers.hpp:69
row_type operator[](std::size_t row) const
Get a row of the the table.
Definition: GridHelpers.hpp:60
Holds the implementation of the CpGrid as a pimple.
Definition: CellQuadrature.hpp:29
face_tag
Connection taxonomy.
Definition: preprocess.h:66
Maps the grid type to the associated type of the cell to faces mapping.
Definition: GridHelpers.hpp:277
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
Definition: IteratorRange.hpp:29
Data structure for an unstructured grid, unstructured meaning that any cell may have an arbitrary num...
Definition: UnstructuredGrid.h:99
int * face_cells
For a face f, face_cells[2*f] and face_cells[2*f + 1] contain the cell indices of the cells adjacent ...
Definition: UnstructuredGrid.h:138