24 #ifndef OPM_UPSCALING_RELPERM_UTILS_HPP 25 #define OPM_UPSCALING_RELPERM_UTILS_HPP 27 #include <opm/parser/eclipse/Parser/ParseContext.hpp> 28 #include <opm/parser/eclipse/Parser/Parser.hpp> 29 #include <opm/parser/eclipse/Deck/Deck.hpp> 31 #include <opm/core/utility/MonotCubicInterpolator.hpp> 33 #include <opm/upscaling/ParserAdditions.hpp> 34 #include <opm/upscaling/SinglePhaseUpscaler.hpp> 68 std::array<std::array<std::vector<MonotCubicInterpolator>,2>,3>
Krfunctions;
86 template <
class String>
103 std::vector<std::vector<double>>
getRelPerm(
int phase)
const;
115 const double minPerm,
116 const double maxPerm,
117 const double minPoro);
160 bool checkCurve(MonotCubicInterpolator& func);
166 double critRelpThresh;
167 std::vector<int> node_vs_pressurepoint;
168 std::array<std::vector<std::vector<double>>,2> PhasePerm;
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;
178 #if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP) 180 #endif // UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP 181 std::vector<double> dP;
183 #if defined(UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP) 185 #endif // UNITTEST_TRESPASS_PRIVATE_PROPERTY_DP 187 std::map<std::string,std::string>& options;
190 template <
class String>
194 auto parser = Parser{};
197 return parser.parseFile(eclipseFileName, ParseContext{});
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