20 #ifndef OPM_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
21 #define OPM_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
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"
39 template <
class ScalarT,
class ParamsT>
42 template <
class ScalarT>
46 typedef ScalarT Scalar;
47 void init(
const Opm::Deck& deck)
49 ParseContext parseContext;
50 EclipseState eclipseState(deck , parseContext);
52 const auto& tables = eclipseState.getTableManager();
53 const auto& swofTable = tables.getSwofTables().getTable<SwofTable>(0);
54 const auto& sgofTable = tables.getSgofTables().getTable<SgofTable>(0);
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();
67 buildUniformMonotoneTable(sw, krw, samples, krw_);
68 buildUniformMonotoneTable(sw, krow, samples, krow_);
69 buildUniformMonotoneTable(sg, krg, samples, krg_);
70 buildUniformMonotoneTable(sg, krog, samples, krog_);
74 buildUniformMonotoneTable(sw, pcow, samples, pcow_);
75 buildUniformMonotoneTable(sg, pcog, samples, pcog_);
79 template <
class S,
class P>
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_;
97 template <
class ScalarT,
class ParamsT = Flu
idMatrixInteractionBlackoilParams<ScalarT> >
101 typedef ParamsT Params;
102 typedef typename Params::Scalar Scalar;
114 template <
class pcContainerT,
class SatContainerT>
115 static void pC(pcContainerT &pc,
116 const Params ¶ms,
117 const SatContainerT &saturations,
120 Scalar sw = saturations[Aqua];
121 Scalar sg = saturations[Vapour];
123 pc[Aqua] = params.pcow_(sw);
124 pc[Vapour] = params.pcog_(sg);
131 template <
class krContainerT,
class SatContainerT>
132 static void kr(krContainerT &
kr,
133 const Params ¶ms,
134 const SatContainerT &saturations,
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_;
147 kr[Liquid] = krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
148 if (kr[Liquid] < 0.0) {
157 template <
class krContainerT,
class SatContainerT>
159 const Params ¶ms,
160 const SatContainerT &saturations,
163 for (
int p1 = 0; p1 < numPhases; ++p1) {
164 for (
int p2 = 0; p2 < numPhases; ++p2) {
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);
192 #endif // OPM_FLUIDMATRIXINTERACTIONBLACKOIL_HEADER_INCLUDED
Definition: BlackoilDefs.hpp:33
static void dkr(krContainerT &dkr, const Params ¶ms, const SatContainerT &saturations, Scalar)
The saturation derivatives of relative permeability of all phases.
Definition: FluidMatrixInteractionBlackoil.hpp:158
static void pC(pcContainerT &pc, const Params ¶ms, const SatContainerT &saturations, Scalar)
The linear capillary pressure-saturation curve.
Definition: FluidMatrixInteractionBlackoil.hpp:115
static void kr(krContainerT &kr, const Params ¶ms, 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