EclEpsScalingPoints.hpp
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef OPM_ECL_EPS_SCALING_POINTS_HPP
28 #define OPM_ECL_EPS_SCALING_POINTS_HPP
29 
30 #include "EclEpsConfig.hpp"
31 
32 #if HAVE_OPM_PARSER
33 #include <opm/parser/eclipse/Deck/Deck.hpp>
34 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
35 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/parser/eclipse/EclipseState/Grid/GridProperty.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/SgfnTable.hpp>
38 #include <opm/parser/eclipse/EclipseState/Tables/SgofTable.hpp>
39 #include <opm/parser/eclipse/EclipseState/Tables/SlgofTable.hpp>
40 #include <opm/parser/eclipse/EclipseState/Tables/Sof3Table.hpp>
41 #include <opm/parser/eclipse/EclipseState/Tables/SwfnTable.hpp>
42 #include <opm/parser/eclipse/EclipseState/Tables/SwofTable.hpp>
43 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
44 #endif
45 
47 
48 #include <array>
49 #include <string>
50 #include <iostream>
51 #include <cassert>
52 #include <algorithm>
53 
54 namespace Opm {
62 {
63  typedef std::vector<int> IntData;
64  typedef std::vector<double> DoubleData;
65 
66 public:
67 #if HAVE_OPM_PARSER
68  void initFromDeck(const Opm::Deck& /* deck */,
69  const Opm::EclipseState& eclState,
70  bool useImbibition)
71  {
72  std::string kwPrefix = useImbibition?"I":"";
73 
74  if (useImbibition)
75  satnum = &eclState.get3DProperties().getIntGridProperty("IMBNUM").getData();
76  else
77  satnum = &eclState.get3DProperties().getIntGridProperty("SATNUM").getData();
78 
79  retrieveGridPropertyData_(&swl, eclState, kwPrefix+"SWL");
80  retrieveGridPropertyData_(&sgl, eclState, kwPrefix+"SGL");
81  retrieveGridPropertyData_(&swcr, eclState, kwPrefix+"SWCR");
82  retrieveGridPropertyData_(&sgcr, eclState, kwPrefix+"SGCR");
83  retrieveGridPropertyData_(&sowcr, eclState, kwPrefix+"SOWCR");
84  retrieveGridPropertyData_(&sogcr, eclState, kwPrefix+"SOGCR");
85  retrieveGridPropertyData_(&swu, eclState, kwPrefix+"SWU");
86  retrieveGridPropertyData_(&sgu, eclState, kwPrefix+"SGU");
87  retrieveGridPropertyData_(&pcw, eclState, kwPrefix+"PCW");
88  retrieveGridPropertyData_(&pcg, eclState, kwPrefix+"PCG");
89  retrieveGridPropertyData_(&krw, eclState, kwPrefix+"KRW");
90  retrieveGridPropertyData_(&kro, eclState, kwPrefix+"KRO");
91  retrieveGridPropertyData_(&krg, eclState, kwPrefix+"KRG");
92 
93  // _may_ be needed to calculate the Leverett capillary pressure scaling factor
94  const auto& ecl3dProps = eclState.get3DProperties();
95  poro = &ecl3dProps.getDoubleGridProperty("PORO").getData();
96 
97  if (ecl3dProps.hasDeckDoubleGridProperty("PERMX")) {
98  permx = &ecl3dProps.getDoubleGridProperty("PERMX").getData();
99  permy = permx;
100  permz = permx;
101  }
102 
103  if (ecl3dProps.hasDeckDoubleGridProperty("PERMY"))
104  permy = &ecl3dProps.getDoubleGridProperty("PERMY").getData();
105 
106  if (ecl3dProps.hasDeckDoubleGridProperty("PERMZ"))
107  permz = &ecl3dProps.getDoubleGridProperty("PERMZ").getData();
108  }
109 #endif
110 
111  const IntData* satnum;
112 
113  const DoubleData* swl;
114  const DoubleData* sgl;
115  const DoubleData* swcr;
116  const DoubleData* sgcr;
117  const DoubleData* sowcr;
118  const DoubleData* sogcr;
119  const DoubleData* swu;
120  const DoubleData* sgu;
121  const DoubleData* pcw;
122  const DoubleData* pcg;
123  const DoubleData* krw;
124  const DoubleData* kro;
125  const DoubleData* krg;
126  const DoubleData* poro;
127  const DoubleData* permx;
128  const DoubleData* permy;
129  const DoubleData* permz;
130 
131 private:
132 #if HAVE_OPM_PARSER
133  // this method makes sure that a grid property is not created if it is not explicitly
134  // mentioned in the deck. (saves memory.)
135  void retrieveGridPropertyData_(const DoubleData **data,
136  const Opm::EclipseState& eclState,
137  const std::string& properyName)
138  {
139  (*data) = 0;
140  if (eclState.get3DProperties().hasDeckDoubleGridProperty(properyName))
141  (*data) = &eclState.get3DProperties().getDoubleGridProperty(properyName).getData();
142  }
143 #endif
144 };
145 
153 template <class Scalar>
155 {
156  // connate saturations
157  Scalar Swl; // oil
158  Scalar Sgl; // gas
159  Scalar Sowl; // oil for the oil-water system
160  Scalar Sogl; // oil for the gas-oil system
161 
162  // critical water and gas saturations
163  Scalar Swcr; // oil
164  Scalar Sgcr; // gas
165  Scalar Sowcr; // oil for the oil-water system
166  Scalar Sogcr; // oil for the gas-oil system
167 
168  // maximum saturations
169  Scalar Swu; // oil
170  Scalar Sgu; // gas
171  Scalar Sowu; // oil for the oil-water system
172  Scalar Sogu; // oil for the gas-oil system
173 
174  // maximum capillary pressures
175  Scalar maxPcow; // maximum capillary pressure of the oil-water system
176  Scalar maxPcgo; // maximum capillary pressure of the gas-oil system
177 
178  // the Leverett capillary pressure scaling factors. (those only make sense for the
179  // scaled points, for the unscaled ones they are 1.0.)
180  Scalar pcowLeverettFactor;
181  Scalar pcgoLeverettFactor;
182 
183  // maximum relative permabilities
184  Scalar maxKrw; // maximum relative permability of water
185  Scalar maxKrow; // maximum relative permability of oil in the oil-water system
186  Scalar maxKrog; // maximum relative permability of oil in the gas-oil system
187  Scalar maxKrg; // maximum relative permability of gas
188 
189  void print() const
190  {
191  std::cout << " Swl: " << Swl << "\n"
192  << " Sgl: " << Sgl << "\n"
193  << " Sowl: " << Sowl << "\n"
194  << " Sogl: " << Sogl << "\n"
195  << " Swcr: " << Swcr << "\n"
196  << " Sgcr: " << Sgcr << "\n"
197  << " Sowcr: " << Sowcr << "\n"
198  << " Sogcr: " << Sogcr << "\n"
199  << " Swu: " << Swu << "\n"
200  << " Sgu: " << Sgu << "\n"
201  << " Sowu: " << Sowu << "\n"
202  << " Sogu: " << Sogu << "\n"
203  << " maxPcow: " << maxPcow << "\n"
204  << " maxPcgo: " << maxPcgo << "\n"
205  << " pcowLeverettFactor: " << pcowLeverettFactor << "\n"
206  << " pcgoLeverettFactor: " << pcgoLeverettFactor << "\n"
207  << " maxKrw: " << maxKrw << "\n"
208  << " maxKrg: " << maxKrg << "\n"
209  << " maxKrow: " << maxKrow << "\n"
210  << " maxKrog: " << maxKrog << "\n";
211  }
212 
213 #if HAVE_OPM_PARSER
214 
220  void extractUnscaled(const Opm::Deck& deck,
221  const Opm::EclipseState& eclState,
222  unsigned satRegionIdx)
223  {
224  // TODO: support for the SOF2/SOF3 keyword family
225  const auto& tables = eclState.getTableManager();
226  const TableContainer& swofTables = tables.getSwofTables();
227  const TableContainer& sgofTables = tables.getSgofTables();
228  const TableContainer& slgofTables = tables.getSlgofTables();
229  const TableContainer& swfnTables = tables.getSwfnTables();
230  const TableContainer& sgfnTables = tables.getSgfnTables();
231  const TableContainer& sof3Tables = tables.getSof3Tables();
232 
233  bool hasWater = deck.hasKeyword("WATER");
234  bool hasGas = deck.hasKeyword("GAS");
235  bool hasOil = deck.hasKeyword("OIL");
236 
237  if (!hasWater) {
238  Swl = 0.0;
239  Swu = 0.0;
240  Swcr = 0.0;
241  if (!sgofTables.empty())
242  extractUnscaledSgof_(sgofTables.getTable<SgofTable>(satRegionIdx));
243  else {
244  assert(!slgofTables.empty());
245  extractUnscaledSlgof_(slgofTables.getTable<SlgofTable>(satRegionIdx));
246  }
247  return;
248  }
249  else if (!hasGas) {
250  assert(!swofTables.empty());
251  Sgl = 0.0;
252  Sgu = 0.0;
253  Sgcr = 0.0;
254  extractUnscaledSwof_(swofTables.getTable<SwofTable>(satRegionIdx));
255  return;
256  }
257 
258  // so far, only water-oil and oil-gas simulations are supported, i.e.,
259  // there's no gas-water yet.
260  if (!hasWater || !hasGas || !hasOil)
261  throw std::domain_error("The specified phase configuration is not suppored");
262 
263  bool family1 = (!sgofTables.empty() || !slgofTables.empty()) && !swofTables.empty();
264  bool family2 = !swfnTables.empty() && !sgfnTables.empty() && !sof3Tables.empty();
265 
266  if (family1) {
267  extractUnscaledSwof_(swofTables.getTable<SwofTable>(satRegionIdx));
268 
269  if (!sgofTables.empty()) {
270  // gas-oil parameters are specified using the SGOF keyword
271  extractUnscaledSgof_(sgofTables.getTable<SgofTable>(satRegionIdx));
272  }
273  else {
274  // gas-oil parameters are specified using the SLGOF keyword
275  assert(!slgofTables.empty());
276 
277  extractUnscaledSlgof_(slgofTables.getTable<SlgofTable>(satRegionIdx));
278  }
279  }
280  else if (family2) {
281  extractUnscaledSwfn_(swfnTables.getTable<SwfnTable>(satRegionIdx));
282  extractUnscaledSgfn_(sgfnTables.getTable<SgfnTable>(satRegionIdx));
283  extractUnscaledSof3_(sof3Tables.getTable<Sof3Table>(satRegionIdx));
284  }
285  else {
286  throw std::domain_error("No valid saturation keyword family specified");
287  }
288 
289  // there are no "unscaled" Leverett factors, so we just set them to 1.0
290  pcowLeverettFactor = 1.0;
291  pcgoLeverettFactor = 1.0;
292  }
293 
299  void extractScaled(const Opm::EclipseState& eclState,
300  const EclEpsGridProperties& epsProperties,
301  unsigned cartesianCellIdx)
302  {
303  // overwrite the unscaled values with the values for the cell if it is
304  // explicitly specified by the corresponding keyword.
305  extractGridPropertyValue_(Swl, epsProperties.swl, cartesianCellIdx);
306  extractGridPropertyValue_(Sgl, epsProperties.sgl, cartesianCellIdx);
307  extractGridPropertyValue_(Swcr, epsProperties.swcr, cartesianCellIdx);
308  extractGridPropertyValue_(Sgcr, epsProperties.sgcr, cartesianCellIdx);
309  extractGridPropertyValue_(Sowcr, epsProperties.sowcr, cartesianCellIdx);
310  extractGridPropertyValue_(Sogcr, epsProperties.sogcr, cartesianCellIdx);
311  extractGridPropertyValue_(Swu, epsProperties.swu, cartesianCellIdx);
312  extractGridPropertyValue_(Sgu, epsProperties.sgu, cartesianCellIdx);
313 
314  extractGridPropertyValue_(maxPcow, epsProperties.pcw, cartesianCellIdx);
315  extractGridPropertyValue_(maxPcgo, epsProperties.pcg, cartesianCellIdx);
316 
317  extractGridPropertyValue_(maxKrw, epsProperties.krw, cartesianCellIdx);
318  extractGridPropertyValue_(maxKrg, epsProperties.krg, cartesianCellIdx);
319 
320  // quite likely that's wrong!
321  extractGridPropertyValue_(maxKrow, epsProperties.kro, cartesianCellIdx);
322  extractGridPropertyValue_(maxKrog, epsProperties.kro, cartesianCellIdx);
323 
324  // compute the Leverett capillary pressure scaling factors if applicable. note
325  // that this needs to be done using non-SI units to make it correspond to the
326  // documentation.
327  pcowLeverettFactor = 1.0;
328  pcgoLeverettFactor = 1.0;
329  if (eclState.getTableManager().useJFunc()) {
330  const auto& jfunc = eclState.getTableManager().getJFunc();
331  const auto& jfuncDir = jfunc.direction();
332 
333  Scalar perm;
334  if (jfuncDir == Opm::JFunc::Direction::X)
335  perm =
336  (*epsProperties.permx)[cartesianCellIdx];
337  else if (jfuncDir == Opm::JFunc::Direction::Y)
338  perm =
339  (*epsProperties.permy)[cartesianCellIdx];
340  else if (jfuncDir == Opm::JFunc::Direction::Z)
341  perm =
342  (*epsProperties.permz)[cartesianCellIdx];
343  else if (jfuncDir == Opm::JFunc::Direction::XY)
344  // TODO: verify that this really is the arithmetic mean. (the
345  // documentation just says that the "average" should be used, IMO the
346  // harmonic mean would be more appropriate because that's what's usually
347  // applied when calculating the fluxes.)
348  perm =
349  Opm::arithmeticMean((*epsProperties.permx)[cartesianCellIdx],
350  (*epsProperties.permy)[cartesianCellIdx]);
351  else
352  OPM_THROW(std::runtime_error, "Illegal direction indicator for the JFUNC "
353  "keyword ("<<static_cast<int>(jfuncDir)<<")");
354 
355  // convert permeability from m^2 to mD
356  perm *= 1.01325e15;
357 
358  Scalar poro = (*epsProperties.poro)[cartesianCellIdx];
359  Scalar alpha = jfunc.alphaFactor();
360  Scalar beta = jfunc.betaFactor();
361 
362  // the part of the Leverett capillary pressure which does not depend on
363  // surface tension.
364  Scalar commonFactor = std::pow(poro, alpha)/std::pow(perm, beta);
365 
366  // multiply the documented constant by 10^5 because we want the pressures
367  // in [Pa], not in [bar]
368  const Scalar Uconst = 0.318316 * 1e5;
369 
370  // compute the oil-water Leverett factor.
371  const auto& jfuncFlag = jfunc.flag();
372  if (jfuncFlag == Opm::JFunc::Flag::WATER || jfuncFlag == Opm::JFunc::Flag::BOTH) {
373  // note that we use the surface tension in terms of [dyn/cm]
374  Scalar gamma =
375  jfunc.owSurfaceTension();
376  pcowLeverettFactor = commonFactor*gamma*Uconst;
377  }
378 
379  // compute the gas-oil Leverett factor.
380  if (jfuncFlag == Opm::JFunc::Flag::GAS || jfuncFlag == Opm::JFunc::Flag::BOTH) {
381  // note that we use the surface tension in terms of [dyn/cm]
382  Scalar gamma =
383  jfunc.goSurfaceTension();
384  pcgoLeverettFactor = commonFactor*gamma*Uconst;
385  }
386  }
387  }
388 #endif
389 
390 private:
391 #if HAVE_OPM_PARSER
392  void extractUnscaledSgof_(const Opm::SgofTable& sgofTable)
393  {
394  // minimum gas and oil-in-gas-oil saturation
395  Sgl = sgofTable.getSgColumn().front();
396  Sogl = 1.0 - sgofTable.getSgColumn().back();
397 
398  // maximum gas and oil-in-gas-oil saturation
399  Sgu = sgofTable.getSgColumn().back();
400  Sogu = 1.0 - sgofTable.getSgColumn().front();
401 
402  // critical gas saturation
403  Sgcr = 0.0;
404  for (size_t rowIdx = 0; rowIdx < sgofTable.numRows(); ++ rowIdx) {
405  if (sgofTable.getKrgColumn()[rowIdx] > 0.0)
406  break;
407 
408  Sgcr = sgofTable.getSgColumn()[rowIdx];
409  }
410 
411  // critical oil saturation of gas-oil system
412  Sogcr = 0.0;
413  for (int rowIdx = static_cast<int>(sgofTable.numRows() - 1);
414  rowIdx >= 0;
415  -- rowIdx)
416  {
417  if (sgofTable.getKrogColumn()[static_cast<size_t>(rowIdx)] > 0.0)
418  break;
419 
420  Sogcr = 1.0 - sgofTable.getSgColumn()[static_cast<size_t>(rowIdx)];
421  }
422 
423  // maximum gas-oil capillary pressure
424  maxPcgo = sgofTable.getPcogColumn().back();
425 
426  // maximum gas-* relperms
427  maxKrg = sgofTable.getKrgColumn().back();
428  maxKrog = sgofTable.getKrogColumn().front();
429  }
430 
431  void extractUnscaledSlgof_(const Opm::SlgofTable& slgofTable)
432  {
433  // minimum gas and oil-in-gas-oil saturation
434  Sgl = 1.0 - slgofTable.getSlColumn().back();
435  Sogl = slgofTable.getSlColumn().front();
436 
437  // maximum gas and oil-in-gas-oil saturation
438  Sgu = 1.0 - slgofTable.getSlColumn().front();
439  Sogu = slgofTable.getSlColumn().back();
440 
441  // critical gas saturation
442  Sgcr = 0.0;
443  for (int rowIdx = static_cast<int>(slgofTable.numRows()) - 1;
444  rowIdx >= 0;
445  -- rowIdx)
446  {
447  if (slgofTable.getKrgColumn()[static_cast<size_t>(rowIdx)] > 0.0)
448  break;
449 
450  Sgcr = 1 - slgofTable.getSlColumn()[static_cast<size_t>(rowIdx)];
451  }
452 
453  // critical oil saturation of gas-oil system
454  Sogcr = 0.0;
455  for (size_t rowIdx = 0; rowIdx < slgofTable.numRows(); ++ rowIdx) {
456  if (slgofTable.getKrogColumn()[rowIdx] > 0.0)
457  break;
458 
459  Sogcr = slgofTable.getSlColumn()[rowIdx];
460  }
461 
462  // maximum gas-oil capillary pressure
463  maxPcgo = slgofTable.getPcogColumn().front();
464 
465  // maximum gas-* relperms
466  maxKrg = slgofTable.getKrgColumn().front();
467  maxKrog = slgofTable.getKrogColumn().back();
468  }
469 
470  void extractUnscaledSwof_(const Opm::SwofTable& swofTable)
471  {
472  // connate saturations
473  Swl = swofTable.getSwColumn().front();
474  Sowl = 1.0 - swofTable.getSwColumn().back();
475 
476  // maximum water and oil-in-oil-water saturations
477  Swu = swofTable.getSwColumn().back();
478  Sowu = 1.0 - swofTable.getSwColumn().front();
479 
480  // critical water saturation
481  Swcr = 0.0;
482  for (size_t rowIdx = 0; rowIdx < swofTable.numRows(); ++ rowIdx) {
483  if (swofTable.getKrwColumn()[rowIdx] > 0.0)
484  break;
485 
486  Swcr = swofTable.getSwColumn()[rowIdx];
487  }
488 
489  // critical oil saturation of oil-water system
490  Sowcr = 0.0;
491  for (int rowIdx = static_cast<int>(swofTable.numRows()) - 1;
492  rowIdx >= 0;
493  -- rowIdx)
494  {
495  if (swofTable.getKrowColumn()[static_cast<size_t>(rowIdx)] > 0.0)
496  break;
497 
498  Sowcr = 1.0 - swofTable.getSwColumn()[static_cast<size_t>(rowIdx)];
499  }
500 
501  // maximum oil-water capillary pressures
502  maxPcow = swofTable.getPcowColumn().front();
503 
504  // maximum water-* relative permeabilities
505  maxKrw = swofTable.getKrwColumn().back();
506  maxKrow = swofTable.getKrowColumn().front();
507  }
508 
509  void extractUnscaledSwfn_(const Opm::SwfnTable& swfnTable)
510  {
511  // connate water saturation
512  Swl = swfnTable.getSwColumn().front();
513 
514  // maximum water saturation
515  Swu = swfnTable.getSwColumn().back();
516 
517  // critical water saturation
518  Swcr = 0.0;
519  for (size_t rowIdx = 0; rowIdx < swfnTable.numRows(); ++ rowIdx) {
520  if (swfnTable.getKrwColumn()[rowIdx] > 0)
521  break;
522 
523  Swcr = swfnTable.getSwColumn()[rowIdx];
524  }
525 
526  // maximum oil-water capillary pressure
527  maxPcow = swfnTable.getPcowColumn().front();
528 
529  // maximum water relative permeability
530  maxKrw = swfnTable.getKrwColumn().back();
531  }
532 
533  void extractUnscaledSgfn_(const Opm::SgfnTable& sgfnTable)
534  {
535  // connate gas saturation
536  Sgl = sgfnTable.getSgColumn().front();
537 
538  // maximum gas saturations
539  Sgu = sgfnTable.getSgColumn().back();
540  Sogu = 1 - sgfnTable.getSgColumn().front();
541 
542  // critical gas saturation
543  Sgcr = 0.0;
544  for (size_t rowIdx = 0; rowIdx < sgfnTable.numRows(); ++ rowIdx) {
545  if (sgfnTable.getKrgColumn()[rowIdx] > 0)
546  break;
547 
548  Sgcr = sgfnTable.getSgColumn()[rowIdx];
549  }
550 
551  // maximum capillary pressure
552  maxPcgo = sgfnTable.getPcogColumn().back();
553 
554  // maximum relative gas permeability
555  maxKrg = sgfnTable.getKrgColumn().back();
556  }
557 
558  void extractUnscaledSof3_(const Opm::Sof3Table& sof3Table)
559  {
560  // connate oil saturations
561  Sowl = sof3Table.getSoColumn().front() + Sgl;
562  Sogl = sof3Table.getSoColumn().front() + Swl;
563 
564  // maximum oil saturations
565  Sowu = sof3Table.getSoColumn().back();
566 
567  // critical oil saturation of oil-water system
568  Sowcr = 0.0;
569  for (size_t rowIdx = 0 ; rowIdx < sof3Table.numRows(); ++ rowIdx) {
570  if (sof3Table.getKrowColumn()[rowIdx] > 0) {
571  break;
572  };
573 
574  Sowcr = sof3Table.getSoColumn()[rowIdx];
575  }
576 
577  // critical oil saturation of gas-oil system
578  Sogcr = 0.0;
579  for (size_t rowIdx = 0 ; rowIdx < sof3Table.numRows(); ++ rowIdx) {
580  if (sof3Table.getKrogColumn()[rowIdx] > 0)
581  break;
582 
583  Sogcr = sof3Table.getSoColumn()[rowIdx];
584  }
585 
586  // maximum relative oil permeabilities
587  maxKrow = sof3Table.getKrowColumn().back();
588  maxKrog = sof3Table.getKrogColumn().back();
589  }
590 #endif // HAVE_OPM_PARSER
591 
592  void extractGridPropertyValue_(Scalar& targetValue,
593  const std::vector<double>* propData,
594  unsigned cartesianCellIdx)
595  {
596  if (!propData)
597  return;
598 
599  targetValue = (*propData)[cartesianCellIdx];
600  }
601 };
602 
609 template <class Scalar>
611 {
612 public:
617  const EclEpsConfig& config,
618  EclTwoPhaseSystemType epsSystemType)
619  {
620  if (epsSystemType == EclOilWaterSystem) {
621  // saturation scaling for capillary pressure
622  saturationPcPoints_[0] = epsInfo.Swl;
623  saturationPcPoints_[1] = epsInfo.Swu;
624 
625  // krw saturation scaling endpoints
626  if (config.enableThreePointKrSatScaling()) {
627  saturationKrwPoints_[0] = epsInfo.Swcr;
628  saturationKrwPoints_[1] = 1.0 - epsInfo.Sowcr - epsInfo.Sgl;
629  saturationKrwPoints_[2] = epsInfo.Swu;
630  }
631  else {
632  saturationKrwPoints_[0] = epsInfo.Swcr;
633  saturationKrwPoints_[1] = epsInfo.Swu;
634  }
635 
636  // krn saturation scaling endpoints (with the non-wetting phase being oil).
637  // because opm-material specifies non-wetting phase relperms in terms of the
638  // wetting phase saturations, the code here uses 1 minus the values specified
639  // by the Eclipse TD and the order of the scaling points is reversed
640  if (config.enableThreePointKrSatScaling()) {
641  saturationKrnPoints_[2] = 1.0 - epsInfo.Sowcr;
642  saturationKrnPoints_[1] = epsInfo.Swcr + epsInfo.Sgl;
643  saturationKrnPoints_[0] = epsInfo.Swl + epsInfo.Sgl;
644  }
645  else {
646  saturationKrnPoints_[1] = 1 - epsInfo.Sowcr;
647  saturationKrnPoints_[0] = epsInfo.Swl + epsInfo.Sgl;
648  }
649 
650  if (config.enableLeverettScaling())
651  maxPcnwOrLeverettFactor_ = epsInfo.pcowLeverettFactor;
652  else
653  maxPcnwOrLeverettFactor_ = epsInfo.maxPcow;
654  maxKrw_ = epsInfo.maxKrw;
655  maxKrn_ = epsInfo.maxKrow;
656  }
657  else {
658  assert(epsSystemType == EclGasOilSystem);
659 
660  // saturation scaling for capillary pressure
661  saturationPcPoints_[0] = 1.0 - epsInfo.Sgu;
662  saturationPcPoints_[1] = 1.0 - epsInfo.Sgl;
663 
664  // krw saturation scaling endpoints
665  if (config.enableThreePointKrSatScaling()) {
666  saturationKrwPoints_[0] = epsInfo.Sogcr;
667  saturationKrwPoints_[1] = 1 - epsInfo.Sgcr - epsInfo.Swl;
668  saturationKrwPoints_[2] = 1 - epsInfo.Swl - epsInfo.Sgl;
669  }
670  else {
671  saturationKrwPoints_[0] = epsInfo.Sogcr;
672  saturationKrwPoints_[1] = 1 - epsInfo.Swl - epsInfo.Sgl;
673  }
674 
675  // krn saturation scaling endpoints (with the non-wetting phase being gas).
676  // because opm-material specifies non-wetting phase relperms in terms of the
677  // wetting phase saturations, the code here uses 1 minus the values specified
678  // by the Eclipse TD and the order of the scaling points is reversed
679  if (config.enableThreePointKrSatScaling()) {
680  saturationKrnPoints_[2] = 1.0 - epsInfo.Sgcr;
681  saturationKrnPoints_[1] = epsInfo.Sogcr + epsInfo.Swl;
682  saturationKrnPoints_[0] = 1.0 - epsInfo.Sgu;
683  }
684  else {
685  saturationKrnPoints_[1] = 1.0 - epsInfo.Sgcr;
686  saturationKrnPoints_[0] = 1.0 - epsInfo.Sgu;
687  }
688 
689  if (config.enableLeverettScaling())
690  maxPcnwOrLeverettFactor_ = epsInfo.pcgoLeverettFactor;
691  else
692  maxPcnwOrLeverettFactor_ = epsInfo.maxPcgo;
693 
694  maxKrw_ = epsInfo.maxKrog;
695  maxKrn_ = epsInfo.maxKrg;
696  }
697  }
698 
702  void setSaturationPcPoint(unsigned pointIdx, Scalar value)
703  { saturationPcPoints_[pointIdx] = value; }
704 
708  const std::array<Scalar, 2>& saturationPcPoints() const
709  { return saturationPcPoints_; }
710 
714  void setSaturationKrwPoint(unsigned pointIdx, Scalar value)
715  { saturationKrwPoints_[pointIdx] = value; }
716 
720  const std::array<Scalar, 3>& saturationKrwPoints() const
721  { return saturationKrwPoints_; }
722 
726  void setSaturationKrnPoint(unsigned pointIdx, Scalar value)
727  { saturationKrnPoints_[pointIdx] = value; }
728 
732  const std::array<Scalar, 3>& saturationKrnPoints() const
733  { return saturationKrnPoints_; }
734 
738  void setMaxPcnw(Scalar value)
739  { maxPcnwOrLeverettFactor_ = value; }
740 
744  Scalar maxPcnw() const
745  { return maxPcnwOrLeverettFactor_; }
746 
750  void setLeverettFactor(Scalar value)
751  { maxPcnwOrLeverettFactor_ = value; }
752 
756  Scalar leverettFactor() const
757  { return maxPcnwOrLeverettFactor_; }
758 
762  void setMaxKrw(Scalar value)
763  { maxKrw_ = value; }
764 
768  Scalar maxKrw() const
769  { return maxKrw_; }
770 
774  void setMaxKrn(Scalar value)
775  { maxKrn_ = value; }
776 
780  Scalar maxKrn() const
781  { return maxKrn_; }
782 
783  void print() const
784  {
785  std::cout << " saturationKrnPoints_[0]: " << saturationKrnPoints_[0] << "\n"
786  << " saturationKrnPoints_[1]: " << saturationKrnPoints_[1] << "\n"
787  << " saturationKrnPoints_[2]: " << saturationKrnPoints_[2] << "\n";
788  }
789 
790 private:
791  // The the points used for the "y-axis" scaling of capillary pressure
792  Scalar maxPcnwOrLeverettFactor_;
793 
794  // The the points used for the "y-axis" scaling of wetting phase relative permability
795  Scalar maxKrw_;
796 
797  // The the points used for the "y-axis" scaling of non-wetting phase relative permability
798  Scalar maxKrn_;
799 
800  // The the points used for saturation ("x-axis") scaling of capillary pressure
801  std::array<Scalar, 2> saturationPcPoints_;
802 
803  // The the points used for saturation ("x-axis") scaling of wetting phase relative permeability
804  std::array<Scalar, 3> saturationKrwPoints_;
805 
806  // The the points used for saturation ("x-axis") scaling of non-wetting phase relative permeability
807  std::array<Scalar, 3> saturationKrnPoints_;
808 };
809 
810 } // namespace Opm
811 
812 #endif
Implements some common averages.
void setMaxKrw(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:762
void setSaturationPcPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:702
void setMaxPcnw(Scalar value)
Sets the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:738
Scalar maxKrw() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:768
void setSaturationKrnPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:726
Collects all grid properties which are relevant for end point scaling.
Definition: EclEpsScalingPoints.hpp:61
void setSaturationKrwPoint(unsigned pointIdx, Scalar value)
Sets an saturation value for wetting-phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:714
Definition: Air_Mesitylene.hpp:33
Specifies the configuration used by the endpoint scaling code.
Definition: EclEpsConfig.hpp:63
Scalar leverettFactor() const
Returns the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:756
This structure represents all values which can be possibly used as scaling points in the endpoint sca...
Definition: EclEpsScalingPoints.hpp:154
Scalar maxKrn() const
Returns the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:780
Represents the points on the X and Y axis to be scaled if endpoint scaling is used.
Definition: EclEpsScalingPoints.hpp:610
void setLeverettFactor(Scalar value)
Sets the Leverett scaling factor for capillary pressure.
Definition: EclEpsScalingPoints.hpp:750
Scalar arithmeticMean(Scalar x, Scalar y)
Computes the arithmetic average of two values.
Definition: Means.hpp:45
bool enableLeverettScaling() const
Returns whether the Leverett capillary pressure scaling is enabled.
Definition: EclEpsConfig.hpp:156
const std::array< Scalar, 3 > & saturationKrwPoints() const
Returns the points used for wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:720
Specifies the configuration used by the endpoint scaling code.
void setMaxKrn(Scalar value)
Sets the maximum wetting phase relative permeability.
Definition: EclEpsScalingPoints.hpp:774
Scalar maxPcnw() const
Returns the maximum capillary pressure.
Definition: EclEpsScalingPoints.hpp:744
bool enableThreePointKrSatScaling() const
Returns whether three point saturation scaling is enabled for the relative permeabilities.
Definition: EclEpsConfig.hpp:100
const std::array< Scalar, 3 > & saturationKrnPoints() const
Returns the points used for non-wetting phase relperm saturation scaling.
Definition: EclEpsScalingPoints.hpp:732
EclTwoPhaseSystemType
Specified which fluids are involved in a given twophase material law for endpoint scaling...
Definition: EclEpsConfig.hpp:51
const std::array< Scalar, 2 > & saturationPcPoints() const
Returns the points used for capillary pressure saturation scaling.
Definition: EclEpsScalingPoints.hpp:708
void init(const EclEpsScalingPointsInfo< Scalar > &epsInfo, const EclEpsConfig &config, EclTwoPhaseSystemType epsSystemType)
Assigns the scaling points which actually ought to be used.
Definition: EclEpsScalingPoints.hpp:616