phaseUsageFromDeck.hpp
1 /*
2  Copyright 2012 SINTEF ICT, Applied Mathematics.
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_PHASEUSAGEFROMDECK_HEADER_INCLUDED
21 #define OPM_PHASEUSAGEFROMDECK_HEADER_INCLUDED
22 
23 #include <opm/core/props/BlackoilPhases.hpp>
24 #include <opm/common/ErrorMacros.hpp>
25 
26 #include <opm/parser/eclipse/Deck/Deck.hpp>
27 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
28 #include <opm/parser/eclipse/EclipseState/Runspec.hpp>
29 
30 
31 namespace Opm
32 {
33 
36  inline PhaseUsage phaseUsageFromDeck(const Opm::EclipseState& eclipseState)
37  {
38  PhaseUsage pu;
39  std::fill(pu.phase_used, pu.phase_used + BlackoilPhases::MaxNumPhases, 0);
40 
41  const auto& phase = eclipseState.runspec().phases();
42  // Discover phase usage.
43  if (phase.active(Phase::WATER)) {
44  pu.phase_used[BlackoilPhases::Aqua] = 1;
45  }
46  if (phase.active(Phase::OIL)) {
47  pu.phase_used[BlackoilPhases::Liquid] = 1;
48  }
49  if (phase.active(Phase::GAS)) {
50  pu.phase_used[BlackoilPhases::Vapour] = 1;
51  }
52  pu.num_phases = 0;
53  for (int i = 0; i < BlackoilPhases::MaxNumPhases; ++i) {
54  pu.phase_pos[i] = pu.num_phases;
55  pu.num_phases += pu.phase_used[i];
56  }
57 
58  // Only 2 or 3 phase systems handled.
59  if (pu.num_phases < 2 || pu.num_phases > 3) {
60  OPM_THROW(std::runtime_error, "Cannot handle cases with " << pu.num_phases << " phases.");
61  }
62 
63  // We need oil systems, since we do not support the keywords needed for
64  // water-gas systems.
65  if (!pu.phase_used[BlackoilPhases::Liquid]) {
66  OPM_THROW(std::runtime_error, "Cannot handle cases with no OIL, i.e. water-gas systems.");
67  }
68 
69  // Add solvent info
70  pu.has_solvent = false;
71  if (phase.active(Phase::SOLVENT)) {
72  pu.has_solvent = true;
73  }
74 
75  // Add polyme info
76  pu.has_polymer = false;
77  if (phase.active(Phase::POLYMER)) {
78  pu.has_polymer = true;
79  }
80 
81  return pu;
82  }
83 
86  inline PhaseUsage phaseUsageFromDeck(const Opm::Deck& deck)
87  {
88  PhaseUsage pu;
89  std::fill(pu.phase_used, pu.phase_used + BlackoilPhases::MaxNumPhases, 0);
90 
91  Runspec runspec( deck );
92  const auto& phase = runspec.phases();
93 
94  // Discover phase usage.
95  if (phase.active( Phase::WATER )) {
96  pu.phase_used[BlackoilPhases::Aqua] = 1;
97  }
98  if (phase.active( Phase::OIL )) {
99  pu.phase_used[BlackoilPhases::Liquid] = 1;
100  }
101  if (phase.active( Phase::GAS )) {
102  pu.phase_used[BlackoilPhases::Vapour] = 1;
103  }
104  pu.num_phases = 0;
105  for (int i = 0; i < BlackoilPhases::MaxNumPhases; ++i) {
106  if (pu.phase_used[i]) {
107  pu.phase_pos[i] = pu.num_phases;
108  pu.num_phases += 1;
109  }
110  else {
111  //Set to ridiculous value on purpose: should never be used
112  pu.phase_pos[i] = 2000000000;
113  }
114  }
115 
116  // Only 2 or 3 phase systems handled.
117  if (pu.num_phases < 2 || pu.num_phases > 3) {
118  OPM_THROW(std::runtime_error, "Cannot handle cases with " << pu.num_phases << " phases.");
119  }
120 
121  // We need oil systems, since we do not support the keywords needed for
122  // water-gas systems.
123  if (!pu.phase_used[BlackoilPhases::Liquid]) {
124  OPM_THROW(std::runtime_error, "Cannot handle cases with no OIL, i.e. water-gas systems.");
125  }
126 
127  // Add solvent info
128  pu.has_solvent = false;
129  if (phase.active(Phase::SOLVENT)) {
130  pu.has_solvent = true;
131  }
132 
133  // Add polyme info
134  pu.has_polymer = false;
135  if (phase.active(Phase::POLYMER)) {
136  pu.has_polymer = true;
137  }
138 
139  return pu;
140  }
141 
142 }
143 
144 #endif // OPM_PHASEUSAGEFROMDECK_HEADER_INCLUDED
Definition: AnisotropicEikonal.cpp:446