RockFromDeck.hpp
1 /*
2  Copyright 2012 SINTEF ICT, Applied Mathematics.
3 
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 3 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 */
19 
20 #ifndef OPM_ROCKFROMDECK_HEADER_INCLUDED
21 #define OPM_ROCKFROMDECK_HEADER_INCLUDED
22 
23 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
24 
25 #include <vector>
26 
27 struct UnstructuredGrid;
28 
29 namespace Opm
30 {
31 
33  {
34  // BlackoilPropsDataHandle needs mutable
35  // access to porosity and permeability
36  friend class BlackoilPropsDataHandle;
37 
38  public:
40  RockFromDeck();
43  explicit RockFromDeck(std::size_t number_of_cells);
50  void init(const Opm::EclipseState& eclState,
51  int number_of_cells, const int* global_cell,
52  const int* cart_dims);
53 
55  int numDimensions() const
56  {
57  return 3;
58  }
59 
61  int numCells() const
62  {
63  return porosity_.size();
64  }
65 
67  const double* porosity() const
68  {
69  return &porosity_[0];
70  }
71 
75  const double* permeability() const
76  {
77  return &permeability_[0];
78  }
79 
89  static
90  void extractInterleavedPermeability(const Opm::EclipseState& eclState,
91  const int number_of_cells,
92  const int* global_cell,
93  const int* cart_dims,
94  const double perm_threshold,
95  std::vector<double>& permeability);
96 
97  private:
98  void assignPorosity(const Opm::EclipseState& eclState,
99  int number_of_cells,
100  const int* global_cell);
101 
102  std::vector<double> porosity_;
103  std::vector<double> permeability_;
104  std::vector<unsigned char> permfield_valid_;
105  };
106 
107 
108 
109 } // namespace Opm
110 
111 
112 #endif // OPM_ROCKFROMDECK_HEADER_INCLUDED
const double * porosity() const
Definition: RockFromDeck.hpp:67
const double * permeability() const
Definition: RockFromDeck.hpp:75
int numCells() const
Definition: RockFromDeck.hpp:61
void init(const Opm::EclipseState &eclState, int number_of_cells, const int *global_cell, const int *cart_dims)
Initialize from deck and cell mapping.
Definition: RockFromDeck.cpp:81
static void extractInterleavedPermeability(const Opm::EclipseState &eclState, const int number_of_cells, const int *global_cell, const int *cart_dims, const double perm_threshold, std::vector< double > &permeability)
Convert the permeabilites for the logically Cartesian grid in EclipseState to an array of size number...
Definition: RockFromDeck.cpp:110
Definition: AnisotropicEikonal.cpp:446
int numDimensions() const
Definition: RockFromDeck.hpp:55
Definition: RockFromDeck.hpp:32
RockFromDeck()
Default constructor.
Definition: RockFromDeck.cpp:70