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>
39template <
class ScalarT,
class ParamsT>
40class FluidMatrixInteractionBlackoil;
42template <
class ScalarT>
47 void init(Opm::DeckConstPtr 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_;
97template <
class ScalarT,
class ParamsT = Flu
idMatrixInteractionBlackoilParams<ScalarT> >
114 template <
class pcContainerT,
class SatContainerT>
115 static void pC(pcContainerT &pc,
117 const SatContainerT &saturations,
123 pc[
Aqua] = params.pcow_(sw);
124 pc[
Vapour] = params.pcog_(sg);
138 template <
class SatContainerT,
class pcContainerT>
139 static void S(SatContainerT &saturations,
141 const pcContainerT &pc,
144 std::cerr <<
"FluidMatrixInteractionBlackoil::S() is not implemented yet\n";
145 throw std::logic_error(
"Not implemented");
152 template <
class krContainerT,
class SatContainerT>
153 static void kr(krContainerT &
kr,
155 const SatContainerT &saturations,
161 Scalar krw = params.krw_(sw);
162 Scalar krg = params.krg_(sg);
163 Scalar krow = params.krow_(sw);
164 Scalar krog = params.krog_(sg);
165 Scalar krocw = params.krocw_;
168 kr[
Liquid] = krocw*((krow/krocw + krw)*(krog/krocw + krg) - krw - krg);
178 template <
class krContainerT,
class SatContainerT>
181 const SatContainerT &saturations,
192 Scalar krw = params.krw_(sw);
193 Scalar dkrww = params.krw_.derivative(sw);
194 Scalar krg = params.krg_(sg);
195 Scalar dkrgg = params.krg_.derivative(sg);
196 Scalar krow = params.krow_(sw);
197 Scalar dkrow = params.krow_.derivative(sw);
198 Scalar krog = params.krog_(sg);
199 Scalar dkrog = params.krog_.derivative(sg);
200 Scalar krocw = params.krocw_;
203 dkr[
Liquid][
Aqua] = krocw*((dkrow/krocw + dkrww)*(krog/krocw + krg) - dkrww);
204 dkr[
Liquid][
Vapour] = krocw*((krow/krocw + krw)*(dkrog/krocw + dkrgg) - dkrgg);
Definition: BlackoilDefs.hpp:31
double Scalar
Definition: BlackoilDefs.hpp:39
@ Aqua
Definition: BlackoilDefs.hpp:37
@ Vapour
Definition: BlackoilDefs.hpp:37
@ Liquid
Definition: BlackoilDefs.hpp:37
@ numPhases
Definition: BlackoilDefs.hpp:34
Definition: FluidMatrixInteractionBlackoil.hpp:44
void init(Opm::DeckConstPtr deck)
Definition: FluidMatrixInteractionBlackoil.hpp:47
ScalarT Scalar
Definition: FluidMatrixInteractionBlackoil.hpp:46
Capillary pressures and relative permeabilities for a black oil system.
Definition: FluidMatrixInteractionBlackoil.hpp:99
static void kr(krContainerT &kr, const Params ¶ms, const SatContainerT &saturations, Scalar)
The relative permeability of all phases.
Definition: FluidMatrixInteractionBlackoil.hpp:153
static void dkr(krContainerT &dkr, const Params ¶ms, const SatContainerT &saturations, Scalar)
The saturation derivatives of relative permeability of all phases.
Definition: FluidMatrixInteractionBlackoil.hpp:179
Params::Scalar Scalar
Definition: FluidMatrixInteractionBlackoil.hpp:102
static void pC(pcContainerT &pc, const Params ¶ms, const SatContainerT &saturations, Scalar)
The linear capillary pressure-saturation curve.
Definition: FluidMatrixInteractionBlackoil.hpp:115
static void S(SatContainerT &saturations, const Params ¶ms, const pcContainerT &pc, Scalar)
The saturation-capillary pressure curve.
Definition: FluidMatrixInteractionBlackoil.hpp:139
ParamsT Params
Definition: FluidMatrixInteractionBlackoil.hpp:101
Definition: BlackoilFluid.hpp:32