All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
FluidMatrixInteractionBlackoil.hpp
1 /*
2  Copyright 2010 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_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
21 #define OPM_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
22 
23 #include <opm/core/utility/UniformTableLinear.hpp>
24 #include <opm/core/utility/buildUniformMonotoneTable.hpp>
25 #include <opm/parser/eclipse/Deck/Deck.hpp>
26 #include <opm/parser/eclipse/Parser/ParseContext.hpp>
27 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
28 #include <opm/parser/eclipse/EclipseState/Tables/SwofTable.hpp>
29 #include <opm/parser/eclipse/EclipseState/Tables/SgofTable.hpp>
30 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
31 #include "BlackoilDefs.hpp"
32 #include <iostream>
33 #include <stdexcept>
34 
35 namespace Opm
36 {
37 
38 // Forward declaration needed by associated parameters class.
39 template <class ScalarT, class ParamsT>
41 
42 template <class ScalarT>
44 {
45 public:
46  typedef ScalarT Scalar;
47  void init(const Opm::Deck& deck)
48  {
49  ParseContext parseContext;
50  EclipseState eclipseState(deck , parseContext);
51  // Extract input data.
52  const auto& tables = eclipseState.getTableManager();
53  const auto& swofTable = tables.getSwofTables().getTable<SwofTable>(0);
54  const auto& sgofTable = tables.getSgofTables().getTable<SgofTable>(0);
55 
56  std::vector<double> sw = swofTable.getSwColumn().vectorCopy();
57  std::vector<double> krw = swofTable.getKrwColumn().vectorCopy();
58  std::vector<double> krow = swofTable.getKrowColumn().vectorCopy();
59  std::vector<double> pcow = swofTable.getPcowColumn().vectorCopy();
60  std::vector<double> sg = sgofTable.getSgColumn().vectorCopy();
61  std::vector<double> krg = sgofTable.getKrgColumn().vectorCopy();
62  std::vector<double> krog = sgofTable.getKrogColumn().vectorCopy();
63  std::vector<double> pcog = sgofTable.getPcogColumn().vectorCopy();
64 
65  // Create tables for krw, krow, krg and krog.
66  int samples = 200;
67  buildUniformMonotoneTable(sw, krw, samples, krw_);
68  buildUniformMonotoneTable(sw, krow, samples, krow_);
69  buildUniformMonotoneTable(sg, krg, samples, krg_);
70  buildUniformMonotoneTable(sg, krog, samples, krog_);
71  krocw_ = krow[0]; // At connate water -> ecl. SWOF
72 
73  // Create tables for pcow and pcog.
74  buildUniformMonotoneTable(sw, pcow, samples, pcow_);
75  buildUniformMonotoneTable(sg, pcog, samples, pcog_);
76  }
77 
78 private:
79  template <class S, class P>
80  friend class FluidMatrixInteractionBlackoil;
81 
82  Opm::utils::UniformTableLinear<Scalar> krw_;
83  Opm::utils::UniformTableLinear<Scalar> krow_;
84  Opm::utils::UniformTableLinear<Scalar> pcow_;
85  Opm::utils::UniformTableLinear<Scalar> krg_;
86  Opm::utils::UniformTableLinear<Scalar> krog_;
87  Opm::utils::UniformTableLinear<Scalar> pcog_;
88  Scalar krocw_; // = krow_(s_wc)
89 };
90 
91 
97 template <class ScalarT, class ParamsT = FluidMatrixInteractionBlackoilParams<ScalarT> >
99 {
100 public:
101  typedef ParamsT Params;
102  typedef typename Params::Scalar Scalar;
103 
114  template <class pcContainerT, class SatContainerT>
115  static void pC(pcContainerT &pc,
116  const Params &params,
117  const SatContainerT &saturations,
118  Scalar /*temperature*/)
119  {
120  Scalar sw = saturations[Aqua];
121  Scalar sg = saturations[Vapour];
122  pc[Liquid] = 0.0;
123  pc[Aqua] = params.pcow_(sw);
124  pc[Vapour] = params.pcog_(sg);
125  }
126 
127 
131  template <class krContainerT, class SatContainerT>
132  static void kr(krContainerT &kr,
133  const Params &params,
134  const SatContainerT &saturations,
135  Scalar /*temperature*/)
136  {
137  // Stone-II relative permeability model.
138  Scalar sw = saturations[Aqua];
139  Scalar sg = saturations[Vapour];
140  Scalar krw = params.krw_(sw);
141  Scalar krg = params.krg_(sg);
142  Scalar krow = params.krow_(sw);
143  Scalar krog = params.krog_(sg);
144  Scalar krocw = params.krocw_;
145  kr[Aqua] = krw;
146  kr[Vapour] = krg;
147  kr[Liquid] = krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
148  if (kr[Liquid] < 0.0) {
149  kr[Liquid] = 0.0;
150  }
151  }
152 
153 
157  template <class krContainerT, class SatContainerT>
158  static void dkr(krContainerT &dkr,
159  const Params &params,
160  const SatContainerT &saturations,
161  Scalar /*temperature*/)
162  {
163  for (int p1 = 0; p1 < numPhases; ++p1) {
164  for (int p2 = 0; p2 < numPhases; ++p2) {
165  dkr[p1][p2] = 0.0;
166  }
167  }
168  // Stone-II relative permeability model.
169  Scalar sw = saturations[Aqua];
170  Scalar sg = saturations[Vapour];
171  Scalar krw = params.krw_(sw);
172  Scalar dkrww = params.krw_.derivative(sw);
173  Scalar krg = params.krg_(sg);
174  Scalar dkrgg = params.krg_.derivative(sg);
175  Scalar krow = params.krow_(sw);
176  Scalar dkrow = params.krow_.derivative(sw);
177  Scalar krog = params.krog_(sg);
178  Scalar dkrog = params.krog_.derivative(sg);
179  Scalar krocw = params.krocw_;
180  dkr[Aqua][Aqua] = dkrww;
181  dkr[Vapour][Vapour] = dkrgg;
182  dkr[Liquid][Aqua] = krocw*((dkrow/krocw + dkrww)*(krog/krocw + krg) - dkrww);
183  dkr[Liquid][Vapour] = krocw*((krow/krocw + krw)*(dkrog/krocw + dkrgg) - dkrgg);
184  }
185 };
186 
187 } // namespace Opm
188 
189 
190 
191 
192 #endif // OPM_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
Definition: BlackoilDefs.hpp:33
static void dkr(krContainerT &dkr, const Params &params, const SatContainerT &saturations, Scalar)
The saturation derivatives of relative permeability of all phases.
Definition: FluidMatrixInteractionBlackoil.hpp:158
static void pC(pcContainerT &pc, const Params &params, const SatContainerT &saturations, Scalar)
The linear capillary pressure-saturation curve.
Definition: FluidMatrixInteractionBlackoil.hpp:115
static void kr(krContainerT &kr, const Params &params, const SatContainerT &saturations, Scalar)
The relative permeability of all phases.
Definition: FluidMatrixInteractionBlackoil.hpp:132
Capillary pressures and relative permeabilities for a black oil system.
Definition: FluidMatrixInteractionBlackoil.hpp:40
Definition: FluidMatrixInteractionBlackoil.hpp:43