27 #ifndef OPM_SOLVENT_PVT_HPP
28 #define OPM_SOLVENT_PVT_HPP
35 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
36 #include <opm/input/eclipse/Schedule/Schedule.hpp>
37 #include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
38 #include <opm/input/eclipse/EclipseState/Tables/PvdsTable.hpp>
48 template <
class Scalar>
51 typedef std::vector<std::pair<Scalar, Scalar> > SamplingPoints;
57 SolventPvt(
const std::vector<Scalar>& solventReferenceDensity,
58 const std::vector<TabulatedOneDFunction>& inverseSolventB,
59 const std::vector<TabulatedOneDFunction>& solventMu,
60 const std::vector<TabulatedOneDFunction>& inverseSolventBMu)
61 : solventReferenceDensity_(solventReferenceDensity)
62 , inverseSolventB_(inverseSolventB)
63 , solventMu_(solventMu)
64 , inverseSolventBMu_(inverseSolventBMu)
74 void initFromState(
const EclipseState& eclState,
const Schedule&)
76 const auto& pvdsTables = eclState.getTableManager().getPvdsTables();
77 const auto& sdensityTables = eclState.getTableManager().getSolventDensityTables();
79 assert(pvdsTables.size() == sdensityTables.size());
84 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
85 Scalar rhoRefS = sdensityTables[regionIdx].getSolventDensityColumn().front();
89 const auto& pvdsTable = pvdsTables.getTable<PvdsTable>(regionIdx);
93 std::vector<Scalar> invB(pvdsTable.numRows());
94 const auto& Bg = pvdsTable.getFormationFactorColumn();
95 for (
unsigned i = 0; i < Bg.size(); ++ i) {
99 size_t numSamples = invB.size();
100 inverseSolventB_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), invB);
101 solventMu_[regionIdx].setXYArrays(numSamples, pvdsTable.getPressureColumn(), pvdsTable.getViscosityColumn());
120 { solventReferenceDensity_[regionIdx] = rhoRefSolvent; }
128 { solventMu_[regionIdx] = mug; }
137 SamplingPoints tmp(samplePoints);
138 auto it = tmp.begin();
139 const auto& endIt = tmp.end();
140 for (; it != endIt; ++ it)
141 std::get<1>(*it) = 1.0/std::get<1>(*it);
143 inverseSolventB_[regionIdx].setContainerOfTuples(tmp);
144 assert(inverseSolventB_[regionIdx].monotonic());
154 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
157 const auto& solventMu = solventMu_[regionIdx];
158 const auto& invSolventB = inverseSolventB_[regionIdx];
159 assert(solventMu.numSamples() == invSolventB.numSamples());
161 std::vector<Scalar> pressureValues(solventMu.numSamples());
162 std::vector<Scalar> invSolventBMuValues(solventMu.numSamples());
163 for (
unsigned pIdx = 0; pIdx < solventMu.numSamples(); ++pIdx) {
164 pressureValues[pIdx] = invSolventB.xAt(pIdx);
165 invSolventBMuValues[pIdx] = invSolventB.valueAt(pIdx) * (1.0/solventMu.valueAt(pIdx));
168 inverseSolventBMu_[regionIdx].setXYContainers(pressureValues, invSolventBMuValues);
176 {
return solventReferenceDensity_.size(); }
181 template <
class Evaluation>
184 const Evaluation& pressure)
const
186 const Evaluation& invBg = inverseSolventB_[regionIdx].eval(pressure,
true);
187 const Evaluation& invMugBg = inverseSolventBMu_[regionIdx].eval(pressure,
true);
189 return invBg/invMugBg;
192 template <
class Evaluation>
193 Evaluation diffusionCoefficient(
const Evaluation& ,
197 throw std::runtime_error(
"Not implemented: The PVT model does not provide a diffusionCoefficient()");
204 {
return solventReferenceDensity_[regionIdx]; }
209 template <
class Evaluation>
212 const Evaluation& pressure)
const
213 {
return inverseSolventB_[regionIdx].eval(pressure,
true); }
215 const std::vector<Scalar>& solventReferenceDensity()
const
216 {
return solventReferenceDensity_; }
218 const std::vector<TabulatedOneDFunction>& inverseSolventB()
const
219 {
return inverseSolventB_; }
221 const std::vector<TabulatedOneDFunction>& solventMu()
const
222 {
return solventMu_; }
224 const std::vector<TabulatedOneDFunction>& inverseSolventBMu()
const
225 {
return inverseSolventBMu_; }
227 bool operator==(
const SolventPvt<Scalar>& data)
const
229 return solventReferenceDensity_ == data.solventReferenceDensity_ &&
230 inverseSolventB_ == data.inverseSolventB_ &&
231 solventMu_ == data.solventMu_ &&
232 inverseSolventBMu_ == data.inverseSolventBMu_;
236 std::vector<Scalar> solventReferenceDensity_;
237 std::vector<TabulatedOneDFunction> inverseSolventB_;
238 std::vector<TabulatedOneDFunction> solventMu_;
239 std::vector<TabulatedOneDFunction> inverseSolventBMu_;
A central place for various physical constants occuring in some equations.
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:50
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:182
void setReferenceDensity(unsigned regionIdx, Scalar rhoRefSolvent)
Initialize the reference density of the solvent gas for a given PVT region.
Definition: SolventPvt.hpp:119
void setSolventViscosity(unsigned regionIdx, const TabulatedOneDFunction &mug)
Initialize the viscosity of the solvent gas phase.
Definition: SolventPvt.hpp:127
Scalar referenceDensity(unsigned regionIdx) const
Return the reference density the solvent phase for a given PVT region.
Definition: SolventPvt.hpp:203
void initEnd()
Finish initializing the oil phase PVT properties.
Definition: SolventPvt.hpp:150
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: SolventPvt.hpp:175
void setSolventFormationVolumeFactor(unsigned regionIdx, const SamplingPoints &samplePoints)
Initialize the function for the formation volume factor of solvent gas.
Definition: SolventPvt.hpp:135
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure) const
Returns the formation volume factor [-] of the fluid phase.
Definition: SolventPvt.hpp:210
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:47