27 #ifndef OPM_SOLVENT_PVT_HPP
28 #define OPM_SOLVENT_PVT_HPP
33 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
34 #include <opm/input/eclipse/Schedule/Schedule.hpp>
35 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
36 #include <opm/input/eclipse/EclipseState/Tables/PvdsTable.hpp>
46 template <
class Scalar>
49 using SamplingPoints = std::vector<std::pair<Scalar, Scalar>>;
55 SolventPvt(
const std::vector<Scalar>& solventReferenceDensity,
56 const std::vector<TabulatedOneDFunction>& inverseSolventB,
57 const std::vector<TabulatedOneDFunction>& solventMu,
58 const std::vector<TabulatedOneDFunction>& inverseSolventBMu)
59 : solventReferenceDensity_(solventReferenceDensity)
60 , inverseSolventB_(inverseSolventB)
61 , solventMu_(solventMu)
62 , inverseSolventBMu_(inverseSolventBMu)
72 void initFromState(
const EclipseState& eclState,
const Schedule&)
74 const auto& pvdsTables = eclState.getTableManager().getPvdsTables();
75 const auto& sdensityTables = eclState.getTableManager().getSolventDensityTables();
77 assert(pvdsTables.size() == sdensityTables.size());
82 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
83 Scalar rhoRefS = sdensityTables[regionIdx].getSolventDensityColumn().front();
87 const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
91 std::vector<Scalar> invB(pvdsTable.numRows());
92 const auto& Bg = pvdsTable.getFormationFactorColumn();
93 for (
unsigned i = 0; i < Bg.size(); ++ i) {
97 size_t numSamples = invB.size();
98 inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
99 solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
118 { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
126 { solventMu_[regionIdx] = mug; }
135 SamplingPoints tmp(samplePoints);
136 auto it = tmp.begin();
137 const auto& endIt = tmp.end();
138 for (; it != endIt; ++ it)
139 std::get<1>(*it) = 1.0/std::get<1>(*it);
141 inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
142 assert(inverseSolventB_[regionIdx].monotonic());
152 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
155 const auto& solventMu = solventMu_[regionIdx];
156 const auto& invSolventB = inverseSolventB_[regionIdx];
157 assert(solventMu.numSamples() == invSolventB.numSamples());
159 std::vector<Scalar> pressureValues(solventMu.numSamples());
160 std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
161 for (
unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
162 pressureValues[pIdx] = invSolventB.xAt(pIdx);
163 invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
166 inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
174 {
return solventReferenceDensity_.size(); }
179 template <
class Evaluation>
182 const Evaluation& pressure)
const
184 const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure,
true);
185 const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure,
true);
187 return invBg/invMugBg;
190 template <
class Evaluation>
191 Evaluation diffusionCoefficient(
const Evaluation& ,
195 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
202 {
return solventReferenceDensity_[regionIdx]; }
207 template <
class Evaluation>
210 const Evaluation& pressure)
const
211 {
return inverseSolventB_[regionIdx].eval(pressure,
true); }
213 const std::vector<Scalar>& solventReferenceDensity()
const
214 {
return solventReferenceDensity_; }
216 const std::vector<TabulatedOneDFunction>& inverseSolventB()
const
217 {
return inverseSolventB_; }
219 const std::vector<TabulatedOneDFunction>& solventMu()
const
220 {
return solventMu_; }
222 const std::vector<TabulatedOneDFunction>& inverseSolventBMu()
const
223 {
return inverseSolventBMu_; }
225 bool operator==(
const SolventPvt<Scalar>& data)
const
227 return solventReferenceDensity_ == data.solventReferenceDensity_ &&
228 inverseSolventB_ == data.inverseSolventB_ &&
229 solventMu_ == data.solventMu_ &&
230 inverseSolventBMu_ == data.inverseSolventBMu_;
234 std::vector<Scalar> solventReferenceDensity_;
235 std::vector<TabulatedOneDFunction> inverseSolventB_;
236 std::vector<TabulatedOneDFunction> solventMu_;
237 std::vector<TabulatedOneDFunction> inverseSolventBMu_;
Implements a linearly interpolated scalar function that depends on one variable.
This class represents the Pressure-Volume-Temperature relations of the "second" gas phase in the of E...
Definition: SolventPvt.hpp:48
Evaluation viscosity(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: SolventPvt.hpp:180
void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
Initialize the reference density of the solvent gas for a given PVT region.
Definition: SolventPvt.hpp:117
void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the solvent gas phase.
Definition: SolventPvt.hpp:125
Scalar referenceDensity(unsigned regionIdx) const
Return the reference density the solvent phase for a given PVT region.
Definition: SolventPvt.hpp:201
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: SolventPvt.hpp:148
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: SolventPvt.hpp:173
void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of solvent gas.
Definition: SolventPvt.hpp:133
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: SolventPvt.hpp:208
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48