24#ifndef OPM_UPSCALING_RELPERM_UTILS_HPP
25#define OPM_UPSCALING_RELPERM_UTILS_HPP
27#include <opm/input/eclipse/Parser/Parser.hpp>
28#include <opm/input/eclipse/Deck/Deck.hpp>
30#include <opm/common/utility/numeric/MonotCubicInterpolator.hpp>
32#include <opm/upscaling/ParserAdditions.hpp>
33#include <opm/upscaling/SinglePhaseUpscaler.hpp>
68 std::array<std::array<std::vector<MonotCubicInterpolator>,2>,3>
Krfunctions;
85 template <
class String>
102 std::vector<std::vector<double>>
getRelPerm(
int phase)
const;
114 const double minPerm,
115 const double maxPerm,
116 const double minPoro);
159 bool checkCurve(MonotCubicInterpolator& func);
165 double critRelpThresh;
166 std::vector<int> node_vs_pressurepoint;
167 std::array<std::vector<std::vector<double>>,2> PhasePerm;
170 std::array<std::vector<double>,3> perms;
171 std::vector<double> poros;
172 std::vector<double> zcorns;
173 double minSinglePhasePerm;
174 std::vector<double> cellPoreVolumes;
175 MonotCubicInterpolator WaterSaturationVsCapPressure;
177#if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP)
180 std::vector<double> dP;
182#if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP)
186 std::map<std::string,std::string>& options;
189 template <
class String>
193 auto parser = Parser{};
196 return parser.parseFile(eclipseFileName);
Helper class for relperm upscaling applications.
Definition RelPermUtils.hpp:55
std::tuple< double, double > upscalePermeability(int mpi_rank)
Upscale permeabilities.
Definition RelPermUtils.cpp:997
SinglePhaseUpscaler::BoundaryConditionType boundaryCondition
Boundary conditions to use.
Definition RelPermUtils.hpp:69
std::vector< int > satnums
Cell satnums.
Definition RelPermUtils.hpp:65
double volume
Total volume.
Definition RelPermUtils.hpp:73
void upscaleSinglePhasePermeability()
Upscale single phase permeability.
Definition RelPermUtils.cpp:355
void upscaleCapillaryPressure()
Upscale capillary pressure.
Definition RelPermUtils.cpp:888
SinglePhaseUpscaler::permtensor_t permTensor
Tensor of upscaled results.
Definition RelPermUtils.hpp:64
static Deck parseEclipseFile(const String &eclipseFileName)
Form Deck object from input file (ECLIPSE format)
Definition RelPermUtils.hpp:191
std::vector< MonotCubicInterpolator > SwPcfunctions
Holds Sw(Pc) for each rocktype.
Definition RelPermUtils.hpp:70
void checkCriticalSaturations()
Check that input relperm curevs specify critical saturations.
Definition RelPermUtils.cpp:618
bool isMaster
Whether this is the master MPI node or not.
Definition RelPermUtils.hpp:57
std::vector< MonotCubicInterpolator > InvJfunctions
Inverse of the loaded J-functions.
Definition RelPermUtils.hpp:66
std::vector< double > pressurePoints
Vector of capillary pressure points between Swor and Swir.
Definition RelPermUtils.hpp:75
void sanityCheckInput(const Opm::Deck &deck, const double minPerm, const double maxPerm, const double minPoro)
Do sanity checks for input file.
Definition RelPermUtils.cpp:392
size_t tesselatedCells
Number of "active" cells (Sintef interpretation of "active")
Definition RelPermUtils.hpp:72
void calculateMinMaxCapillaryPressure()
Calculate minimum and maximum capillary pressures.
Definition RelPermUtils.cpp:733
int tensorElementCount
Number of independent elements in resulting tensor.
Definition RelPermUtils.hpp:61
std::array< std::array< std::vector< MonotCubicInterpolator >, 2 >, 3 > Krfunctions
Relperm-curves for each (component->phase->stone type)
Definition RelPermUtils.hpp:68
int points
Number of saturation points to upscale for.
Definition RelPermUtils.hpp:60
bool doEclipseCheck
Whether to check that input relperm curves include relperm at critical saturation points.
Definition RelPermUtils.hpp:58
bool anisotropic_input
Whether input eclipse file has diagonal anisotrophy.
Definition RelPermUtils.hpp:59
void collectResults()
Collect results from all MPI nodes.
Definition RelPermUtils.cpp:239
std::vector< double > WaterSaturation
Re-upscaled water saturation for the computed pressure points.
Definition RelPermUtils.hpp:63
double tesselateGrid(const Opm::Deck &deck)
Tesselate grid.
Definition RelPermUtils.cpp:653
bool upscaleBothPhases
Whether to upscale both phases.
Definition RelPermUtils.hpp:62
std::string saturationstring
Fluid system type.
Definition RelPermUtils.hpp:71
std::vector< std::vector< double > > getRelPerm(int phase) const
Calculate relperm values from phase permeabilities.
Definition RelPermUtils.cpp:311
double poreVolume
Total pore volume.
Definition RelPermUtils.hpp:74
void setupBoundaryConditions()
Setup requested boundary conditions.
Definition RelPermUtils.cpp:638
void calculateCellPressureGradients()
Find cell center pressure gradient for every cell.
Definition RelPermUtils.cpp:690
A class for doing single phase (permeability) upscaling.
Definition SinglePhaseUpscaler.hpp:51
ResProp::MutablePermTensor permtensor_t
A type for the upscaled permeability.
Definition UpscalerBase.hpp:66
Inverting small matrices.
Definition ImplicitAssembly.hpp:43
void addNonStandardUpscalingKeywords(Parser &parser)
This function registers non-standard keywords used by opm-upscaling in a parser object.
Definition ParserAdditions.cpp:41