27#ifndef OPM_ECL_EPS_SCALING_POINTS_HPP
28#define OPM_ECL_EPS_SCALING_POINTS_HPP
34#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
35#include <opm/input/eclipse/EclipseState/Runspec.hpp>
36#include <opm/input/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
37#include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
58template <
class Scalar>
122 std::cout <<
" Swl: " <<
Swl <<
'\n'
123 <<
" Sgl: " <<
Sgl <<
'\n'
124 <<
" Swcr: " <<
Swcr <<
'\n'
125 <<
" Sgcr: " <<
Sgcr <<
'\n'
126 <<
" Sowcr: " <<
Sowcr <<
'\n'
127 <<
" Sogcr: " <<
Sogcr <<
'\n'
128 <<
" Swu: " <<
Swu <<
'\n'
129 <<
" Sgu: " <<
Sgu <<
'\n'
130 <<
" maxPcow: " <<
maxPcow <<
'\n'
131 <<
" maxPcgo: " <<
maxPcgo <<
'\n'
134 <<
" Krwr: " <<
Krwr <<
'\n'
135 <<
" Krgr: " <<
Krgr <<
'\n'
136 <<
" Krorw: " <<
Krorw <<
'\n'
137 <<
" Krorg: " <<
Krorg <<
'\n'
138 <<
" maxKrw: " <<
maxKrw <<
'\n'
139 <<
" maxKrg: " <<
maxKrg <<
'\n'
140 <<
" maxKrow: " <<
maxKrow <<
'\n'
141 <<
" maxKrog: " <<
maxKrog <<
'\n';
151 void extractUnscaled(
const satfunc::RawTableEndPoints& rtep,
152 const satfunc::RawFunctionValues& rfunc,
153 const std::vector<double>::size_type satRegionIdx)
155 this->Swl = rtep.connate.water[satRegionIdx];
156 this->Sgl = rtep.connate.gas [satRegionIdx];
158 this->Swcr = rtep.critical.water [satRegionIdx];
159 this->Sgcr = rtep.critical.gas [satRegionIdx];
160 this->Sowcr = rtep.critical.oil_in_water[satRegionIdx];
161 this->Sogcr = rtep.critical.oil_in_gas [satRegionIdx];
163 this->Swu = rtep.maximum.water[satRegionIdx];
164 this->Sgu = rtep.maximum.gas [satRegionIdx];
166 this->maxPcgo = rfunc.pc.g[satRegionIdx];
167 this->maxPcow = rfunc.pc.w[satRegionIdx];
170 this->pcowLeverettFactor = 1.0;
171 this->pcgoLeverettFactor = 1.0;
173 this->Krwr = rfunc.krw.r [satRegionIdx];
174 this->Krgr = rfunc.krg.r [satRegionIdx];
175 this->Krorw = rfunc.kro.rw[satRegionIdx];
176 this->Krorg = rfunc.kro.rg[satRegionIdx];
178 this->maxKrw = rfunc.krw.max[satRegionIdx];
179 this->maxKrow = rfunc.kro.max[satRegionIdx];
180 this->maxKrog = rfunc.kro.max[satRegionIdx];
181 this->maxKrg = rfunc.krg.max[satRegionIdx];
184 void update(Scalar& targetValue,
const double * value_ptr) {
186 targetValue = *value_ptr;
194 void extractScaled(
const EclipseState& eclState,
195 const EclEpsGridProperties& epsProperties,
196 unsigned activeIndex)
200 update(
Swl, epsProperties.swl(activeIndex));
201 update(
Sgl, epsProperties.sgl(activeIndex));
202 update(
Swcr, epsProperties.swcr(activeIndex));
203 update(
Sgcr, epsProperties.sgcr(activeIndex));
205 update(
Sowcr, epsProperties.sowcr(activeIndex));
206 update(
Sogcr, epsProperties.sogcr(activeIndex));
207 update(
Swu, epsProperties.swu(activeIndex));
208 update(
Sgu, epsProperties.sgu(activeIndex));
209 update(
maxPcow, epsProperties.pcw(activeIndex));
210 update(
maxPcgo, epsProperties.pcg(activeIndex));
212 update(this->Krwr, epsProperties.krwr(activeIndex));
213 update(this->Krgr, epsProperties.krgr(activeIndex));
214 update(this->Krorw, epsProperties.krorw(activeIndex));
215 update(this->Krorg, epsProperties.krorg(activeIndex));
217 update(
maxKrw, epsProperties.krw(activeIndex));
218 update(
maxKrg, epsProperties.krg(activeIndex));
219 update(
maxKrow, epsProperties.kro(activeIndex));
220 update(
maxKrog, epsProperties.kro(activeIndex));
227 if (eclState.getTableManager().useJFunc()) {
228 const auto& jfunc = eclState.getTableManager().getJFunc();
229 const auto& jfuncDir = jfunc.direction();
232 if (jfuncDir == JFunc::Direction::X)
233 perm = epsProperties.permx(activeIndex);
234 else if (jfuncDir == JFunc::Direction::Y)
235 perm = epsProperties.permy(activeIndex);
236 else if (jfuncDir == JFunc::Direction::Z)
237 perm = epsProperties.permz(activeIndex);
238 else if (jfuncDir == JFunc::Direction::XY)
244 double permx = epsProperties.permx(activeIndex);
245 double permy = epsProperties.permy(activeIndex);
248 throw std::runtime_error(
"Illegal direction indicator for the JFUNC "
249 "keyword ("+std::to_string(
int(jfuncDir))+
")");
254 Scalar poro = epsProperties.poro(activeIndex);
255 Scalar alpha = jfunc.alphaFactor();
256 Scalar beta = jfunc.betaFactor();
264 const Scalar Uconst = 0.318316 * 1e5;
267 const auto& jfuncFlag = jfunc.flag();
268 if (jfuncFlag == JFunc::Flag::WATER || jfuncFlag == JFunc::Flag::BOTH) {
271 jfunc.owSurfaceTension();
276 if (jfuncFlag == JFunc::Flag::GAS || jfuncFlag == JFunc::Flag::BOTH) {
279 jfunc.goSurfaceTension();
287 void extractGridPropertyValue_(Scalar& targetValue,
288 const std::vector<double>* propData,
289 unsigned cartesianCellIdx)
294 targetValue = (*propData)[cartesianCellIdx];
304template <
class Scalar>
317 saturationPcPoints_[0] = epsInfo.
Swl;
318 saturationPcPoints_[2] = saturationPcPoints_[1] = epsInfo.
Swu;
321 saturationKrwPoints_[0] = epsInfo.
Swcr;
322 saturationKrwPoints_[1] = 1.0 - epsInfo.
Sowcr - epsInfo.
Sgl;
323 saturationKrwPoints_[2] = epsInfo.
Swu;
329 saturationKrnPoints_[2] = 1.0 - epsInfo.
Sowcr;
330 saturationKrnPoints_[1] = epsInfo.
Swcr + epsInfo.
Sgl;
331 saturationKrnPoints_[0] = epsInfo.
Swl + epsInfo.
Sgl;
336 maxPcnwOrLeverettFactor_ = epsInfo.
maxPcow;
338 Krwr_ = epsInfo.
Krwr;
339 Krnr_ = epsInfo.
Krorw;
348 saturationPcPoints_[0] = 1.0 - epsInfo.
Swl - epsInfo.
Sgu;
349 saturationPcPoints_[2] = saturationPcPoints_[1] = 1.0 - epsInfo.
Swl - epsInfo.
Sgl;
352 saturationKrwPoints_[0] = epsInfo.
Sogcr;
353 saturationKrwPoints_[1] = 1.0 - epsInfo.
Sgcr - epsInfo.
Swl;
354 saturationKrwPoints_[2] = 1.0 - epsInfo.
Swl - epsInfo.
Sgl;
362 saturationKrnPoints_[2] = 1.0 - epsInfo.
Swl - epsInfo.
Sgcr;
363 saturationKrnPoints_[1] = epsInfo.
Sogcr;
364 saturationKrnPoints_[0] = 1.0 - epsInfo.
Swl - epsInfo.
Sgu;
369 maxPcnwOrLeverettFactor_ = epsInfo.
maxPcgo;
371 Krwr_ = epsInfo.
Krorg;
372 Krnr_ = epsInfo.
Krgr;
383 { saturationPcPoints_[pointIdx] = value; }
389 {
return saturationPcPoints_; }
395 { saturationKrwPoints_[pointIdx] = value; }
401 {
return saturationKrwPoints_; }
407 { saturationKrnPoints_[pointIdx] = value; }
413 {
return saturationKrnPoints_; }
419 { maxPcnwOrLeverettFactor_ = value; }
425 {
return maxPcnwOrLeverettFactor_; }
431 { maxPcnwOrLeverettFactor_ = value; }
437 {
return maxPcnwOrLeverettFactor_; }
444 { this->Krwr_ = value; }
451 {
return this->Krwr_; }
470 { this->Krnr_ = value; }
477 {
return this->Krnr_; }
493 std::cout <<
" saturationKrnPoints_[0]: " << saturationKrnPoints_[0] <<
"\n"
494 <<
" saturationKrnPoints_[1]: " << saturationKrnPoints_[1] <<
"\n"
495 <<
" saturationKrnPoints_[2]: " << saturationKrnPoints_[2] <<
"\n";
500 Scalar maxPcnwOrLeverettFactor_;
517 std::array<Scalar, 3> saturationPcPoints_;
520 std::array<Scalar, 3> saturationKrwPoints_;
523 std::array<Scalar, 3> saturationKrnPoints_;
Implements some common averages.
Specifies the configuration used by the endpoint scaling code.
Definition: EclEpsConfig.hpp:57
bool enableLeverettScaling() const
Returns whether the Leverett capillary pressure scaling is enabled.
Definition: EclEpsConfig.hpp:177
Represents the points on the X and Y axis to be scaled if endpoint scaling is used.
Definition: EclEpsScalingPoints.hpp:306
const std::array< Scalar, 3 > & saturationKrwPoints() const
Returns the points used for wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:400
void setMaxKrn(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:482
void print() const
Definition: EclEpsScalingPoints.hpp:491
Scalar krnr() const
Returns non-wetting phase relative permeability at residual saturation of wetting phase.
Definition: EclEpsScalingPoints.hpp:476
Scalar maxKrn() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:488
void setKrnr(Scalar value)
Set non-wetting phase relative permeability at residual saturation of wetting phase.
Definition: EclEpsScalingPoints.hpp:469
void setMaxKrw(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:456
Scalar maxKrw() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:462
void setSaturationKrwPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for wetting-phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:394
void setMaxPcnw(Scalar value)
Sets the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:418
Scalar krwr() const
Returns wetting-phase relative permeability at residual saturation of non-wetting phase.
Definition: EclEpsScalingPoints.hpp:450
Scalar maxPcnw() const
Returns the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:424
void setKrwr(Scalar value)
Set wetting-phase relative permeability at residual saturation of non-wetting phase.
Definition: EclEpsScalingPoints.hpp:443
void setSaturationPcPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:382
const std::array< Scalar, 3 > & saturationKrnPoints() const
Returns the points used for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:412
void setLeverettFactor(Scalar value)
Sets the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:430
const std::array< Scalar, 3 > & saturationPcPoints() const
Returns the points used for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:388
void init(const EclEpsScalingPointsInfo< Scalar > &epsInfo, const EclEpsConfig &config, EclTwoPhaseSystemType epsSystemType)
Assigns the scaling points which actually ought to be used.
Definition: EclEpsScalingPoints.hpp:311
void setSaturationKrnPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:406
Scalar leverettFactor() const
Returns the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:436
Definition: Air_Mesitylene.hpp:34
EclTwoPhaseSystemType
Specified which fluids are involved in a given twophase material law for endpoint scaling.
Definition: EclEpsConfig.hpp:43
@ EclGasOilSystem
Definition: EclEpsConfig.hpp:44
@ EclGasWaterSystem
Definition: EclEpsConfig.hpp:46
@ EclOilWaterSystem
Definition: EclEpsConfig.hpp:45
Scalar arithmeticMean(Scalar x, Scalar y)
Computes the arithmetic average of two values.
Definition: Means.hpp:45
ReturnEval_< Evaluation1, Evaluation2 >::type pow(const Evaluation1 &base, const Evaluation2 &exp)
Definition: MathToolbox.hpp:416
This structure represents all values which can be possibly used as scaling points in the endpoint sca...
Definition: EclEpsScalingPoints.hpp:60
Scalar Swu
Definition: EclEpsScalingPoints.hpp:72
Scalar Krwr
Definition: EclEpsScalingPoints.hpp:85
void print() const
Definition: EclEpsScalingPoints.hpp:120
Scalar pcgoLeverettFactor
Definition: EclEpsScalingPoints.hpp:82
Scalar Sogcr
Definition: EclEpsScalingPoints.hpp:69
bool operator==(const EclEpsScalingPointsInfo< Scalar > &data) const
Definition: EclEpsScalingPoints.hpp:96
Scalar pcowLeverettFactor
Definition: EclEpsScalingPoints.hpp:81
Scalar Krgr
Definition: EclEpsScalingPoints.hpp:86
Scalar Swl
Definition: EclEpsScalingPoints.hpp:62
Scalar Krorg
Definition: EclEpsScalingPoints.hpp:88
Scalar Sgu
Definition: EclEpsScalingPoints.hpp:73
Scalar maxKrow
Definition: EclEpsScalingPoints.hpp:92
Scalar Krorw
Definition: EclEpsScalingPoints.hpp:87
Scalar Sgl
Definition: EclEpsScalingPoints.hpp:63
Scalar maxKrw
Definition: EclEpsScalingPoints.hpp:91
Scalar maxKrg
Definition: EclEpsScalingPoints.hpp:94
Scalar maxPcgo
Definition: EclEpsScalingPoints.hpp:77
Scalar maxKrog
Definition: EclEpsScalingPoints.hpp:93
Scalar maxPcow
Definition: EclEpsScalingPoints.hpp:76
Scalar Sgcr
Definition: EclEpsScalingPoints.hpp:67
Scalar Sowcr
Definition: EclEpsScalingPoints.hpp:68
Scalar Swcr
Definition: EclEpsScalingPoints.hpp:66