All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Pages
RelpermDiagnostics_impl.hpp
1 /*
2  Copyright 2016 Statoil ASA.
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 3 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 
20 #ifndef OPM_RELPERMDIAGNOSTICS_IMPL_HEADER_INCLUDED
21 #define OPM_RELPERMDIAGNOSTICS_IMPL_HEADER_INCLUDED
22 
23 #include <vector>
24 #include <utility>
25 
26 #include <opm/core/props/satfunc/RelpermDiagnostics.hpp>
27 #include <opm/core/utility/compressedToCartesian.hpp>
28 
29 namespace Opm {
30 
31  template <class GridT>
32  void RelpermDiagnostics::diagnosis(const Opm::EclipseState& eclState,
33  const Opm::Deck& deck,
34  const GridT& grid)
35  {
36  OpmLog::info("\n===============Saturation Functions Diagnostics===============\n");
37  phaseCheck_(eclState);
38  satFamilyCheck_(eclState);
39  tableCheck_(eclState);
40  unscaledEndPointsCheck_(deck, eclState);
41  scaledEndPointsCheck_(deck, eclState, grid);
42  }
43 
44  template <class GridT>
45  void RelpermDiagnostics::scaledEndPointsCheck_(const Deck& deck,
46  const EclipseState& eclState,
47  const GridT& grid)
48  {
49  // All end points are subject to round-off errors, checks should account for it
50  const float tolerance = 1e-6;
51  const int nc = Opm::UgGridHelpers::numCells(grid);
52  const auto& global_cell = Opm::UgGridHelpers::globalCell(grid);
53  const auto dims = Opm::UgGridHelpers::cartDims(grid);
54  const auto& compressedToCartesianIdx = Opm::compressedToCartesian(nc, global_cell);
55  scaledEpsInfo_.resize(nc);
56  EclEpsGridProperties epsGridProperties;
57  epsGridProperties.initFromDeck(deck, eclState, /*imbibition=*/false);
58  const auto& satnum = eclState.get3DProperties().getIntGridProperty("SATNUM");
59 
60  const std::string tag = "Scaled endpoints";
61  for (int c = 0; c < nc; ++c) {
62  const int cartIdx = compressedToCartesianIdx[c];
63  const std::string satnumIdx = std::to_string(satnum.iget(cartIdx));
64  std::array<int, 3> ijk;
65  ijk[0] = cartIdx % dims[0];
66  ijk[1] = (cartIdx / dims[0]) % dims[1];
67  ijk[2] = cartIdx / dims[0] / dims[1];
68  const std::string cellIdx = "(" + std::to_string(ijk[0]) + ", " +
69  std::to_string(ijk[1]) + ", " +
70  std::to_string(ijk[2]) + ")";
71  scaledEpsInfo_[c].extractScaled(eclState, epsGridProperties, cartIdx);
72 
73  // SGU <= 1.0 - SWL
74  if (scaledEpsInfo_[c].Sgu > (1.0 - scaledEpsInfo_[c].Swl + tolerance)) {
75  const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGU exceed 1.0 - SWL";
76  OpmLog::warning(tag, msg);
77  }
78 
79  // SGL <= 1.0 - SWU
80  if (scaledEpsInfo_[c].Sgl > (1.0 - scaledEpsInfo_[c].Swu + tolerance)) {
81  const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SGL exceed 1.0 - SWU";
82  OpmLog::warning(tag, msg);
83  }
84 
85  if (deck.hasKeyword("SCALECRS") && fluidSystem_ == FluidSystem::BlackOil) {
86  // Mobilility check.
87  if ((scaledEpsInfo_[c].Sowcr + scaledEpsInfo_[c].Swcr) >= (1.0 + tolerance)) {
88  const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOWCR + SWCR exceed 1.0";
89  OpmLog::warning(tag, msg);
90  }
91 
92  if ((scaledEpsInfo_[c].Sogcr + scaledEpsInfo_[c].Sgcr + scaledEpsInfo_[c].Swl) >= (1.0 + tolerance)) {
93  const std::string msg = "For scaled endpoints input, cell" + cellIdx + " SATNUM = " + satnumIdx + ", SOGCR + SGCR + SWL exceed 1.0";
94  OpmLog::warning(tag, msg);
95  }
96  }
97  }
98  }
99 
100 } //namespace Opm
101 
102 #endif // OPM_RELPERMDIAGNOSTICS_IMPL_HEADER_INCLUDED
void diagnosis(const EclipseState &eclState, const Deck &deck, const GridT &grid)
This function is used to diagnosis relperm in eclipse data file.