00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_ECL_EPS_SCALING_POINTS_HPP
00028 #define OPM_ECL_EPS_SCALING_POINTS_HPP
00029
00030 #include "EclEpsConfig.hpp"
00031
00032 #if HAVE_OPM_PARSER
00033 #include <opm/parser/eclipse/Deck/Deck.hpp>
00034 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
00035 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
00036 #include <opm/parser/eclipse/EclipseState/Grid/GridProperty.hpp>
00037 #include <opm/parser/eclipse/EclipseState/Tables/SgfnTable.hpp>
00038 #include <opm/parser/eclipse/EclipseState/Tables/SgofTable.hpp>
00039 #include <opm/parser/eclipse/EclipseState/Tables/SlgofTable.hpp>
00040 #include <opm/parser/eclipse/EclipseState/Tables/Sof3Table.hpp>
00041 #include <opm/parser/eclipse/EclipseState/Tables/SwfnTable.hpp>
00042 #include <opm/parser/eclipse/EclipseState/Tables/SwofTable.hpp>
00043 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
00044 #endif
00045
00046 #include <opm/material/common/Means.hpp>
00047
00048 #include <array>
00049 #include <string>
00050 #include <iostream>
00051 #include <cassert>
00052 #include <algorithm>
00053
00054 namespace Opm {
00061 class EclEpsGridProperties
00062 {
00063 typedef std::vector<int> IntData;
00064 typedef std::vector<double> DoubleData;
00065
00066 public:
00067 #if HAVE_OPM_PARSER
00068 void initFromDeck(const Opm::Deck& ,
00069 const Opm::EclipseState& eclState,
00070 bool useImbibition)
00071 {
00072 std::string kwPrefix = useImbibition?"I":"";
00073
00074 if (useImbibition)
00075 satnum = &eclState.get3DProperties().getIntGridProperty("IMBNUM").getData();
00076 else
00077 satnum = &eclState.get3DProperties().getIntGridProperty("SATNUM").getData();
00078
00079 retrieveGridPropertyData_(&swl, eclState, kwPrefix+"SWL");
00080 retrieveGridPropertyData_(&sgl, eclState, kwPrefix+"SGL");
00081 retrieveGridPropertyData_(&swcr, eclState, kwPrefix+"SWCR");
00082 retrieveGridPropertyData_(&sgcr, eclState, kwPrefix+"SGCR");
00083 retrieveGridPropertyData_(&sowcr, eclState, kwPrefix+"SOWCR");
00084 retrieveGridPropertyData_(&sogcr, eclState, kwPrefix+"SOGCR");
00085 retrieveGridPropertyData_(&swu, eclState, kwPrefix+"SWU");
00086 retrieveGridPropertyData_(&sgu, eclState, kwPrefix+"SGU");
00087 retrieveGridPropertyData_(&pcw, eclState, kwPrefix+"PCW");
00088 retrieveGridPropertyData_(&pcg, eclState, kwPrefix+"PCG");
00089 retrieveGridPropertyData_(&krw, eclState, kwPrefix+"KRW");
00090 retrieveGridPropertyData_(&kro, eclState, kwPrefix+"KRO");
00091 retrieveGridPropertyData_(&krg, eclState, kwPrefix+"KRG");
00092
00093
00094 const auto& ecl3dProps = eclState.get3DProperties();
00095 poro = &ecl3dProps.getDoubleGridProperty("PORO").getData();
00096
00097 if (ecl3dProps.hasDeckDoubleGridProperty("PERMX")) {
00098 permx = &ecl3dProps.getDoubleGridProperty("PERMX").getData();
00099 permy = permx;
00100 permz = permx;
00101 }
00102
00103 if (ecl3dProps.hasDeckDoubleGridProperty("PERMY"))
00104 permy = &ecl3dProps.getDoubleGridProperty("PERMY").getData();
00105
00106 if (ecl3dProps.hasDeckDoubleGridProperty("PERMZ"))
00107 permz = &ecl3dProps.getDoubleGridProperty("PERMZ").getData();
00108 }
00109 #endif
00110
00111 const IntData* satnum;
00112
00113 const DoubleData* swl;
00114 const DoubleData* sgl;
00115 const DoubleData* swcr;
00116 const DoubleData* sgcr;
00117 const DoubleData* sowcr;
00118 const DoubleData* sogcr;
00119 const DoubleData* swu;
00120 const DoubleData* sgu;
00121 const DoubleData* pcw;
00122 const DoubleData* pcg;
00123 const DoubleData* krw;
00124 const DoubleData* kro;
00125 const DoubleData* krg;
00126 const DoubleData* poro;
00127 const DoubleData* permx;
00128 const DoubleData* permy;
00129 const DoubleData* permz;
00130
00131 private:
00132 #if HAVE_OPM_PARSER
00133
00134
00135 void retrieveGridPropertyData_(const DoubleData **data,
00136 const Opm::EclipseState& eclState,
00137 const std::string& properyName)
00138 {
00139 (*data) = 0;
00140 if (eclState.get3DProperties().hasDeckDoubleGridProperty(properyName))
00141 (*data) = &eclState.get3DProperties().getDoubleGridProperty(properyName).getData();
00142 }
00143 #endif
00144 };
00145
00153 template <class Scalar>
00154 struct EclEpsScalingPointsInfo
00155 {
00156
00157 Scalar Swl;
00158 Scalar Sgl;
00159 Scalar Sowl;
00160 Scalar Sogl;
00161
00162
00163 Scalar Swcr;
00164 Scalar Sgcr;
00165 Scalar Sowcr;
00166 Scalar Sogcr;
00167
00168
00169 Scalar Swu;
00170 Scalar Sgu;
00171 Scalar Sowu;
00172 Scalar Sogu;
00173
00174
00175 Scalar maxPcow;
00176 Scalar maxPcgo;
00177
00178
00179
00180 Scalar pcowLeverettFactor;
00181 Scalar pcgoLeverettFactor;
00182
00183
00184 Scalar maxKrw;
00185 Scalar maxKrow;
00186 Scalar maxKrog;
00187 Scalar maxKrg;
00188
00189 void print() const
00190 {
00191 std::cout << " Swl: " << Swl << "\n"
00192 << " Sgl: " << Sgl << "\n"
00193 << " Sowl: " << Sowl << "\n"
00194 << " Sogl: " << Sogl << "\n"
00195 << " Swcr: " << Swcr << "\n"
00196 << " Sgcr: " << Sgcr << "\n"
00197 << " Sowcr: " << Sowcr << "\n"
00198 << " Sogcr: " << Sogcr << "\n"
00199 << " Swu: " << Swu << "\n"
00200 << " Sgu: " << Sgu << "\n"
00201 << " Sowu: " << Sowu << "\n"
00202 << " Sogu: " << Sogu << "\n"
00203 << " maxPcow: " << maxPcow << "\n"
00204 << " maxPcgo: " << maxPcgo << "\n"
00205 << " pcowLeverettFactor: " << pcowLeverettFactor << "\n"
00206 << " pcgoLeverettFactor: " << pcgoLeverettFactor << "\n"
00207 << " maxKrw: " << maxKrw << "\n"
00208 << " maxKrg: " << maxKrg << "\n"
00209 << " maxKrow: " << maxKrow << "\n"
00210 << " maxKrog: " << maxKrog << "\n";
00211 }
00212
00213 #if HAVE_OPM_PARSER
00214
00220 void extractUnscaled(const Opm::Deck& deck,
00221 const Opm::EclipseState& eclState,
00222 unsigned satRegionIdx)
00223 {
00224
00225 const auto& tables = eclState.getTableManager();
00226 const TableContainer& swofTables = tables.getSwofTables();
00227 const TableContainer& sgofTables = tables.getSgofTables();
00228 const TableContainer& slgofTables = tables.getSlgofTables();
00229 const TableContainer& swfnTables = tables.getSwfnTables();
00230 const TableContainer& sgfnTables = tables.getSgfnTables();
00231 const TableContainer& sof3Tables = tables.getSof3Tables();
00232
00233 bool hasWater = deck.hasKeyword("WATER");
00234 bool hasGas = deck.hasKeyword("GAS");
00235 bool hasOil = deck.hasKeyword("OIL");
00236
00237 if (!hasWater) {
00238 Swl = 0.0;
00239 Swu = 0.0;
00240 Swcr = 0.0;
00241 if (!sgofTables.empty())
00242 extractUnscaledSgof_(sgofTables.getTable<SgofTable>(satRegionIdx));
00243 else {
00244 assert(!slgofTables.empty());
00245 extractUnscaledSlgof_(slgofTables.getTable<SlgofTable>(satRegionIdx));
00246 }
00247 return;
00248 }
00249 else if (!hasGas) {
00250 assert(!swofTables.empty());
00251 Sgl = 0.0;
00252 Sgu = 0.0;
00253 Sgcr = 0.0;
00254 extractUnscaledSwof_(swofTables.getTable<SwofTable>(satRegionIdx));
00255 return;
00256 }
00257
00258
00259
00260 if (!hasWater || !hasGas || !hasOil)
00261 throw std::domain_error("The specified phase configuration is not suppored");
00262
00263 bool family1 = (!sgofTables.empty() || !slgofTables.empty()) && !swofTables.empty();
00264 bool family2 = !swfnTables.empty() && !sgfnTables.empty() && !sof3Tables.empty();
00265
00266 if (family1) {
00267 extractUnscaledSwof_(swofTables.getTable<SwofTable>(satRegionIdx));
00268
00269 if (!sgofTables.empty()) {
00270
00271 extractUnscaledSgof_(sgofTables.getTable<SgofTable>(satRegionIdx));
00272 }
00273 else {
00274
00275 assert(!slgofTables.empty());
00276
00277 extractUnscaledSlgof_(slgofTables.getTable<SlgofTable>(satRegionIdx));
00278 }
00279 }
00280 else if (family2) {
00281 extractUnscaledSwfn_(swfnTables.getTable<SwfnTable>(satRegionIdx));
00282 extractUnscaledSgfn_(sgfnTables.getTable<SgfnTable>(satRegionIdx));
00283 extractUnscaledSof3_(sof3Tables.getTable<Sof3Table>(satRegionIdx));
00284 }
00285 else {
00286 throw std::domain_error("No valid saturation keyword family specified");
00287 }
00288
00289
00290 pcowLeverettFactor = 1.0;
00291 pcgoLeverettFactor = 1.0;
00292 }
00293
00299 void extractScaled(const Opm::EclipseState& eclState,
00300 const EclEpsGridProperties& epsProperties,
00301 unsigned cartesianCellIdx)
00302 {
00303
00304
00305 extractGridPropertyValue_(Swl, epsProperties.swl, cartesianCellIdx);
00306 extractGridPropertyValue_(Sgl, epsProperties.sgl, cartesianCellIdx);
00307 extractGridPropertyValue_(Swcr, epsProperties.swcr, cartesianCellIdx);
00308 extractGridPropertyValue_(Sgcr, epsProperties.sgcr, cartesianCellIdx);
00309 extractGridPropertyValue_(Sowcr, epsProperties.sowcr, cartesianCellIdx);
00310 extractGridPropertyValue_(Sogcr, epsProperties.sogcr, cartesianCellIdx);
00311 extractGridPropertyValue_(Swu, epsProperties.swu, cartesianCellIdx);
00312 extractGridPropertyValue_(Sgu, epsProperties.sgu, cartesianCellIdx);
00313
00314 extractGridPropertyValue_(maxPcow, epsProperties.pcw, cartesianCellIdx);
00315 extractGridPropertyValue_(maxPcgo, epsProperties.pcg, cartesianCellIdx);
00316
00317 extractGridPropertyValue_(maxKrw, epsProperties.krw, cartesianCellIdx);
00318 extractGridPropertyValue_(maxKrg, epsProperties.krg, cartesianCellIdx);
00319
00320
00321 extractGridPropertyValue_(maxKrow, epsProperties.kro, cartesianCellIdx);
00322 extractGridPropertyValue_(maxKrog, epsProperties.kro, cartesianCellIdx);
00323
00324
00325
00326
00327 pcowLeverettFactor = 1.0;
00328 pcgoLeverettFactor = 1.0;
00329 if (eclState.getTableManager().useJFunc()) {
00330 const auto& jfunc = eclState.getTableManager().getJFunc();
00331 const auto& jfuncDir = jfunc.direction();
00332
00333 Scalar perm;
00334 if (jfuncDir == Opm::JFunc::Direction::X)
00335 perm =
00336 (*epsProperties.permx)[cartesianCellIdx];
00337 else if (jfuncDir == Opm::JFunc::Direction::Y)
00338 perm =
00339 (*epsProperties.permy)[cartesianCellIdx];
00340 else if (jfuncDir == Opm::JFunc::Direction::Z)
00341 perm =
00342 (*epsProperties.permz)[cartesianCellIdx];
00343 else if (jfuncDir == Opm::JFunc::Direction::XY)
00344
00345
00346
00347
00348 perm =
00349 Opm::arithmeticMean((*epsProperties.permx)[cartesianCellIdx],
00350 (*epsProperties.permy)[cartesianCellIdx]);
00351 else
00352 OPM_THROW(std::runtime_error, "Illegal direction indicator for the JFUNC "
00353 "keyword ("<<static_cast<int>(jfuncDir)<<")");
00354
00355
00356 perm *= 1.01325e15;
00357
00358 Scalar poro = (*epsProperties.poro)[cartesianCellIdx];
00359 Scalar alpha = jfunc.alphaFactor();
00360 Scalar beta = jfunc.betaFactor();
00361
00362
00363
00364 Scalar commonFactor = std::pow(poro, alpha)/std::pow(perm, beta);
00365
00366
00367
00368 const Scalar Uconst = 0.318316 * 1e5;
00369
00370
00371 const auto& jfuncFlag = jfunc.flag();
00372 if (jfuncFlag == Opm::JFunc::Flag::WATER || jfuncFlag == Opm::JFunc::Flag::BOTH) {
00373
00374 Scalar gamma =
00375 jfunc.owSurfaceTension();
00376 pcowLeverettFactor = commonFactor*gamma*Uconst;
00377 }
00378
00379
00380 if (jfuncFlag == Opm::JFunc::Flag::GAS || jfuncFlag == Opm::JFunc::Flag::BOTH) {
00381
00382 Scalar gamma =
00383 jfunc.goSurfaceTension();
00384 pcgoLeverettFactor = commonFactor*gamma*Uconst;
00385 }
00386 }
00387 }
00388 #endif
00389
00390 private:
00391 #if HAVE_OPM_PARSER
00392 void extractUnscaledSgof_(const Opm::SgofTable& sgofTable)
00393 {
00394
00395 Sgl = sgofTable.getSgColumn().front();
00396 Sogl = 1.0 - sgofTable.getSgColumn().back();
00397
00398
00399 Sgu = sgofTable.getSgColumn().back();
00400 Sogu = 1.0 - sgofTable.getSgColumn().front();
00401
00402
00403 Sgcr = 0.0;
00404 for (size_t rowIdx = 0; rowIdx < sgofTable.numRows(); ++ rowIdx) {
00405 if (sgofTable.getKrgColumn()[rowIdx] > 0.0)
00406 break;
00407
00408 Sgcr = sgofTable.getSgColumn()[rowIdx];
00409 }
00410
00411
00412 Sogcr = 0.0;
00413 for (int rowIdx = static_cast<int>(sgofTable.numRows() - 1);
00414 rowIdx >= 0;
00415 -- rowIdx)
00416 {
00417 if (sgofTable.getKrogColumn()[static_cast<size_t>(rowIdx)] > 0.0)
00418 break;
00419
00420 Sogcr = 1.0 - sgofTable.getSgColumn()[static_cast<size_t>(rowIdx)];
00421 }
00422
00423
00424 maxPcgo = sgofTable.getPcogColumn().back();
00425
00426
00427 maxKrg = sgofTable.getKrgColumn().back();
00428 maxKrog = sgofTable.getKrogColumn().front();
00429 }
00430
00431 void extractUnscaledSlgof_(const Opm::SlgofTable& slgofTable)
00432 {
00433
00434 Sgl = 1.0 - slgofTable.getSlColumn().back();
00435 Sogl = slgofTable.getSlColumn().front();
00436
00437
00438 Sgu = 1.0 - slgofTable.getSlColumn().front();
00439 Sogu = slgofTable.getSlColumn().back();
00440
00441
00442 Sgcr = 0.0;
00443 for (int rowIdx = static_cast<int>(slgofTable.numRows()) - 1;
00444 rowIdx >= 0;
00445 -- rowIdx)
00446 {
00447 if (slgofTable.getKrgColumn()[static_cast<size_t>(rowIdx)] > 0.0)
00448 break;
00449
00450 Sgcr = 1 - slgofTable.getSlColumn()[static_cast<size_t>(rowIdx)];
00451 }
00452
00453
00454 Sogcr = 0.0;
00455 for (size_t rowIdx = 0; rowIdx < slgofTable.numRows(); ++ rowIdx) {
00456 if (slgofTable.getKrogColumn()[rowIdx] > 0.0)
00457 break;
00458
00459 Sogcr = slgofTable.getSlColumn()[rowIdx];
00460 }
00461
00462
00463 maxPcgo = slgofTable.getPcogColumn().front();
00464
00465
00466 maxKrg = slgofTable.getKrgColumn().front();
00467 maxKrog = slgofTable.getKrogColumn().back();
00468 }
00469
00470 void extractUnscaledSwof_(const Opm::SwofTable& swofTable)
00471 {
00472
00473 Swl = swofTable.getSwColumn().front();
00474 Sowl = 1.0 - swofTable.getSwColumn().back();
00475
00476
00477 Swu = swofTable.getSwColumn().back();
00478 Sowu = 1.0 - swofTable.getSwColumn().front();
00479
00480
00481 Swcr = 0.0;
00482 for (size_t rowIdx = 0; rowIdx < swofTable.numRows(); ++ rowIdx) {
00483 if (swofTable.getKrwColumn()[rowIdx] > 0.0)
00484 break;
00485
00486 Swcr = swofTable.getSwColumn()[rowIdx];
00487 }
00488
00489
00490 Sowcr = 0.0;
00491 for (int rowIdx = static_cast<int>(swofTable.numRows()) - 1;
00492 rowIdx >= 0;
00493 -- rowIdx)
00494 {
00495 if (swofTable.getKrowColumn()[static_cast<size_t>(rowIdx)] > 0.0)
00496 break;
00497
00498 Sowcr = 1.0 - swofTable.getSwColumn()[static_cast<size_t>(rowIdx)];
00499 }
00500
00501
00502 maxPcow = swofTable.getPcowColumn().front();
00503
00504
00505 maxKrw = swofTable.getKrwColumn().back();
00506 maxKrow = swofTable.getKrowColumn().front();
00507 }
00508
00509 void extractUnscaledSwfn_(const Opm::SwfnTable& swfnTable)
00510 {
00511
00512 Swl = swfnTable.getSwColumn().front();
00513
00514
00515 Swu = swfnTable.getSwColumn().back();
00516
00517
00518 Swcr = 0.0;
00519 for (size_t rowIdx = 0; rowIdx < swfnTable.numRows(); ++ rowIdx) {
00520 if (swfnTable.getKrwColumn()[rowIdx] > 0)
00521 break;
00522
00523 Swcr = swfnTable.getSwColumn()[rowIdx];
00524 }
00525
00526
00527 maxPcow = swfnTable.getPcowColumn().front();
00528
00529
00530 maxKrw = swfnTable.getKrwColumn().back();
00531 }
00532
00533 void extractUnscaledSgfn_(const Opm::SgfnTable& sgfnTable)
00534 {
00535
00536 Sgl = sgfnTable.getSgColumn().front();
00537
00538
00539 Sgu = sgfnTable.getSgColumn().back();
00540 Sogu = 1 - sgfnTable.getSgColumn().front();
00541
00542
00543 Sgcr = 0.0;
00544 for (size_t rowIdx = 0; rowIdx < sgfnTable.numRows(); ++ rowIdx) {
00545 if (sgfnTable.getKrgColumn()[rowIdx] > 0)
00546 break;
00547
00548 Sgcr = sgfnTable.getSgColumn()[rowIdx];
00549 }
00550
00551
00552 maxPcgo = sgfnTable.getPcogColumn().back();
00553
00554
00555 maxKrg = sgfnTable.getKrgColumn().back();
00556 }
00557
00558 void extractUnscaledSof3_(const Opm::Sof3Table& sof3Table)
00559 {
00560
00561 Sowl = sof3Table.getSoColumn().front() + Sgl;
00562 Sogl = sof3Table.getSoColumn().front() + Swl;
00563
00564
00565 Sowu = sof3Table.getSoColumn().back();
00566
00567
00568 Sowcr = 0.0;
00569 for (size_t rowIdx = 0 ; rowIdx < sof3Table.numRows(); ++ rowIdx) {
00570 if (sof3Table.getKrowColumn()[rowIdx] > 0) {
00571 break;
00572 };
00573
00574 Sowcr = sof3Table.getSoColumn()[rowIdx];
00575 }
00576
00577
00578 Sogcr = 0.0;
00579 for (size_t rowIdx = 0 ; rowIdx < sof3Table.numRows(); ++ rowIdx) {
00580 if (sof3Table.getKrogColumn()[rowIdx] > 0)
00581 break;
00582
00583 Sogcr = sof3Table.getSoColumn()[rowIdx];
00584 }
00585
00586
00587 maxKrow = sof3Table.getKrowColumn().back();
00588 maxKrog = sof3Table.getKrogColumn().back();
00589 }
00590 #endif // HAVE_OPM_PARSER
00591
00592 void extractGridPropertyValue_(Scalar& targetValue,
00593 const std::vector<double>* propData,
00594 unsigned cartesianCellIdx)
00595 {
00596 if (!propData)
00597 return;
00598
00599 targetValue = (*propData)[cartesianCellIdx];
00600 }
00601 };
00602
00609 template <class Scalar>
00610 class EclEpsScalingPoints
00611 {
00612 public:
00616 void init(const EclEpsScalingPointsInfo<Scalar>& epsInfo,
00617 const EclEpsConfig& config,
00618 EclTwoPhaseSystemType epsSystemType)
00619 {
00620 if (epsSystemType == EclOilWaterSystem) {
00621
00622 saturationPcPoints_[0] = epsInfo.Swl;
00623 saturationPcPoints_[1] = epsInfo.Swu;
00624
00625
00626 if (config.enableThreePointKrSatScaling()) {
00627 saturationKrwPoints_[0] = epsInfo.Swcr;
00628 saturationKrwPoints_[1] = 1.0 - epsInfo.Sowcr - epsInfo.Sgl;
00629 saturationKrwPoints_[2] = epsInfo.Swu;
00630 }
00631 else {
00632 saturationKrwPoints_[0] = epsInfo.Swcr;
00633 saturationKrwPoints_[1] = epsInfo.Swu;
00634 }
00635
00636
00637
00638
00639
00640 if (config.enableThreePointKrSatScaling()) {
00641 saturationKrnPoints_[2] = 1.0 - epsInfo.Sowcr;
00642 saturationKrnPoints_[1] = epsInfo.Swcr + epsInfo.Sgl;
00643 saturationKrnPoints_[0] = epsInfo.Swl + epsInfo.Sgl;
00644 }
00645 else {
00646 saturationKrnPoints_[1] = 1 - epsInfo.Sowcr;
00647 saturationKrnPoints_[0] = epsInfo.Swl + epsInfo.Sgl;
00648 }
00649
00650 if (config.enableLeverettScaling())
00651 maxPcnwOrLeverettFactor_ = epsInfo.pcowLeverettFactor;
00652 else
00653 maxPcnwOrLeverettFactor_ = epsInfo.maxPcow;
00654 maxKrw_ = epsInfo.maxKrw;
00655 maxKrn_ = epsInfo.maxKrow;
00656 }
00657 else {
00658 assert(epsSystemType == EclGasOilSystem);
00659
00660
00661 saturationPcPoints_[0] = 1.0 - epsInfo.Sgu;
00662 saturationPcPoints_[1] = 1.0 - epsInfo.Sgl;
00663
00664
00665 if (config.enableThreePointKrSatScaling()) {
00666 saturationKrwPoints_[0] = epsInfo.Sogcr;
00667 saturationKrwPoints_[1] = 1 - epsInfo.Sgcr - epsInfo.Swl;
00668 saturationKrwPoints_[2] = 1 - epsInfo.Swl - epsInfo.Sgl;
00669 }
00670 else {
00671 saturationKrwPoints_[0] = epsInfo.Sogcr;
00672 saturationKrwPoints_[1] = 1 - epsInfo.Swl - epsInfo.Sgl;
00673 }
00674
00675
00676
00677
00678
00679 if (config.enableThreePointKrSatScaling()) {
00680 saturationKrnPoints_[2] = 1.0 - epsInfo.Sgcr;
00681 saturationKrnPoints_[1] = epsInfo.Sogcr + epsInfo.Swl;
00682 saturationKrnPoints_[0] = 1.0 - epsInfo.Sgu;
00683 }
00684 else {
00685 saturationKrnPoints_[1] = 1.0 - epsInfo.Sgcr;
00686 saturationKrnPoints_[0] = 1.0 - epsInfo.Sgu;
00687 }
00688
00689 if (config.enableLeverettScaling())
00690 maxPcnwOrLeverettFactor_ = epsInfo.pcgoLeverettFactor;
00691 else
00692 maxPcnwOrLeverettFactor_ = epsInfo.maxPcgo;
00693
00694 maxKrw_ = epsInfo.maxKrog;
00695 maxKrn_ = epsInfo.maxKrg;
00696 }
00697 }
00698
00702 void setSaturationPcPoint(unsigned pointIdx, Scalar value)
00703 { saturationPcPoints_[pointIdx] = value; }
00704
00708 const std::array<Scalar, 2>& saturationPcPoints() const
00709 { return saturationPcPoints_; }
00710
00714 void setSaturationKrwPoint(unsigned pointIdx, Scalar value)
00715 { saturationKrwPoints_[pointIdx] = value; }
00716
00720 const std::array<Scalar, 3>& saturationKrwPoints() const
00721 { return saturationKrwPoints_; }
00722
00726 void setSaturationKrnPoint(unsigned pointIdx, Scalar value)
00727 { saturationKrnPoints_[pointIdx] = value; }
00728
00732 const std::array<Scalar, 3>& saturationKrnPoints() const
00733 { return saturationKrnPoints_; }
00734
00738 void setMaxPcnw(Scalar value)
00739 { maxPcnwOrLeverettFactor_ = value; }
00740
00744 Scalar maxPcnw() const
00745 { return maxPcnwOrLeverettFactor_; }
00746
00750 void setLeverettFactor(Scalar value)
00751 { maxPcnwOrLeverettFactor_ = value; }
00752
00756 Scalar leverettFactor() const
00757 { return maxPcnwOrLeverettFactor_; }
00758
00762 void setMaxKrw(Scalar value)
00763 { maxKrw_ = value; }
00764
00768 Scalar maxKrw() const
00769 { return maxKrw_; }
00770
00774 void setMaxKrn(Scalar value)
00775 { maxKrn_ = value; }
00776
00780 Scalar maxKrn() const
00781 { return maxKrn_; }
00782
00783 void print() const
00784 {
00785 std::cout << " saturationKrnPoints_[0]: " << saturationKrnPoints_[0] << "\n"
00786 << " saturationKrnPoints_[1]: " << saturationKrnPoints_[1] << "\n"
00787 << " saturationKrnPoints_[2]: " << saturationKrnPoints_[2] << "\n";
00788 }
00789
00790 private:
00791
00792 Scalar maxPcnwOrLeverettFactor_;
00793
00794
00795 Scalar maxKrw_;
00796
00797
00798 Scalar maxKrn_;
00799
00800
00801 std::array<Scalar, 2> saturationPcPoints_;
00802
00803
00804 std::array<Scalar, 3> saturationKrwPoints_;
00805
00806
00807 std::array<Scalar, 3> saturationKrnPoints_;
00808 };
00809
00810 }
00811
00812 #endif