27 #ifndef OPM_CONSTANT_COMPRESSIBILITY_BRINE_PVT_HPP
28 #define OPM_CONSTANT_COMPRESSIBILITY_BRINE_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/PvtwsaltTable.hpp>
41 template <
class Scalar,
bool enableThermal,
bool enableBrine>
42 class WaterPvtMultiplexer;
47 template <
class Scalar>
55 const std::vector<Scalar>& referencePressure,
56 const std::vector<TabulatedFunction> formationVolumeTables,
57 const std::vector<TabulatedFunction> compressibilityTables,
58 const std::vector<TabulatedFunction> viscosityTables,
59 const std::vector<TabulatedFunction> viscosibilityTables)
60 : waterReferenceDensity_(waterReferenceDensity)
61 , referencePressure_(referencePressure)
62 , formationVolumeTables_(formationVolumeTables)
63 , compressibilityTables_(compressibilityTables)
64 , viscosityTables_(viscosityTables)
65 , viscosibilityTables_(viscosibilityTables)
72 void initFromState(
const EclipseState& eclState,
const Schedule&)
74 const auto& tableManager = eclState.getTableManager();
75 size_t numRegions = tableManager.getTabdims().getNumPVTTables();
76 const auto& densityTable = tableManager.getDensityTable();
84 const auto& pvtwsaltTables = tableManager.getPvtwSaltTables();
85 if(!pvtwsaltTables.empty()){
87 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++ regionIdx) {
88 const auto& pvtwsaltTable = pvtwsaltTables[regionIdx];
89 const auto& c = pvtwsaltTable.getSaltConcentrationColumn();
91 const auto& B = pvtwsaltTable.getFormationVolumeFactorColumn();
92 formationVolumeTables_[regionIdx].setXYContainers(c, B);
94 const auto& compressibility = pvtwsaltTable.getCompressibilityColumn();
95 compressibilityTables_[regionIdx].setXYContainers(c, compressibility);
97 const auto& viscositytable = pvtwsaltTable.getViscosityColumn();
98 viscosityTables_[regionIdx].setXYContainers(c, viscositytable);
100 const auto& viscosibility = pvtwsaltTable.getViscosibilityColumn();
101 viscosibilityTables_[regionIdx].setXYContainers(c, viscosibility);
102 referencePressure_[regionIdx] = pvtwsaltTable.getReferencePressureValue();
106 throw std::runtime_error(
"PVTWSALT must be specified in BRINE runs\n");
111 setNumRegions(numPvtwRegions);
113 for (
unsigned regionIdx = 0; regionIdx < numPvtwRegions; ++ regionIdx) {
115 waterReferenceDensity_[regionIdx] = densityTable[regionIdx].water;
126 for (
unsigned regionIdx = 0; regionIdx <
numRegions; ++regionIdx) {
138 { waterReferenceDensity_[regionIdx] = rhoRefWater; }
151 {
return waterReferenceDensity_.size(); }
156 template <
class Evaluation>
160 const Evaluation&)
const
162 throw std::runtime_error(
"Requested the enthalpy of water but the thermal option is not enabled");
168 template <
class Evaluation>
170 const Evaluation& temperature,
171 const Evaluation& pressure,
172 const Evaluation& saltconcentration)
const
175 Scalar pRef = referencePressure_[regionIdx];
176 const Evaluation C = compressibilityTables_[regionIdx].eval(saltconcentration,
true);
177 const Evaluation Cv = viscosibilityTables_[regionIdx].eval(saltconcentration,
true);
178 const Evaluation BwRef = formationVolumeTables_[regionIdx].eval(saltconcentration,
true);
179 const Evaluation Y = (C-Cv)* (pressure - pRef);
180 Evaluation MuwRef = viscosityTables_[regionIdx].eval(saltconcentration,
true);
184 return MuwRef*BwRef*bw/(1 + Y*(1 + Y/2));
190 template <
class Evaluation>
193 const Evaluation& pressure,
194 const Evaluation& saltconcentration)
const
196 Scalar pRef = referencePressure_[regionIdx];
198 const Evaluation BwRef = formationVolumeTables_[regionIdx].eval(saltconcentration,
true);
199 const Evaluation C = compressibilityTables_[regionIdx].eval(saltconcentration,
true);
200 const Evaluation X = C * (pressure - pRef);
202 return (1.0 + X*(1.0 + X/2.0))/BwRef;
206 const Scalar waterReferenceDensity(
unsigned regionIdx)
const
207 {
return waterReferenceDensity_[regionIdx]; }
209 const std::vector<Scalar>& referencePressure()
const
210 {
return referencePressure_; }
212 const std::vector<TabulatedFunction>& formationVolumeTables()
const
213 {
return formationVolumeTables_; }
215 const std::vector<TabulatedFunction>& compressibilityTables()
const
216 {
return compressibilityTables_; }
218 const std::vector<TabulatedFunction>& viscosityTables()
const
219 {
return viscosityTables_; }
221 const std::vector<TabulatedFunction>& viscosibilityTables()
const
222 {
return viscosibilityTables_; }
224 bool operator==(
const ConstantCompressibilityBrinePvt<Scalar>& data)
const
226 return this->waterReferenceDensity_ == data.waterReferenceDensity_ &&
227 this->referencePressure() == data.referencePressure() &&
228 this->formationVolumeTables() == data.formationVolumeTables() &&
229 this->compressibilityTables() == data.compressibilityTables() &&
230 this->viscosityTables() == data.viscosityTables() &&
231 this->viscosibilityTables() == data.viscosibilityTables();
235 std::vector<Scalar> waterReferenceDensity_;
236 std::vector<Scalar> referencePressure_;
237 std::vector<TabulatedFunction> formationVolumeTables_;
238 std::vector<TabulatedFunction> compressibilityTables_;
239 std::vector<TabulatedFunction> viscosityTables_;
240 std::vector<TabulatedFunction> viscosibilityTables_;
Implements a linearly interpolated scalar function that depends on one variable.
This class represents the Pressure-Volume-Temperature relations of the gas phase without vaporized oi...
Definition: ConstantCompressibilityBrinePvt.hpp:49
unsigned numRegions() const
Return the number of PVT regions which are considered by this PVT-object.
Definition: ConstantCompressibilityBrinePvt.hpp:150
Evaluation inverseFormationVolumeFactor(unsigned regionIdx, const Evaluation &, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the formation volume factor [-] of the fluid phase.
Definition: ConstantCompressibilityBrinePvt.hpp:191
Evaluation internalEnergy(unsigned, const Evaluation &, const Evaluation &, const Evaluation &) const
Returns the specific enthalpy [J/kg] of water given a set of parameters.
Definition: ConstantCompressibilityBrinePvt.hpp:157
void setReferenceDensities(unsigned regionIdx, Scalar, Scalar, Scalar rhoRefWater)
Set the water reference density [kg / m^3].
Definition: ConstantCompressibilityBrinePvt.hpp:134
void initEnd()
Finish initializing the water phase PVT properties.
Definition: ConstantCompressibilityBrinePvt.hpp:143
Evaluation viscosity(unsigned regionIdx, const Evaluation &temperature, const Evaluation &pressure, const Evaluation &saltconcentration) const
Returns the dynamic viscosity [Pa s] of the fluid phase given a set of parameters.
Definition: ConstantCompressibilityBrinePvt.hpp:169
Implements a linearly interpolated scalar function that depends on one variable.
Definition: Tabulated1DFunction.hpp:48