27 #ifndef OPM_ECL_EPS_SCALING_POINTS_HPP
28 #define OPM_ECL_EPS_SCALING_POINTS_HPP
34 #include <opm/parser/eclipse/Deck/Deck.hpp>
35 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37 #include <opm/parser/eclipse/EclipseState/Runspec.hpp>
38 #include <opm/parser/eclipse/EclipseState/Grid/SatfuncPropertyInitializers.hpp>
39 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
60 template <
class Scalar>
83 Scalar pcowLeverettFactor;
84 Scalar pcgoLeverettFactor;
100 return Swl == data.Swl &&
104 Sowcr == data.Sowcr &&
105 Sogcr == data.Sogcr &&
108 maxPcow == data.maxPcow &&
109 maxPcgo == data.maxPcgo &&
110 pcowLeverettFactor == data.pcowLeverettFactor &&
111 pcgoLeverettFactor == data.pcgoLeverettFactor &&
114 Krorw == data.Krorw &&
115 Krorg == data.Krorg &&
116 maxKrw == data.maxKrw &&
117 maxKrow == data.maxKrow &&
118 maxKrog == data.maxKrog &&
119 maxKrg == data.maxKrg;
124 std::cout <<
" Swl: " << Swl <<
'\n'
125 <<
" Sgl: " << Sgl <<
'\n'
126 <<
" Swcr: " << Swcr <<
'\n'
127 <<
" Sgcr: " << Sgcr <<
'\n'
128 <<
" Sowcr: " << Sowcr <<
'\n'
129 <<
" Sogcr: " << Sogcr <<
'\n'
130 <<
" Swu: " << Swu <<
'\n'
131 <<
" Sgu: " << Sgu <<
'\n'
132 <<
" maxPcow: " << maxPcow <<
'\n'
133 <<
" maxPcgo: " << maxPcgo <<
'\n'
134 <<
" pcowLeverettFactor: " << pcowLeverettFactor <<
'\n'
135 <<
" pcgoLeverettFactor: " << pcgoLeverettFactor <<
'\n'
136 <<
" Krwr: " << Krwr <<
'\n'
137 <<
" Krgr: " << Krgr <<
'\n'
138 <<
" Krorw: " << Krorw <<
'\n'
139 <<
" Krorg: " << Krorg <<
'\n'
140 <<
" maxKrw: " << maxKrw <<
'\n'
141 <<
" maxKrg: " << maxKrg <<
'\n'
142 <<
" maxKrow: " << maxKrow <<
'\n'
143 <<
" maxKrog: " << maxKrog <<
'\n';
153 void extractUnscaled(
const satfunc::RawTableEndPoints& rtep,
154 const satfunc::RawFunctionValues& rfunc,
155 const std::vector<double>::size_type satRegionIdx)
157 this->Swl = rtep.connate.water[satRegionIdx];
158 this->Sgl = rtep.connate.gas [satRegionIdx];
160 this->Swcr = rtep.critical.water [satRegionIdx];
161 this->Sgcr = rtep.critical.gas [satRegionIdx];
162 this->Sowcr = rtep.critical.oil_in_water[satRegionIdx];
163 this->Sogcr = rtep.critical.oil_in_gas [satRegionIdx];
165 this->Swu = rtep.maximum.water[satRegionIdx];
166 this->Sgu = rtep.maximum.gas [satRegionIdx];
168 this->maxPcgo = rfunc.pc.g[satRegionIdx];
169 this->maxPcow = rfunc.pc.w[satRegionIdx];
172 this->pcowLeverettFactor = 1.0;
173 this->pcgoLeverettFactor = 1.0;
175 this->Krwr = rfunc.krw.r [satRegionIdx];
176 this->Krgr = rfunc.krg.r [satRegionIdx];
177 this->Krorw = rfunc.kro.rw[satRegionIdx];
178 this->Krorg = rfunc.kro.rg[satRegionIdx];
180 this->maxKrw = rfunc.krw.max[satRegionIdx];
181 this->maxKrow = rfunc.kro.max[satRegionIdx];
182 this->maxKrog = rfunc.kro.max[satRegionIdx];
183 this->maxKrg = rfunc.krg.max[satRegionIdx];
186 void update(Scalar& targetValue,
const double * value_ptr) {
188 targetValue = *value_ptr;
196 void extractScaled(
const EclipseState& eclState,
198 unsigned activeIndex)
202 update(Swl, epsProperties.swl(activeIndex));
203 update(Sgl, epsProperties.sgl(activeIndex));
204 update(Swcr, epsProperties.swcr(activeIndex));
205 update(Sgcr, epsProperties.sgcr(activeIndex));
207 update(Sowcr, epsProperties.sowcr(activeIndex));
208 update(Sogcr, epsProperties.sogcr(activeIndex));
209 update(Swu, epsProperties.swu(activeIndex));
210 update(Sgu, epsProperties.sgu(activeIndex));
211 update(maxPcow, epsProperties.pcw(activeIndex));
212 update(maxPcgo, epsProperties.pcg(activeIndex));
214 update(this->Krwr, epsProperties.krwr(activeIndex));
215 update(this->Krgr, epsProperties.krgr(activeIndex));
216 update(this->Krorw, epsProperties.krorw(activeIndex));
217 update(this->Krorg, epsProperties.krorg(activeIndex));
219 update(maxKrw, epsProperties.krw(activeIndex));
220 update(maxKrg, epsProperties.krg(activeIndex));
221 update(maxKrow, epsProperties.kro(activeIndex));
222 update(maxKrog, epsProperties.kro(activeIndex));
227 pcowLeverettFactor = 1.0;
228 pcgoLeverettFactor = 1.0;
229 if (eclState.getTableManager().useJFunc()) {
230 const auto& jfunc = eclState.getTableManager().getJFunc();
231 const auto& jfuncDir = jfunc.direction();
234 if (jfuncDir == JFunc::Direction::X)
235 perm = epsProperties.permx(activeIndex);
236 else if (jfuncDir == JFunc::Direction::Y)
237 perm = epsProperties.permy(activeIndex);
238 else if (jfuncDir == JFunc::Direction::Z)
239 perm = epsProperties.permz(activeIndex);
240 else if (jfuncDir == JFunc::Direction::XY)
246 double permx = epsProperties.permx(activeIndex);
247 double permy = epsProperties.permy(activeIndex);
250 throw std::runtime_error(
"Illegal direction indicator for the JFUNC "
251 "keyword ("+std::to_string(
int(jfuncDir))+
")");
256 Scalar poro = epsProperties.poro(activeIndex);
257 Scalar alpha = jfunc.alphaFactor();
258 Scalar beta = jfunc.betaFactor();
262 Scalar commonFactor = std::pow(poro, alpha)/std::pow(perm, beta);
266 const Scalar Uconst = 0.318316 * 1e5;
269 const auto& jfuncFlag = jfunc.flag();
270 if (jfuncFlag == JFunc::Flag::WATER || jfuncFlag == JFunc::Flag::BOTH) {
273 jfunc.owSurfaceTension();
274 pcowLeverettFactor = commonFactor*gamma*Uconst;
278 if (jfuncFlag == JFunc::Flag::GAS || jfuncFlag == JFunc::Flag::BOTH) {
281 jfunc.goSurfaceTension();
282 pcgoLeverettFactor = commonFactor*gamma*Uconst;
289 void extractGridPropertyValue_(Scalar& targetValue,
290 const std::vector<double>* propData,
291 unsigned cartesianCellIdx)
296 targetValue = (*propData)[cartesianCellIdx];
306 template <
class Scalar>
317 if (epsSystemType == EclOilWaterSystem) {
319 saturationPcPoints_[0] = epsInfo.Swl;
320 saturationPcPoints_[2] = saturationPcPoints_[1] = epsInfo.Swu;
323 saturationKrwPoints_[0] = epsInfo.Swcr;
324 saturationKrwPoints_[1] = 1.0 - epsInfo.Sowcr - epsInfo.Sgl;
325 saturationKrwPoints_[2] = epsInfo.Swu;
331 saturationKrnPoints_[2] = 1.0 - epsInfo.Sowcr;
332 saturationKrnPoints_[1] = epsInfo.Swcr + epsInfo.Sgl;
333 saturationKrnPoints_[0] = epsInfo.Swl + epsInfo.Sgl;
336 maxPcnwOrLeverettFactor_ = epsInfo.pcowLeverettFactor;
338 maxPcnwOrLeverettFactor_ = epsInfo.maxPcow;
340 Krwr_ = epsInfo.Krwr;
341 Krnr_ = epsInfo.Krorw;
343 maxKrw_ = epsInfo.maxKrw;
344 maxKrn_ = epsInfo.maxKrow;
347 assert((epsSystemType == EclGasOilSystem) ||(epsSystemType == EclGasWaterSystem) );
350 saturationPcPoints_[0] = 1.0 - epsInfo.Swl - epsInfo.Sgu;
351 saturationPcPoints_[2] = saturationPcPoints_[1] = 1.0 - epsInfo.Swl - epsInfo.Sgl;
354 saturationKrwPoints_[0] = epsInfo.Sogcr;
355 saturationKrwPoints_[1] = 1.0 - epsInfo.Sgcr - epsInfo.Swl;
356 saturationKrwPoints_[2] = 1.0 - epsInfo.Swl - epsInfo.Sgl;
364 saturationKrnPoints_[2] = 1.0 - epsInfo.Swl - epsInfo.Sgcr;
365 saturationKrnPoints_[1] = epsInfo.Sogcr;
366 saturationKrnPoints_[0] = 1.0 - epsInfo.Swl - epsInfo.Sgu;
369 maxPcnwOrLeverettFactor_ = epsInfo.pcgoLeverettFactor;
371 maxPcnwOrLeverettFactor_ = epsInfo.maxPcgo;
373 Krwr_ = epsInfo.Krorg;
374 Krnr_ = epsInfo.Krgr;
376 maxKrw_ = epsInfo.maxKrog;
377 maxKrn_ = epsInfo.maxKrg;
385 { saturationPcPoints_[pointIdx] = value; }
391 {
return saturationPcPoints_; }
397 { saturationKrwPoints_[pointIdx] = value; }
403 {
return saturationKrwPoints_; }
409 { saturationKrnPoints_[pointIdx] = value; }
415 {
return saturationKrnPoints_; }
421 { maxPcnwOrLeverettFactor_ = value; }
427 {
return maxPcnwOrLeverettFactor_; }
433 { maxPcnwOrLeverettFactor_ = value; }
439 {
return maxPcnwOrLeverettFactor_; }
446 { this->Krwr_ = value; }
453 {
return this->Krwr_; }
472 { this->Krnr_ = value; }
479 {
return this->Krnr_; }
495 std::cout <<
" saturationKrnPoints_[0]: " << saturationKrnPoints_[0] <<
"\n"
496 <<
" saturationKrnPoints_[1]: " << saturationKrnPoints_[1] <<
"\n"
497 <<
" saturationKrnPoints_[2]: " << saturationKrnPoints_[2] <<
"\n";
502 Scalar maxPcnwOrLeverettFactor_;
519 std::array<Scalar, 3> saturationPcPoints_;
522 std::array<Scalar, 3> saturationKrwPoints_;
525 std::array<Scalar, 3> saturationKrnPoints_;
Specifies the configuration used by the endpoint scaling code.
EclTwoPhaseSystemType
Specified which fluids are involved in a given twophase material law for endpoint scaling.
Definition: EclEpsConfig.hpp:49
Implements some common averages.
Scalar arithmeticMean(Scalar x, Scalar y)
Computes the arithmetic average of two values.
Definition: Means.hpp:45
Specifies the configuration used by the endpoint scaling code.
Definition: EclEpsConfig.hpp:63
bool enableLeverettScaling() const
Returns whether the Leverett capillary pressure scaling is enabled.
Definition: EclEpsConfig.hpp:183
Definition: EclEpsGridProperties.hpp:76
Represents the points on the X and Y axis to be scaled if endpoint scaling is used.
Definition: EclEpsScalingPoints.hpp:308
const std::array< Scalar, 3 > & saturationPcPoints() const
Returns the points used for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:390
void setMaxKrn(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:484
Scalar krnr() const
Returns non-wetting phase relative permeability at residual saturation of wetting phase.
Definition: EclEpsScalingPoints.hpp:478
Scalar maxKrn() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:490
void setKrnr(Scalar value)
Set non-wetting phase relative permeability at residual saturation of wetting phase.
Definition: EclEpsScalingPoints.hpp:471
const std::array< Scalar, 3 > & saturationKrnPoints() const
Returns the points used for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:414
void setMaxKrw(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:458
Scalar maxKrw() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:464
void setSaturationKrwPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for wetting-phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:396
void setMaxPcnw(Scalar value)
Sets the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:420
Scalar krwr() const
Returns wetting-phase relative permeability at residual saturation of non-wetting phase.
Definition: EclEpsScalingPoints.hpp:452
Scalar maxPcnw() const
Returns the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:426
void setKrwr(Scalar value)
Set wetting-phase relative permeability at residual saturation of non-wetting phase.
Definition: EclEpsScalingPoints.hpp:445
void setSaturationPcPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:384
void setLeverettFactor(Scalar value)
Sets the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:432
const std::array< Scalar, 3 > & saturationKrwPoints() const
Returns the points used for wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:402
void init(const EclEpsScalingPointsInfo< Scalar > &epsInfo, const EclEpsConfig &config, EclTwoPhaseSystemType epsSystemType)
Assigns the scaling points which actually ought to be used.
Definition: EclEpsScalingPoints.hpp:313
void setSaturationKrnPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:408
Scalar leverettFactor() const
Returns the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:438
This structure represents all values which can be possibly used as scaling points in the endpoint sca...
Definition: EclEpsScalingPoints.hpp:62