My Project
Loading...
Searching...
No Matches
ReservoirPropertyCommon.hpp
1//===========================================================================
2//
3// File: ReservoirPropertyCommon.hpp
4//
5// Created: Mon Oct 26 08:23:31 2009
6//
7// Author(s): Atgeirr F Rasmussen <atgeirr@sintef.no>
8// B�rd Skaflestad <bard.skaflestad@sintef.no>
9//
10// $Date$
11//
12// $Revision$
13//
14//===========================================================================
15
16/*
17 Copyright 2009, 2010 SINTEF ICT, Applied Mathematics.
18 Copyright 2009, 2010 Statoil ASA.
19
20 This file is part of The Open Reservoir Simulator Project (OpenRS).
21
22 OpenRS is free software: you can redistribute it and/or modify
23 it under the terms of the GNU General Public License as published by
24 the Free Software Foundation, either version 3 of the License, or
25 (at your option) any later version.
26
27 OpenRS is distributed in the hope that it will be useful,
28 but WITHOUT ANY WARRANTY; without even the implied warranty of
29 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 GNU General Public License for more details.
31
32 You should have received a copy of the GNU General Public License
33 along with OpenRS. If not, see <http://www.gnu.org/licenses/>.
34*/
35
36#ifndef OPENRS_RESERVOIRPROPERTYCOMMON_HEADER
37#define OPENRS_RESERVOIRPROPERTYCOMMON_HEADER
38
39#include <opm/input/eclipse/Units/Units.hpp>
40#include <opm/porsol/common/Matrix.hpp>
41
42#include <opm/input/eclipse/Deck/Deck.hpp>
43
44namespace Opm
45{
46
47
48
50 enum PermeabilityKind { ScalarPerm, DiagonalPerm, TensorPerm, None, Invalid };
51
52
53
54
57 template <int dim, class RPImpl, class RockType>
59 {
60 public:
62 typedef ImmutableCMatrix PermTensor;
66 typedef SharedCMatrix SharedPermTensor;
67
69 enum { NumberOfPhases = 2 };
70
73
85 void init(const Opm::Deck&,
86 const std::vector<int>& global_cell,
87 const double perm_threshold = 0.0,
88 const std::string* rock_list_filename = 0,
89 const bool use_jfunction_scaling = true,
90 const double sigma = 1.0,
91 const double theta = 0.0);
92
97 void init(const int num_cells,
98 const double uniform_poro = 0.2,
99 const double uniform_perm = 100.0*Opm::prefix::milli*Opm::unit::darcy);
100
104 void setViscosities(double v1, double v2);
105
109 void setDensities(double d1, double d2);
110
113 double viscosityFirstPhase() const;
114
117 double viscositySecondPhase() const;
118
121 double densityFirstPhase() const;
122
125 double densitySecondPhase() const;
126
130 double porosity(int cell_index) const;
131
135 double ntg(int cell_index) const;
136
140 double swcr(int cell_index) const;
141
145 double sowcr(int cell_index) const;
146
150 PermTensor permeability(int cell_index) const;
151
156
162 template<class Vector>
163 void phaseDensities(int /*cell_index*/, Vector& density) const;
164
167 double densityDifference() const;
168
171 double cflFactor() const;
172
175 double cflFactorGravity() const;
176
179 double cflFactorCapillary() const;
180
185 double capillaryPressure(int cell_index, double saturation) const;
190 double capillaryPressureDeriv(int c, double s) const;
191
192 // @brief Minimum of saturation in rock table.
195 double s_min(int c) const;
196
197 // @brief Maximum of saturation in rock table.
200 double s_max(int c) const;
201
206 double saturationFromCapillaryPressure(int cell_index, double cap_press) const;
207
210 void writeSintefLegacyFormat(const std::string& grid_prefix) const;
211
212 protected:
213 // Methods
214 void assignPorosity(const Opm::Deck& deck,
215 const std::vector<int>& global_cell);
216 void assignNTG(const Opm::Deck& deck,
217 const std::vector<int>& global_cell);
218 void assignSWCR(const Opm::Deck& deck,
219 const std::vector<int>& global_cell);
220 void assignSOWCR(const Opm::Deck& deck,
221 const std::vector<int>& global_cell);
222 void assignPermeability(const Opm::Deck& deck,
223 const std::vector<int>& global_cell,
224 const double perm_threshold);
225 void assignRockTable(const Opm::Deck& deck,
226 const std::vector<int>& global_cell);
227 void readRocks(const std::string& rock_list_file);
228
229 // Supporting Barton/Nackman trick (also known as the curiously recurring template pattern).
230 RPImpl& asImpl();
231
232 // Data members.
233 std::vector<double> porosity_;
234 std::vector<double> ntg_;
235 std::vector<double> swcr_;
236 std::vector<double> sowcr_;
237 std::vector<double> permeability_;
238 std::vector<unsigned char> permfield_valid_;
239 double density1_;
240 double density2_;
241 double viscosity1_;
242 double viscosity2_;
243 double cfl_factor_;
244 double cfl_factor_gravity_;
245 double cfl_factor_capillary_;
246 std::vector<RockType> rock_;
247 std::vector<int> cell_to_rock_;
248 PermeabilityKind permeability_kind_;
249 };
250
251
252} // namespace Opm
253
254#include "ReservoirPropertyCommon_impl.hpp"
255
256#endif // OPENRS_RESERVOIRPROPERTYCOMMON_HEADER
A property class for incompressible two-phase flow.
Definition ReservoirPropertyCommon.hpp:59
double porosity(int cell_index) const
Read-access to porosity.
Definition ReservoirPropertyCommon_impl.hpp:380
ReservoirPropertyCommon()
Default constructor.
Definition ReservoirPropertyCommon_impl.hpp:258
double densitySecondPhase() const
Density of second (oil) phase.
Definition ReservoirPropertyCommon_impl.hpp:373
double cflFactorCapillary() const
A factor useful in gravity cfl computations.
Definition ReservoirPropertyCommon_impl.hpp:461
void setViscosities(double v1, double v2)
Set viscosities of both faces.
Definition ReservoirPropertyCommon_impl.hpp:338
double viscosityFirstPhase() const
Viscosity of first (water) phase.
Definition ReservoirPropertyCommon_impl.hpp:352
double viscositySecondPhase() const
Viscosity of second (oil) phase.
Definition ReservoirPropertyCommon_impl.hpp:359
double saturationFromCapillaryPressure(int cell_index, double cap_press) const
Inverse of the capillary pressure function.
Definition ReservoirPropertyCommon_impl.hpp:519
double densityDifference() const
Difference of densities.
Definition ReservoirPropertyCommon_impl.hpp:440
PermTensor permeability(int cell_index) const
Read-access to permeability.
Definition ReservoirPropertyCommon_impl.hpp:406
double densityFirstPhase() const
Density of first (water) phase.
Definition ReservoirPropertyCommon_impl.hpp:366
double capillaryPressure(int cell_index, double saturation) const
Capillary pressure.
Definition ReservoirPropertyCommon_impl.hpp:467
double capillaryPressureDeriv(int c, double s) const
Derivative of Capillary pressure.
Definition ReservoirPropertyCommon_impl.hpp:480
ImmutableCMatrix PermTensor
Tensor type for read-only access to permeability.
Definition ReservoirPropertyCommon.hpp:62
double swcr(int cell_index) const
Read-access to swcr.
Definition ReservoirPropertyCommon_impl.hpp:392
double cflFactorGravity() const
A factor useful in gravity cfl computations.
Definition ReservoirPropertyCommon_impl.hpp:454
void setDensities(double d1, double d2)
Set densitities of both faces.
Definition ReservoirPropertyCommon_impl.hpp:345
double s_min(int c) const
Definition ReservoirPropertyCommon_impl.hpp:493
OwnCMatrix MutablePermTensor
Tensor type to be used for holding copies of permeability tensors.
Definition ReservoirPropertyCommon.hpp:64
SharedCMatrix SharedPermTensor
Tensor type for read and write access to permeability.
Definition ReservoirPropertyCommon.hpp:66
double ntg(int cell_index) const
Read-access to ntg.
Definition ReservoirPropertyCommon_impl.hpp:386
void writeSintefLegacyFormat(const std::string &grid_prefix) const
Write permeability and porosity in the Sintef legacy format.
Definition ReservoirPropertyCommon_impl.hpp:533
double cflFactor() const
A factor useful in cfl computations.
Definition ReservoirPropertyCommon_impl.hpp:447
void init(const Opm::Deck &, const std::vector< int > &global_cell, const double perm_threshold=0.0, const std::string *rock_list_filename=0, const bool use_jfunction_scaling=true, const double sigma=1.0, const double theta=0.0)
Initialize from a grdecl file.
Definition ReservoirPropertyCommon_impl.hpp:276
void phaseDensities(int, Vector &density) const
Densities for both phases.
Definition ReservoirPropertyCommon_impl.hpp:431
double s_max(int c) const
Definition ReservoirPropertyCommon_impl.hpp:506
SharedPermTensor permeabilityModifiable(int cell_index)
Read- and write-access to permeability.
Definition ReservoirPropertyCommon_impl.hpp:417
double sowcr(int cell_index) const
Read-access to sowcr.
Definition ReservoirPropertyCommon_impl.hpp:398
Inverting small matrices.
Definition ImplicitAssembly.hpp:43
PermeabilityKind
Enum for the kind of permeability field originally retrieved.
Definition ReservoirPropertyCommon.hpp:50
FullMatrix< double, OwnData, COrdering > OwnCMatrix
Convenience typedefs for C-ordered.
Definition Matrix.hpp:579