RelPermUtils.hpp
1 /*
2  Copyright 2010 Statoil ASA.
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 
24 #ifndef OPM_UPSCALING_RELPERM_UTILS_HPP
25 #define OPM_UPSCALING_RELPERM_UTILS_HPP
26 
27 #include <opm/parser/eclipse/Parser/ParseContext.hpp>
28 #include <opm/parser/eclipse/Parser/Parser.hpp>
29 #include <opm/parser/eclipse/Deck/Deck.hpp>
30 
31 #include <opm/core/utility/MonotCubicInterpolator.hpp>
32 
33 #include <opm/upscaling/ParserAdditions.hpp>
34 #include <opm/upscaling/SinglePhaseUpscaler.hpp>
35 
36 #include <array>
37 #include <map>
38 #include <memory>
39 #include <tuple>
40 #include <vector>
41 
42 namespace Opm {
47  double getVoigtValue(const SinglePhaseUpscaler::permtensor_t& K, int voigt_idx);
48 
53  void setVoigtValue(SinglePhaseUpscaler::permtensor_t& K, int voigt_idx, double val);
54 
57  public:
58  bool isMaster;
61  int points;
64  std::vector<double> WaterSaturation;
66  std::vector<int> satnums;
67  std::vector<MonotCubicInterpolator> InvJfunctions;
68  std::array<std::array<std::vector<MonotCubicInterpolator>,2>,3> Krfunctions;
70  SinglePhaseUpscaler::BoundaryConditionType boundaryCondition;
71  std::vector<MonotCubicInterpolator> SwPcfunctions;
72  std::string saturationstring;
73  size_t tesselatedCells;
74  double volume;
75  double poreVolume;
76  std::vector<double> pressurePoints;
77 
86  template <class String>
87  static Deck
88  parseEclipseFile(const String& eclipseFileName);
89 
94  RelPermUpscaleHelper(int mpi_rank, std::map<std::string,std::string>& options_);
95 
97  void collectResults();
98 
103  std::vector<std::vector<double>> getRelPerm(int phase) const;
104 
107 
114  void sanityCheckInput(const Opm::Deck& deck,
115  const double minPerm,
116  const double maxPerm,
117  const double minPoro);
118 
121 
125 
133  double tesselateGrid(const Opm::Deck& deck);
134 
138 
143 
147 
152  std::tuple<double,double> upscalePermeability(int mpi_rank);
153 
154  private:
160  bool checkCurve(MonotCubicInterpolator& func);
161 
162  double Swir;
163  double Swor;
164  double Pcmin;
165  double Pcmax;
166  double critRelpThresh;
167  std::vector<int> node_vs_pressurepoint;
168  std::array<std::vector<std::vector<double>>,2> PhasePerm;
169  SinglePhaseUpscaler::permtensor_t permTensorInv;
170  SinglePhaseUpscaler upscaler;
171  std::array<std::vector<double>,3> perms;
172  std::vector<double> poros;
173  std::vector<double> zcorns;
174  double minSinglePhasePerm;
175  std::vector<double> cellPoreVolumes;
176  MonotCubicInterpolator WaterSaturationVsCapPressure;
177 
178 #if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP)
179  public: // Intrusive unit testing of calculcateCellPressureGradients()
180 #endif // UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP
181  std::vector<double> dP;
182 
183 #if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP)
184  private:
185 #endif // UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP
186 
187  std::map<std::string,std::string>& options;
188  };
189 
190  template <class String>
191  Deck
192  RelPermUpscaleHelper::parseEclipseFile(const String& eclipseFileName)
193  {
194  auto parser = Parser{};
196 
197  return parser.parseFile(eclipseFileName, ParseContext{});
198  }
199 }
200 
201 #endif // OPM_UPSCALING_RELPERM_UTILS_HPP
double volume
Total volume.
Definition: RelPermUtils.hpp:74
void calculateCellPressureGradients()
Find cell center pressure gradient for every cell.
Definition: RelPermUtils.cpp:687
bool doEclipseCheck
Whether to check that input relperm curves include relperm at critical saturation points...
Definition: RelPermUtils.hpp:59
double tesselateGrid(const Opm::Deck &deck)
Tesselate grid.
Definition: RelPermUtils.cpp:650
std::array< std::array< std::vector< MonotCubicInterpolator >, 2 >, 3 > Krfunctions
Relperm-curves for each (component->phase->stone type)
Definition: RelPermUtils.hpp:69
std::vector< std::vector< double > > getRelPerm(int phase) const
Calculate relperm values from phase permeabilities.
Definition: RelPermUtils.cpp:309
void collectResults()
Collect results from all MPI nodes.
Definition: RelPermUtils.cpp:237
int tensorElementCount
Number of independent elements in resulting tensor.
Definition: RelPermUtils.hpp:62
static Deck parseEclipseFile(const String &eclipseFileName)
Form Deck object from input file (ECLIPSE format)
Definition: RelPermUtils.hpp:192
bool anisotropic_input
Whether input eclipse file has diagonal anisotrophy.
Definition: RelPermUtils.hpp:60
void setupBoundaryConditions()
Setup requested boundary conditions.
Definition: RelPermUtils.cpp:635
Class for immiscible dead oil and dry gas.
Definition: applier.hpp:18
bool isMaster
Whether this is the master MPI node or not.
Definition: RelPermUtils.hpp:58
void upscaleCapillaryPressure()
Upscale capillary pressure.
Definition: RelPermUtils.cpp:885
void addNonStandardUpscalingKeywords(Parser &parser)
This function registers non-standard keywords used by opm-upscaling in a parser object.
Definition: ParserAdditions.cpp:41
std::tuple< double, double > upscalePermeability(int mpi_rank)
Upscale permeabilities.
Definition: RelPermUtils.cpp:994
SinglePhaseUpscaler::permtensor_t permTensor
Tensor of upscaled results.
Definition: RelPermUtils.hpp:65
std::vector< int > satnums
Cell satnums.
Definition: RelPermUtils.hpp:66
double poreVolume
Total pore volume.
Definition: RelPermUtils.hpp:75
void upscaleSinglePhasePermeability()
Upscale single phase permeability.
Definition: RelPermUtils.cpp:353
SinglePhaseUpscaler::BoundaryConditionType boundaryCondition
Boundary conditions to use.
Definition: RelPermUtils.hpp:70
std::vector< MonotCubicInterpolator > InvJfunctions
Inverse of the loaded J-functions.
Definition: RelPermUtils.hpp:67
void calculateMinMaxCapillaryPressure()
Calculate minimum and maximum capillary pressures.
Definition: RelPermUtils.cpp:730
std::vector< double > pressurePoints
Vector of capillary pressure points between Swor and Swir.
Definition: RelPermUtils.hpp:76
ResProp::MutablePermTensor permtensor_t
A type for the upscaled permeability.
Definition: UpscalerBase.hpp:66
RelPermUpscaleHelper(int mpi_rank, std::map< std::string, std::string > &options_)
Default constructor.
Definition: RelPermUtils.cpp:206
double getVoigtValue(const SinglePhaseUpscaler::permtensor_t &K, int voigt_idx)
Get value from tensor.
Definition: RelPermUtils.cpp:173
bool upscaleBothPhases
Whether to upscale both phases.
Definition: RelPermUtils.hpp:63
void checkCriticalSaturations()
Check that input relperm curevs specify critical saturations.
Definition: RelPermUtils.cpp:615
void setVoigtValue(SinglePhaseUpscaler::permtensor_t &K, int voigt_idx, double val)
Set value in tensor.
Definition: RelPermUtils.cpp:190
int points
Number of saturation points to upscale for.
Definition: RelPermUtils.hpp:61
std::vector< MonotCubicInterpolator > SwPcfunctions
Holds Sw(Pc) for each rocktype.
Definition: RelPermUtils.hpp:71
void sanityCheckInput(const Opm::Deck &deck, const double minPerm, const double maxPerm, const double minPoro)
Do sanity checks for input file.
Definition: RelPermUtils.cpp:390
size_t tesselatedCells
Number of "active" cells (Sintef interpretation of "active")
Definition: RelPermUtils.hpp:73
Helper class for relperm upscaling applications.
Definition: RelPermUtils.hpp:56
std::string saturationstring
Fluid system type.
Definition: RelPermUtils.hpp:72
A class for doing single phase (permeability) upscaling.
Definition: SinglePhaseUpscaler.hpp:50
std::vector< double > WaterSaturation
Re-upscaled water saturation for the computed pressure points.
Definition: RelPermUtils.hpp:64