27 #ifndef OPM_SOLVENT_PVT_HPP
28 #define OPM_SOLVENT_PVT_HPP
35 #include <opm/parser/eclipse/Deck/Deck.hpp>
36 #include <opm/parser/eclipse/EclipseState/EclipseState.hpp>
37 #include <opm/parser/eclipse/EclipseState/Tables/TableManager.hpp>
38 #include <opm/parser/eclipse/EclipseState/Tables/PvdsTable.hpp>
39 #include <opm/parser/eclipse/Deck/DeckKeyword.hpp>
40 #include <opm/parser/eclipse/Deck/DeckRecord.hpp>
50 template <
class Scalar>
54 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
63 void initFromDeck(
const Deck& deck,
const EclipseState& eclState)
65 const auto& pvdsTables = eclState.getTableManager().getPvdsTables();
66 const auto& sdensityKeyword = deck.getKeyword(
"SDENSITY");
68 assert(pvdsTables.size() == sdensityKeyword.size());
71 setNumRegions(numRegions);
73 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
74 Scalar rhoRefS = sdensityKeyword.getRecord(regionIdx).getItem(
"SOLVENT_DENSITY").getSIDouble(0);
78 const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
82 std::vector<Scalar> invB(pvdsTable.numRows());
83 const auto& Bg = pvdsTable.getFormationFactorColumn();
84 for (
unsigned i = 0; i < Bg.size(); ++ i) {
88 size_t numSamples = invB.size();
89 inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
90 solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
99 solventReferenceDensity_.resize(numRegions);
100 inverseSolventB_.resize(numRegions);
101 inverseSolventBMu_.resize(numRegions);
102 solventMu_.resize(numRegions);
109 { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
117 { solventMu_[regionIdx] = mug; }
126 SamplingPoints tmp(samplePoints);
127 auto it = tmp.begin();
128 const auto& endIt = tmp.end();
129 for (; it != endIt; ++ it)
130 std::get<1>(*it) = 1.0/std::get<1>(*it);
132 inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
133 assert(inverseSolventB_[regionIdx].monotonic());
142 size_t numRegions = solventMu_.size();
143 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
146 const auto& solventMu = solventMu_[regionIdx];
147 const auto& invSolventB = inverseSolventB_[regionIdx];
148 assert(solventMu.numSamples() == invSolventB.numSamples());
150 std::vector<Scalar> pressureValues(solventMu.numSamples());
151 std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
152 for (
unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
153 pressureValues[pIdx] = invSolventB.xAt(pIdx);
154 invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
157 inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
165 {
return solventReferenceDensity_.size(); }
170 template <
class Evaluation>
172 const Evaluation& temperature OPM_UNUSED,
173 const Evaluation& pressure)
const
175 const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure,
true);
176 const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure,
true);
178 return invBg/invMugBg;
185 {
return solventReferenceDensity_[regionIdx]; }
190 template <
class Evaluation>
192 const Evaluation& temperature OPM_UNUSED,
193 const Evaluation& pressure)
const
194 {
return inverseSolventB_[regionIdx].eval(pressure,
true); }
197 std::vector<Scalar> solventReferenceDensity_;
198 std::vector<TabulatedOneDFunction> inverseSolventB_;
199 std::vector<TabulatedOneDFunction> solventMu_;
200 std::vector<TabulatedOneDFunction> inverseSolventBMu_;
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: SolventPvt.hpp:164
This class represents the Pressure-Volume-Temperature relations of the "second" gas phase in the of E...
Definition: SolventPvt.hpp:51
void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
Initialize the reference density of the solvent gas for a given PVT region.
Definition: SolventPvt.hpp:108
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: SolventPvt.hpp:139
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature OPM_UNUSED, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: SolventPvt.hpp:171
Scalar referenceDensity(unsigned regionIdx) const
Return the reference density the solvent phase for a given PVT region.
Definition: SolventPvt.hpp:184
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &temperature OPM_UNUSED, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: SolventPvt.hpp:191
void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the solvent gas phase.
Definition: SolventPvt.hpp:116
Implements a linearly interpolated scalar function that depends on one variable.
A central place for various physical constants occuring in some equations.
void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of solvent gas.
Definition: SolventPvt.hpp:124
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47