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
30namespace Opm
31{
32
33class EclipseGrid;
34
35namespace UgGridHelpers
36{
37
44{
45public:
48
54 SparseTableView(int* data, int *offset, std::size_t size_arg)
55 : data_(data), offset_(offset), size_(size_arg)
56 {}
57
61 row_type operator[](std::size_t row) const
62 {
63 assert(row<=size());
64 return row_type{data_ + offset_[row],
65 data_ + offset_[row+1]};
66 }
67
70 std::size_t size() const
71 {
72 return size_;
73 }
74
76 std::size_t noEntries() const
77 {
78 return offset_[size_];
79 }
80
81private:
83 const int* data_;
87 const int* offset_;
89 std::size_t size_;
90};
91
93int numCells(const UnstructuredGrid& grid);
94
96int numFaces(const UnstructuredGrid& grid);
97
99int dimensions(const UnstructuredGrid& grid);
100
102int numCellFaces(const UnstructuredGrid& grid);
103
105const int* cartDims(const UnstructuredGrid& grid);
106
111const int* globalCell(const UnstructuredGrid& grid);
112
113#if HAVE_ECL_INPUT
118std::vector<int> createACTNUM(const UnstructuredGrid& grid);
119#endif
120
121
127template<class G>
129{
130};
131
132template<>
134{
135 typedef const double* IteratorType;
136 typedef const double* ValueType;
137};
138
144beginCellCentroids(const UnstructuredGrid& grid);
145
146
150double cellCenterDepth(const UnstructuredGrid& grid, int cell_index);
151
157Dune::FieldVector<double,3> faceCenterEcl(const UnstructuredGrid& grid, int cell_index, int face_tag);
158
165Dune::FieldVector<double,3> faceAreaNormalEcl(const UnstructuredGrid& grid, int face_index);
166
167
172double cellCentroidCoordinate(const UnstructuredGrid& grid, int cell_index,
173 int coordinate);
174
175
179const double* cellCentroid(const UnstructuredGrid& grid, int cell_index);
180
181
185double cellVolume(const UnstructuredGrid& grid, int cell_index);
186
192template<class T>
194{
195};
196
197template<>
199{
200 typedef const double* IteratorType;
201};
202
223#if HAVE_ECL_INPUT
226Opm::EclipseGrid createEclipseGrid(const UnstructuredGrid& grid, const Opm::EclipseGrid& inputGrid );
227#endif
228
230const double* beginCellVolumes(const UnstructuredGrid& grid);
231
233const double* endCellVolumes(const UnstructuredGrid& grid);
234
235
236
242template<class G>
244{
245};
246
247template<>
249{
250 typedef const double* IteratorType;
251 typedef const double* ValueType;
252};
253
256beginFaceCentroids(const UnstructuredGrid& grid);
257
262faceCentroid(const UnstructuredGrid& grid, int face_index);
263
267const double* faceNormal(const UnstructuredGrid& grid, int face_index);
268
272double faceArea(const UnstructuredGrid& grid, int face_index);
273
278template<class T>
280{
281};
282
283template<>
285{
286 typedef SparseTableView Type;
287};
288
293template<class T>
295{
296};
297
298template<>
300{
301 typedef SparseTableView Type;
302};
303
306cell2Faces(const UnstructuredGrid& grid);
307
310face2Vertices(const UnstructuredGrid& grid);
311
315const double* vertexCoordinates(const UnstructuredGrid& grid, int index);
316
318{
319public:
321 : face_cells_(grid.face_cells)
322 {}
323 int operator()(int face_index, int local_index) const
324 {
325 return face_cells_[2*face_index+local_index];
326 }
327private:
328 const int* face_cells_;
329};
330
335template<class T>
337{};
338
339template<>
341{
342 typedef FaceCellsProxy Type;
343};
344
347
353template<class T>
354T* increment(T* cc, int i, int dim)
355{
356 return cc+(i*dim);
357}
362template<class T>
363T increment(const T& t, int i, int)
364{
365 return t+i;
366}
367
372template<class T>
373double getCoordinate(T* cc, int i)
374{
375 return cc[i];
376}
377
383template<class T>
384double getCoordinate(T t, int i)
385{
386 return (*t)[i];
387}
388
389
390
391
392
393} // end namespace UGGridHelpers
394} // end namespace OPM
395#endif
Main OPM-Core grid data structure along with helper functions for construction, destruction and readi...
Definition: GridHelpers.hpp:318
Allows viewing a sparse table consisting out of C-array.
Definition: GridHelpers.hpp:44
std::size_t noEntries() const
Get the number of non-zero entries.
Definition: GridHelpers.hpp:76
SparseTableView(int *data, int *offset, std::size_t size_arg)
Creates a sparse table view.
Definition: GridHelpers.hpp:54
std::size_t size() const
Get the size of the table.
Definition: GridHelpers.hpp:70
row_type operator[](std::size_t row) const
Get a row of the the table.
Definition: GridHelpers.hpp:61
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:280
Traits of the cell centroids of a grid.
Definition: GridHelpers.hpp:129
The mapping of the grid type to type of the iterator over the cell volumes.
Definition: GridHelpers.hpp:194
Maps the grid type to the associated type of the face to vertices mapping.
Definition: GridHelpers.hpp:295
Traits of the face to attached cell mappping of a grid.
Definition: GridHelpers.hpp:337
Traits of the face centroids of a grid.
Definition: GridHelpers.hpp:244
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